FORTRAN Routines Research Articles

PSGen, a generator of phase space parameterizations for the multichannel Monte Carlo integration

PSGen is a new general purpose Fortran program which has been written to facilitate the Monte Carlo phase space integration of the S matrix element of any 2→n scattering process, with n=2,...,9, provided by the user. The program is written in Fortran 90/95. It uses a new very fast algorithm that automatically generates calls to Fortran subroutines containing different phase space parameterizations of the considered class of processes. The parameterizations take into account mappings of poles due to the Feynman propagators of unstable heavy particles decaying into 2 or 3 on shell final state particles according to predefined patterns, possible single or double t-channel poles and peaks due to one on shell photon or gluon radiation. The individual subroutines are organized in a single multichannel kinematics subroutine which can be easily called while computing the phase space integral of the S matrix element as a function of generated particle four momenta, in either the leading or higher orders of the perturbation series. The particle four momenta can be used in a quadruple precision version, if necessary. Program summaryProgram Title:PSGenCPC Library link to program files:https://doi.org/10.17632/y7nprysy8g.1Developer's repository link:https://www.kk.us.edu.pl/psgen.htmlLicensing provisions: GPLv3Programming language: Fortran 90/95Nature of problem: Automatic generation of Fortran routines which allow for the multichannel Monte Carlo integration of user's provided S matrix element as a function of generated particle four momenta of any 2→n, with n=2,...,9, scattering process. The S matrix can be defined in either the leading or higher orders.Solution method: A new very fast algorithm written in Fortran 90/95 generates calls to different handwritten subroutines containing phase space parameterizations which, among others, take into account poles due to the Feynman propagators of unstable heavy particles decaying into 2 or 3 on shell final state particles according to predefined patterns, possible single or double t-channel poles and peaks due to one on shell photon or gluon radiation. The generated routines are moved to a directory, where they are combined in a single properly normalized multichannel phase space parameterization. The corresponding kinematics subroutine is self-consistent and can be readily combined with any program which calculates the S matrix element in either the leading or higher orders of the perturbation series. The particle four momenta can be generated in a quadruple precision version, if necessary. The automatically generated subroutine containing calls to different kinematics subroutines can be easily supplemented with calls to user's own made kinematics subroutines.

Efficient and reproducible research in quantum open systems using Maxima with a Fortran interface

Abstract maxqosim is a software package that facilitates the analysis and simulation of quantum open systems, providing a framework where theoretical calculations, efficient numerical calculations and exposition can be coordinated and automated. maxqosim consists of Maxima and Fortran routines, revolving around the symbolic manipulation of equations relevant to quantum open systems theory and their efficient numerical simulation. maxqosim package aims to be a tool for reproducible, extendable and efficient research in quantum open systems.

Prototyping a numerical model coupled with remote sensing for tracking harmful algal blooms in shallow lakes

In this paper, we present the design, development and test cases of LakeM, a software package for modeling the transport pathways of harmful algal blooms in shallow lakes. LakeM integrates implementation of numerical methods and data visualization for the purpose of facilitating the communication between modelers and decision makers. By making use of Finite Difference Method, LakeM is able to numerically solve the Shallow Water Equations on a staggered grid. The core of LakeM consists of two modules, the hydrodynamic and tracer module, which are run separately. The hydrodynamic module simulates the surface level and water movement, and the tracer module simulates the movement of harmful algal blooms subject to advection and dispersion processes. Most of the LakeM programs are developed based on the Python language in an object-oriented fashion, which can be readily adapted to changing requirements. Computationally intensive programs are written in Fortran routines for interoperating with Python. Model validation is performed by comparing model predictions to analytical solutions. Additionally, the capability and effectiveness of LakeM are illustrated by a case study of modeling algal bloom transport in Taihu Lake. The results have shown good agreements with analytical solutions and available remote sensing images. The model is efficient, robust, and can be used to aid environmental managers in developing emergency response actions for HABs. Besides that, the development of LakeM offers a design pattern for future model development efforts in tracking HABs.

Derivation of Entropy and Exergy Transport Equations, and Application to Second Law Analysis of Sugarcane Bagasse Gasification in Bubbling Fluidized Beds

Abstract In the present work, the transport equations for mass, momentum, energy, and chemical species as given by the Euler–Euler formulation for multiphase flows are used together with the second law of thermodynamics to derive the entropy and exergy transport equations, suitable to the study of gas-particle reactive flows, such as those observed during pyrolysis, gasification, and combustion of biomass particles. The terms of the derived equations are discussed, and the exergy destruction contributions are identified. Subsequently, a kinetic model is implemented in a computational fluid dynamics (CFD) open source code for the sugarcane bagasse gasification. Then, the derived exergy destruction terms are implemented numerically through user-defined Fortran routines. Next, the second law analysis of the gasification process of sugarcane bagasse in bubbling fluidized beds is carried out. Detailed results are obtained for the local destructions of exergy along the reactor. This information is important to help improve environmental and sustainable practices and should be of interest to both designers and operators of fluidized bed equipment.

4D multi-physical stress modelling during shale gas production: A case study of Sichuan Basin shale gas reservoir, China

To investigate the time-lapse, three-dimensional (4D) stress during shale gas production, a multi-physical modeling method is proposed. This model enfranchises multiple years' scale dynamic production data. A finite difference method (FDM) reservoir simulator is used to couple thermal-hydrological-chemical (THC) processes, while a finite element method (FEM) geomechanical simulator takes on the role of a thermal-hydrological-mechanical (THM) coupling calculator. In order to couple the flow model with the geomechanical model, an interface (coupling) Python code is provided to communicate data between the finite difference (FD) and finite element (FE) grids. A set of FORTRAN routines are assigned to interpolate initial stress field and mechanical parameters into the geomechanical model based on field and experimental data. Ultimately, a 4D multi-physical model of Sichuan Basin shale gas reservoir (China) is proposed and analyzed. Based on the simulation results, the continuous reduction in reservoir pore pressure leads to evolving changes in the stress state and reservoir properties beyond the pre-consolidation pressure. Transversely isotropic characteristics have a considerable effect on the changes in the reservoir stresses and poro-elastic parameters. The results provide some fundamental insights for optimizing well patterns, ensuring wellbore stability, re-fracturing design and locating, drilling and completing infill wells.

Study of model application for drying of pulp fruit in spouted bed with intermittent feeding and accumulation

ABSTRACTIn the present study, we applied and improved a model to describe the behavior of a spouted bed dryer with intermittent feeding for suspension drying, considering suspension accumulation inside the device. This model is important in macroscopic heat and mass balances, represented by a system of ordinary differential balances, implemented and resolved numerically using Fortran routines, in which the influence of process variables on drying dynamics were assessed. The results obtained in modeling and mathematical simulations were compared with experimental data, indicating that the model fits the process well, based on outlet air temperature and moisture.

Analytical series representing the DE431 ephemerides of the outer planets

This paper presents a new development of orbital elements for the four outer planets in analytical series. These series complete the series DEA431, for which the terms for the terrestrial planets are given in an earlier paper. As the source of planetary coordinates, the development uses the latest long-term numerical ephemerides DE431, which are valid over 30 000 yr [−13 000, 17 000]. The author's modification of the frequency analysis method is used as the expansion tool. The obtained results are compared with those given by recent analytical theories of planetary motion, namely VSOP2013 and TOP2013. The number of terms in DEA431 is about one order of magnitude lower than that in either VSOP2013 or TOP2013. The extensive analytical theories of the motion of the outer planets are more precise than the compact DEA431 series over time intervals of several hundreds of years; the three compared analytical series are of about the same accuracy over time intervals of a few thousand years; and DEA431 has a higher accuracy than both VSOP2013 and TOP2013 when representing the coordinates of the outer planets over longer time intervals, such as tens of thousands of years. Over 24 000 yr centred at epoch J2000, the maximum difference between the mean orbital longitudes of the outer planets given by the DE431 numerical ephemerides and the DEA431 analytical series is 0|${^{\prime\prime}_{.}}$|08 for Jupiter, 0|${^{\prime\prime}_{.}}$|22 for Saturn, 0|${^{\prime\prime}_{.}}$|10 for Uranus, and 0|${^{\prime\prime}_{.}}$|05 for Neptune. The coefficients of the DEA431 series for all the outer planets and the corresponding FORTRAN routines for calculating the orbital elements and rectangular coordinates of the planets are available on the SAI web-server.

Mixed Non-Parametric and Parametric Estimation Techniques in R Package etasFLP for Earthquakes' Description

etasFLP is an R package which fits an epidemic type aftershock sequence (ETAS) model to an earthquake catalog; non-parametric background seismicity can be estimated through a forward predictive likelihood approach, while parametric components of triggered seismicity are estimated through maximum likelihood; estimation steps are alternated until convergence is obtained and for each event the probability of being a background event is estimated. The package includes options which allow its wide use. Methods for plot, summary and profile are defined for the main output class object. The paper provides examples of the package's use with description of the underlying R and Fortran routines.

VESUVIO Data Analysis Goes MANTID

This paper describes ongoing efforts to implement the reduction and analysis of neutron Compton scattering data within the MANTID framework. Recently, extensive work has been carried out to integrate the bespoke data reduction and analysis routines written for VESUVIO with the MANTID framework. While the programs described in this document are designed to replicate the functionality of the Fortran and Genie routines already in use, most of them have been written from scratch and are not based on the original code base.

Comparison of three Stark problem solution techniques for the bounded case

Three methods of obtaining solutions to the Stark problem—one developed by Lantoine and Russell using Jacobi elliptic and related functions, one developed by Biscani and Izzo using Weierstrass elliptic and related functions, and one developed by Pellegrini, Russell, and Vittaldev using $$F$$ and $$G$$ Taylor series extended to the Stark problem—are compared qualitatively and quantitatively for the bounded motion case. For consistency with existing available code for the series solution, Fortran routines of the Lantoine method and Biscani method are newly implemented and made available. For these implementations, the Lantoine formulation is found to be more efficient than the Biscani formulation in the propagation of a single trajectory segment. However, for applications for which acceptable accuracy may be achieved by orders up to 16, the Pellegrini series solution is shown to be more efficient than either analytical method. The three methods are also compared in the propagation of sequentially connected trajectory segments in a low-thrust orbital transfer maneuver. Separate tests are conducted for discretizations between 8 and 96 segments per orbit. For the series solution, the interaction between order and step size leads to computation times that are nearly invariable to discretization for a given truncation error tolerance over the tested range of discretizations. This finding makes the series solution particularly attractive for mission design applications where problems may require both coarse and fine discretizations. Example applications include the modeling of low-thrust propulsion and time-varying perturbations—problems for which the efficient propagation of relatively short Stark segments is paramount because the disturbing acceleration generally varies continuously.

Genome-wide regression and prediction with the BGLR statistical package.

Many modern genomic data analyses require implementing regressions where the number of parameters (p, e.g., the number of marker effects) exceeds sample size (n). Implementing these large-p-with-small-n regressions poses several statistical and computational challenges, some of which can be confronted using Bayesian methods. This approach allows integrating various parametric and nonparametric shrinkage and variable selection procedures in a unified and consistent manner. The BGLR R-package implements a large collection of Bayesian regression models, including parametric variable selection and shrinkage methods and semiparametric procedures (Bayesian reproducing kernel Hilbert spaces regressions, RKHS). The software was originally developed for genomic applications; however, the methods implemented are useful for many nongenomic applications as well. The response can be continuous (censored or not) or categorical (either binary or ordinal). The algorithm is based on a Gibbs sampler with scalar updates and the implementation takes advantage of efficient compiled C and Fortran routines. In this article we describe the methods implemented in BGLR, present examples of the use of the package, and discuss practical issues emerging in real-data analysis.

PySSM: APythonModule for Bayesian Inference of Linear Gaussian State Space Models

PySSM is a Python package that has been developed for the analysis of time series using linear Gaussian state space models. PySSM is easy to use; models can be set up quickly and efficiently and a variety of different settings are available to the user. It also takes advantage of scientific libraries NumPy and SciPy and other high level features of the Python language. PySSM is also used as a platform for interfacing between optimized and parallelized Fortran routines. These Fortran routines heavily utilize basic linear algebra and linear algebra Package functions for maximum performance. PySSM contains classes for filtering, classical smoothing as well as simulation smoothing.

Cts: AnRPackage for Continuous Time Autoregressive Models via Kalman Filter

We describe an R package cts for fitting a modified form of continuous time autoregressive model, which can be particularly useful with unequally sampled time series. The estimation is based on the application of the Kalman filter. The paper provides the methods and algorithms implemented in the package, including parameter estimation, spectral analysis, forecasting, model checking and Kalman smoothing. The package contains R functions which interface underlying Fortran routines. The package is applied to geophysical and medical data for illustration.

On the Calculation of Displacement, Stress, and Strain Induced by Triangular Dislocations

Abstract Integral equation‐based methods in elasticity require the calculation of integrals that involve a singular (or weakly singular) matrix Green’s function and a traction or displacement vector field defined on a surface. A variety of numerical methods based on analytic quadratures have been developed for such calculations, using either rectangular or triangular surface patches. Unfortunately, without correction, these analytic rules are subject to numerical instabilities in certain parameter regimes. In this paper, we present a stable, semi‐analytic collection of quadrature rules that can be applied to both infinite medium and half‐space simulations. We describe our underlying approach and illustrate its performance with numerical examples. Online Material: Fortran routines.

A Mechanical Characterization of the Porcine Atria at the Healthy Stage and After Ventricular Tachypacing

Atrial fibrillation (AF) is a cardiac arrhythmia that highly increases the risk of stroke and is associated with significant but still unexplored changes in the mechanical behavior of the tissue. Planar biaxial tests were performed on tissue specimens from pigs at the healthy stage and after ventricular tachypacing (VTP), a procedure applied to reproduce the relevant features of AF. The local arrangement of the fiber bundles in the tissue was investigated on specimens from rabbit atria by means of circularly polarized light. Based on this, mechanical data were fitted to two anisotropic constitutive relationships, including a four-parameter Fung-type model and a microstructurally-motivated model. Accounting for the fiber-induced anisotropy brought average R(2) = 0.807 for the microstructurally-motivated model and average R(2) = 0.949 for the Fung model. Validation of the fitted constitutive relationships was performed by means of FEM simulations coupled to FORTRAN routines. The performances of the two material models in predicting the second Piola-Kirchhoff stress were comparable, with average errors <3.1%. However, the Fung model outperformed the other in the prediction of the Green-Lagrange strain, with 9.2% maximum average error. To increase model generality, a proper averaging procedure accounting for nonlinearities was used to obtain average material parameters. In general, a stiffer behavior after VTP was noted.

ANALYTICAL EXPRESSIONS FOR THE ENVELOPE BINDING ENERGY OF GIANTS AS A FUNCTION OF BASIC STELLAR PARAMETERS

The common-envelope (CE) phase is an important stage in the evolution of binary stellar populations. The most common way to compute the change in orbital period during a CE is to relate the binding energy of the envelope of the Roche-lobe filling giant to the change in orbital energy. Especially in population-synthesis codes, where the evolution of millions of stars must be computed and detailed evolutionary models are too expensive computationally, simple approximations are made for the envelope binding energy. In this study, we present accurate analytic prescriptions based on detailed stellar-evolution models that provide the envelope binding energy for giants with metallicities between Z = 10-4 and Z = 0.03 and masses between 0.8 Msun and 100 Msun, as a function of the metallicity, mass, radius and evolutionary phase of the star. Our results are also presented in the form of electronic data tables and Fortran routines that use them. We find that the accuracy of our fits is better than 15% for 90% of our model data points in all cases, and better than 10% for 90% of our data points in all cases except the asymptotic giant branches for three of the six metallicities we consider. For very massive stars (M > 50 Msun), when stars lose more than ~20% of their initial mass due to stellar winds, our fits do not describe the models as accurately. Our results are more widely applicable - covering wider ranges of metallicity and mass - and are of higher accuracy than those of previous studies.

Biaxial mechanical modeling of the small intestine

Capsule endoscopes are pill-size devices provided with a camera that capture images of the small intestine from inside the body after being ingested by a patient. The interaction between intestinal tissue and capsule endoscopes needs to be investigated to optimize capsule design while preventing tissue damage. To that purpose, a constitutive model that can reliably predict the mechanical response of the intestinal tissue under complex mechanical loading is required. This paper describes the development and numerical validation of a phenomenological constitutive model for the porcine duodenum, jejunum and ileum. Parameters characterizing the mechanical behavior of the material were estimated from planar biaxial test data, where intestinal tissue specimens were simultaneously loaded along the circumferential and longitudinal directions. Specimen-specific Fung constitutive models were able to accurately predict the planar stress–strain behavior of the tested samples under a wide range of loading conditions. To increase model generality, average anisotropic constitutive relationships were also generated for each tissue region by fitting average stress–strain curves to the Fung potential. Due to the observed variability in the direction of maximum stiffness, the average Fung models were less anisotropic than the specimen-specific models. Hence, average isotropic models in the Neo–Hookean and Mooney–Rivlin forms were attempted, but they could not adequately describe the degree of nonlinearity in the tissue. Values of the R2 for the nonlinear regressions were 0.17, 0.44 and 0.93 for the average Neo–Hookean, Mooney–Rivlin and Fung models, respectively. Average models were successfully implemented into FORTRAN routines and used to simulate capsule deployment with a finite element method analysis.

TRIP

TRIP is an interactive computer algebra system that is devoted to perturbation series computations, and specially adapted to celestial mechanics. Its development started in 1988, as an upgrade of the special purpose FORTRAN routines elaborated by J. Laskar for the demonstration of the chaotic behavior of the Solar System. TRIP is now a mature and efficient tool for handling multivariate generalized power series, and contains also a numerical kernel. A limited version of TRIP is freely available on the WEB at www.imcce.fr/trip.

An efficient analytical approach for calculating line mixing in atmospheric remote sensing applications

Collisional coupling between energy states in a molecule undergoing an optical transition can alter the line shape associated with the transition, an effect known as line mixing. Accounting for this effect in the analysis of remote sensing measurements of Earth's atmosphere by the Atmospheric Chemistry Experiment (ACE) yields reduced residuals, which leads to improved performance in the volume mixing ratio retrievals for some molecules. Analytical expressions are presented for the imaginary components of the polynomial ratios from the Humlicek algorithm, which provides approximate solutions to the complex probability function. These imaginary components are employed in the calculation of line mixing using the Rosenkranz first order approximation. Examples of line mixing in ACE measurements are presented, including a set of CH 4 lines that exhibit both line mixing and speed dependence. An efficient, analytical approach is proposed for calculating line shapes with a combination of line mixing and speed dependence. FORTRAN routines for calculating line mixing effects are provided as a supplement to the paper.

Impact of Solid Particulate on Brittle Materials

This paper deals with the erosion of brittle materials due to the impact of solid particulate. Erosion tests are characterized by a great uncertainty in the results and high costs. These difficulties led the author to develop alternative methods based on numerical simulations for calculating the erosion rate of brittle and ductile materials. Erosion criteria, based on fracture energy, were developed and implemented as Fortran routines in a commercial FEM code. In this paper the proposed method was applied to calculate the erosion rate of glass and glass ceramics. Fracture energy was obtained using four points bending tests. The numerical results were compared with a developed theoretical model and with the experimental results, available in the literature. In spite of several simplifying assumptions, there was a good agreement between the numerical and experimental results. This method could be a powerful tool for assessing the erosion resistance of components in real operating conditions.