Equilibrium Solver Research Articles

Sensitivity analysis of separable traffic equilibrium equilibria with application to bilevel optimization in network design

We provide a sensitivity analysis of separable traffic equilibrium models with travel cost and demand parameters. We establish that while equilibrium link flows may not always be directionally differentiable (even when the link travel costs are strictly increasing), travel demands and link costs are; this improves the general results of Patriksson [Patriksson, M., 2004. Sensitivity analysis of traffic equilibria. Transportation Science 37, 258–281]. The new results contradict common belief that equilibrium cost and demand sensitivities hinge on that of equilibrium flows. The paper by Tobin and Friesz [Tobin, R.L., Friesz, T.L., 1988. Sensitivity analysis for equilibrium network flow. Transportation Science 22, 242–250] brought the classic non-linear programming subject of sensitivity analysis to transportation science. Theirs is still the most widely used device by which “gradients” of traffic equilibrium solutions are calculated, for use in bilevel transportation planning applications such as network design, origin–destination (OD) matrix estimation and problems where link tolls are imposed on the users in order to reach a traffic management objective. However, it is not widely understood that the regularity conditions proposed by them are stronger than necessary. Also, users of their method sometimes misunderstand its limitations and are not aware of the computational advantages offered by more recent methods. In fact, a more often applicable formula was proposed already by Qiu and Magnanti [Qiu, Y., Magnanti, T.L., 1989. Sensitivity analysis for variational inequalities defined on polyhedral sets. Mathematics of Operations Research 14, 410–432], and Bell and Iida [Bell, M.G.H., Iida, Y., 1997. Transportation Network Analysis. John Wiley & Sons, Chichester, UK] describe one of the cases in practice in which the formula by Tobin and Friesz would not be able to generate sensitivity information, because one of their regularity conditions fails to hold. This paper provides an overview of this formula, and illustrates by means of examples that there are several cases where it is not applicable. Our findings are illustrated with small numerical examples, as are our own analysis. The findings of this paper are hoped to motivate replacing the previous approach with the more often applicable one, not only because of this fact but equally importantly because it is intuitive and also can be much more efficiently utilized: the sensitivity problem that provides the directional derivative is a linearized traffic equilibrium problem, and the sensitivity information can be generated efficiently by only slightly modifying a state-of-the-art traffic equilibrium solver. This is essential for bringing the use of sensitivity analysis in transportation planning beyond the solution of only toy problems. We finally utilize a new sensitivity solver in the preliminary testing of a simple heuristic for bilevel optimization in continuous traffic network design, and compare it favourably to previous heuristics on known small-scale problems.

A computationally efficient Multicomponent Equilibrium Solver for Aerosols (MESA)

Development and application of a new Multicomponent Equilibrium Solver for Aerosols (MESA) is described for systems containing H+, NH4+, Na+, Ca2+, SO42−, HSO4−, NO3−, and Cl− ions. The equilibrium solution is obtained by integrating a set of pseudo‐transient ordinary differential equations describing the precipitation and dissolution reactions for all the possible salts to steady state. A comprehensive temperature dependent mutual deliquescence relative humidity (MDRH) parameterization is developed for all the possible salt mixtures, thereby eliminating the need for a rigorous numerical solution when ambient RH is less than MDRH(T). The solver is unconditionally stable, mass conserving, and shows robust convergence. Performance of MESA was evaluated against the Web‐based AIM Model III, which served as a benchmark for accuracy, and the EQUISOLV II solver for speed. Important differences in the convergence and thermodynamic errors in MESA and EQUISOLV II are discussed. The average ratios of speeds of MESA over EQUISOLV II ranged between 1.4 and 5.8, with minimum and maximum ratios of 0.6 and 17, respectively. Because MESA directly diagnoses MDRH, it is significantly more efficient when RH < MDRH. MESA's superior performance is partially due to its “hard‐wired” code for the present system as opposed to EQUISOLV II, which has a more generalized structure for solving any number and type of reactions at temperatures down to 190 K. These considerations suggest that MESA is highly attractive for use in 3‐D aerosol/air‐quality models for lower tropospheric applications (T > 240 K) in which both accuracy and computational efficiency are critical.

Solving the 3D MHD equilibrium equations in toroidal geometry by Newton’s method

We describe a novel form of Newton’s method for computing 3D MHD equilibria. The method has been implemented as an extension to the hybrid spectral/finite-difference Princeton Iterative Equilibrium Solver (PIES) which normally uses Picard iteration on the full nonlinear MHD equilibrium equations. Computing the Newton functional derivative numerically is not feasible in a code of this type but we are able to do the calculation analytically in magnetic coordinates by considering the response of the plasma’s Pfirsch–Schlüter currents to small changes in the magnetic field. Results demonstrate a significant advantage over Picard iteration in many cases, including simple finite- β stellarator equilibria. The method shows promise in cases that are difficult for Picard iteration, although it is sensitive to resolution and imperfections in the magnetic coordinates, and further work is required to adapt it to the presence of magnetic islands and stochastic regions.

PIES free boundary stellarator equilibria with improved initial conditions

The MFBE procedure developed by Strumberger (1997 Nucl. Fusion 37 19) is used to provide an improved starting point for free boundary equilibrium computations in the case of W7-X (Nührenberg and Zille 1986 Phys. Lett. A 114 129) using the Princeton iterative equilibrium solver (PIES) code (Reiman and Greenside 1986 Comput. Phys. Commun. 43 157). Transferring the consistent field found by the variational moments equilibrium code (VMEC) (Hirshmann and Whitson 1983 Phys. Fluids 26 3553) to an extended coordinate system using the VMORPH code, a safe margin between plasma boundary and PIES domain is established. The new EXTENDER_P code implements a generalization of the virtual casing principle, which allows field extension both for VMEC and PIES equilibria. This facilitates analysis of the 5/5 islands of the W7-X standard case without including them in the original PIES computation.

Linear tearing modes calculation for diverted tokamak configurations

Writing the expression of the potential energy in terms of the perturbation of the flux function, and performing an Euler minimization, one obtains a system of ordinary differential equations in that perturbation. For a diverted configuration, the usual vanishing boundary conditions for the perturbed flux function at the magnetic axis and at infinity can no longer be used. An approach to fix “natural” boundary conditions for the perturbed flux function just at the plasma boundary has been developed; this replaces the vanishing boundary conditions at infinity. To obtain the necessary metric coefficients from an equilibrium solver, a function exhibiting the same singularity at the X point has been introduced, and thus the separatrix contour was described by a reduced number of moments.

Constructing integrable high-pressure full-current free-boundary stellarator magnetohydrodynamic equilibrium solutions

For the (non-axisymmetric) stellarator class of plasma confinement devices to be feasible candidates for fusion power stations it is essential that, to a good approximation, the magnetic field lines lie on nested flux surfaces; however, the inherent lack of a continuous symmetry implies that magnetic islands responsible for breaking the smooth topology of the flux surfaces are guaranteed to exist. Thus, the suppression of magnetic islands is a critical issue for stellarator design, particularly for small aspect ratio devices. Pfirsch–Schlüter currents, diamagnetic currents and resonant coil fields contribute to the formation of magnetic islands, and the challenge is to design the plasma and coils such that these effects cancel.Magnetic islands in free-boundary high-pressure full-current stellarator magnetohydrodynamic equilibria are suppressed using a procedure based on the Princeton Iterative Equilibrium Solver (Reiman and Greenside 1986 Comput. Phys. Commun. 43 157) which iterates the equilibrium equations to obtain the plasma equilibrium. At each iteration, changes to a Fourier representation of the coil geometry are made to cancel resonant fields produced by the plasma. The changes are constrained to preserve certain measures of engineering acceptability and to preserve the stability of ideal kink modes. As the iterations continue, the coil geometry and the plasma simultaneously converge to an equilibrium in which the island content is negligible, the plasma is stable to ideal kink modes, and the coils satisfy engineering constraints. The method is applied to a candidate plasma and coil design for the National Compact Stellarator eXperiment (Reiman et al 2001 Phys. Plasma 8 2083).

Local equilibrium of nonrotating plasmas

What is the minimal information that must be supplied on a constant pressure surface of a toroidal plasma to completely specify a local magnetohydrodyamic equilibrium? Clearly the shape of the surface must be specified, but other information such as the separation between neighboring surfaces must also be given. This paper explores the conditions required for three-dimensional equilibria and how local equilibrium constraints are related to global constraints. The resulting conditions provide useful checks on the accuracy of equilibrium solvers, particularly for solvers for three-dimensional equilibria, and allow local refinements in the accuracy of equilibria.

Reduction of islands in full-pressure stellarator equilibria

The control of magnetic islands is a crucial issue in designing stellarators. Islands are associated with resonant radial magnetic fields at rational rotational-transform surfaces and can lead to chaos and poor plasma confinement. In this article it is shown that variations in the resonant fields of a full pressure stellarator equilibrium can be related to variations in the boundary via a coupling matrix, and that inversion of this matrix determines a boundary modification for which the island content is significantly reduced. The numerical procedure is described and the results of island optimization are presented. Equilibria with islands are computed using the Princeton Iterative Equilibrium Solver, and resonant radial fields are calculated via construction of quadratic-flux-minimizing surfaces. A design candidate for the National Compact Stellarator Experiment [A. Reimann, L. Ku, D. Monticello et al., Phys. Plasmas 8, 2083 (2001)], which has a large island, is used to illustrate the technique. Small variations in the boundary shape are used to reduce island size and to reverse the phase of a major island chain.

Magnetic flux evolution in highly shaped plasmas

The resistive evolution of magnetic flux in toroidal devices is studied. The formulation is applicable to general nonaxisymmetric (three-dimensional) toroidal configurations. In particular, it can treat highly shaped, high β, three-dimensional stellarator configurations, as well as two-dimensional (axisymmetric) tokamak plasmas. The time evolution of the poloidal magnetic flux is posed in terms of the rotational transform, ι̸, and allows for a transparent inclusion of stellarator specific current-free contributions to ι̸. Strong diamagnetic and paramagnetic contributions to toroidal magnetic flux, as evident in spherical tokamaks and similar concepts, are calculated by direct iteration with an equilibrium solver. The nonlinear evolution equation is derived using a susceptance matrix formulation originally introduced by Grad and co-workers [Bateman, Nucl. Fusion 13, 227 (1973)]. Here, it is extended to general, nonstraight field line coordinate systems. The basic equations are described, explicit expressions for the susceptance matrix are given, and example applications using the stand-alone code, THRIFT (THRee-dimensional Inductive Flux evolution in Toroidal devices), are discussed.

A parameterization of an aerosol physical chemistry model for the NH3/H2SO4/HNO3/H2O system at cold temperatures

Simple expressions are fitted to the results obtained from ion interaction thermodynamic models for calculating HNO3 and H2O vapor pressures over the NH3/H2SO4/HNO3/H2O system at cold temperatures. The vapor pressure expressions are incorporated into a mass conserving equilibrium solver for computing aerosol compositions in the lower stratosphere and upper troposphere. The compositions calculated from the aerosol physical chemistry model (APCM) are compared against previous parameterizations. The APCM compositions are in better agreement with the compositions obtained from ion interaction models than from other previous formulations of the NH3/H2SO4/HNO3/H2O system. The only advantage of the APCM over the ion interaction approach is that the numerical scheme used in the model is fast and efficient for incorporation into large‐scale models. The APCM is used to calculate HNO3 solubility in aminoniated aerosols as a function of HNO3, H2SO4 and NH3 mass loadings in the lower stratosphere and upper troposphere. While the uptake of HNO3 by ammoniated aerosols is strongly dependent upon the solution neutrality (or pH), we find that in both the lower stratosphere and upper troposphere a significant fraction of HNO3 will exist in aerosol solutions near and below the ice frost point irrespective of solution neutrality.

Studying the effects of calcium and magnesium on size-distributed nitrate and ammonium with EQUISOLV II

A chemical equilibrium code was improved and used to show that calcium and magnesium have a large yet different effect on the aerosol size distribution in different regions of Los Angeles. In the code, a new technique of solving individual equilibrium equations was developed. The technique, the analytical equilibrium iteration (AEI) method, gives the same solutions (to at least 7 decimal places) as the previous technique used, the mass-flux iteration (MFI) method, but consumes 13–48 times less computer time. The model was also updated to include treatment of potassium, calcium, magnesium, and carbonate. Previously, it treated only nitrate, ammonium, chloride, sulfate, and sodium. Predictions from the updated code, EQUISOLV II, were compared with data from an eight-stage Berner impactor at Long Beach, Claremont, and Riverside during the Southern California Air Quality Study. When any equilibrium solver is applied between the gas phase and multiple aerosol size bins, unique solutions are possible only when solids (e.g., NH 4NO 3) that form from two gas-phase species are absent. For this study, unique solutions were possible only when the relative humidity exceeded 62%, and only cases in this regime are discussed. Base-case predictions of nitrate and ammonium matched observations well in most size bins of every case. When Ca and Mg were removed from calculations, coarse-mode nitrate decreased at Long Beach, as expected, to maintain charge balance. At Riverside, removing Ca and Mg had the opposite effect, increasing coarse-mode nitrate, shifting it from the accumulation mode. The reason is explained in terms of mean mixed activity coefficients. At Claremont, the charge-balance and activity-coefficient effects nearly canceled each other.

Finite-beta equilibria for Wendelstein 7-X configurations using the Princeton Iterative Equilibrium Solver code

Fixed- and free-boundary equilibria for Wendelstein 7-X (W7-X) [W. Lotz et al., Plasma Physics and Controlled Nuclear Fusion Research 1990 (Proc. 13th Int. Conf. Washington, DC, 1990), (International Atomic Energy Agency, Vienna, 1991), Vol. 2, p. 603] configurations are calculated using the Princeton Iterative Equilibrium Solver (PIES) [A. H. Reiman et al., Comput. Phys. Commun., 43, 157 (1986)] to deal with magnetic islands and stochastic regions. Usually, these W7-X configurations require a large number of iterations for PIES convergence. Here, two methods have been successfully tested in an attempt to decrease the number of iterations needed for convergence. First, periodic sequences of different blending parameters are used. Second, the initial guess is vastly improved by using results of the Variational Moments Equilibrium Code (VMEC) [S. P. Hirshmann et al., Phys. Fluids 26, 3553 (1983)]. Use of these two methods have allowed verification of the Hamada condition and tendency of “self-healing” of islands has been observed.

Mode coupling in two-dimensional magnetohydrodynamic flows

The spectrum of incompressible waves and instabilities of two-dimensional plasma geometries with background flow is calculated. The equilibrium is solved numerically by the recently developed program FLow Equilibrium Solver (FLES). The spectra of the equilibria are computed by means of another new program, the INcompressible 2-dimensional FLow Eigenvalue Solver (IN2FLES). Magnetic instabilities and instabilities driven by the two-dimensionality and the flow are found. For linear equilibria, the eigenvalues for elliptical geometries remain close to the curves on which the eigenvalues for circular geometries lie. These curves may be found for unbounded domains by a calculation in Fourier space [see Lifschitz, A. In: Proceedings of 1995 International Workshop on Operator Theory and Applications (ed. R. Mennicken and C. Tretter), pp. 97–117, Birkhäuser, Boston, 1998]. Here the relation between a new continuous spectrum of unbounded domains and the discrete spectrum of bounded domains is investigated. Finally, it is found that the two-dimensionality and the background flow may lead to an overstable cluster-point.

Two-dimensional equilibrium in coronal magnetostatic flux tubes: an accurate equilibrium solver

To study linearized magnetohydrodynamic (MHD) waves, continuous spectra, and instabilities in coronal magnetic flux tubes that are anchored in dense chromospheric and photospheric regions, a two-dimensional numerical code, called PARIS, has been developed. PARIS solves the pertinent nonlinear Grad-Shafranov type, partial differential equation for the magnetic flux on a flux coordinate grid. Both a straight field line coordinate system and on orthogonal flux coordinate system are exploited. Isoparametric bicubic Hermite finite elements have been adopted to solve the Grad-Shafranov-like equation. These elements allow for a continuous representation of the flux and the gradient of the flux throughout the tube and can be aligned conveniently along the boundary of the tube. These properties are important to obtain an accurate representation of the solution on flux coordinate grids. An analytical test case is used to show that accurate solutions have been obtained, even for a small number of grid points. The equilibria calculated by PARIS are used to study the continuous spectra of two-dimensional magnetic flux tubes. One illustrative example is given here; extensive results are presented elsewhere (A.J.C. Beliën, S. Poedts and J.P. Goedbloed, Astron. Astrophys. 322 (1997) 995). The equilibria obtained by PARIS are also well suited to calculate the stability and the normal mode MHD spectrum.

Simulating equilibrium within aerosols and nonequilibrium between gases and aerosols

A numerical method to solve chemical equilibrium equations and a method of coupling the equilibrium calculations to nonequilibrium growth and evaporation are discussed. The equilibrium program solves any number of equations for gas, aqueous, ionic, and solid equilibrium concentrations over large spatial grids and particle size grids. It also simultaneously computes electrolyte mean mixed activity coefficients and aerosol liquid water content. Mean mixed activity coefficient calculations require mean binary activity coefficient information. Temperature‐dependent mean binary activity coefficient polynomials were constructed using mean binary activity coefficient data at 298 K, apparent molal enthalpy data, and apparent molal heat capacity data. The equilibrium solver is mole and charge conserving, requires iteration, but always converges. Solutions to the equilibrium equations are used for two purposes. The first is to estimate surface vapor pressures over particles containing a solution and/or a solid phase. Such vapor pressures are then applied in gas‐aerosol transfer equations. The second is to estimate intraparticle composition and size immediately after gas‐aerosol transfer.

Combustion By-Products Destruction by Arc Plasma

The parameters affecting complete thermal decomposition of methane, and 1,1–dichloroethane in a dc argon plasma were investigated. The 1,1–dichloroethane was chosen as a representative reactant from chlorinated–compounds and chemical warfare agents. Theoretical studies were conducted using STANJAN thermochemical equilibrium solver–code on stoichiometric mixtures of methane–air and dichloroethane–air systems. In addition, the effect of 10% excess oxygen in the stoichiometric methane–air system was studied in order to understand the role of oxygen at a plasma condition. The results from these studies give valuable information in predicting the species concentrations and evaluating the quenching process. To conduct experimental investigations, a 10 kW laboratory scale dc plasma generator was built and operated to produce a stable argon plasma. The representative organic species: methane, and dichloroethane were included in the argon stream for decomposition in the plasma. The decomposition of the species in the plasma was monitored by optical emission spectroscopy. The extent of thermal decomposition was measured by the emission spectrum obtained from the data acquisition system. The gas–flow rate of the reactants and the power dissipated in the plasma were found to be important parameters for achieving complete decomposition of a particular chemical species. The argon–methane plasma was used as a reference, because methane is the simplest saturated hydrocarbon and we could use it to simulate the effluent from an incinerator. The plasma destruction of hazardous wastes was thus found to be very promising.

Neural network differential equation and plasma equilibrium solver.

A new generally applicable method to solve differential equations, based on neural networks, is proposed. Straightforward to implement, finite differences and coordinate transformations are not used. The neural network provides a flexible and compact base for representing the solution, found through the global minimization of an error functional. As a proof of principle, a two-dimensional ideal magnetohydrodynamic plasma equilibrium is solved. Since no particular topology is assumed, the technique is especially promising for the three-dimensional plasma equilibrium problem.

Oxygen Based Claus Process for Recovery of Sulfur from H2S Gases

Novel oxygen‐based modifications of the Claus process for sulfur recovery from gases containing hydrogen sulfide are analyzed and compared with the conventional air‐based technology. Oxygen‐based modifications include systems with air enrichment by oxygen, with oxygen‐air‐water oxidizers, with recycle of cooled products from boilers back into the Claus furnace, and systems with consecutive or parallel‐consecutive Claus furnaces. The analysis used in this work involves the process thermodynamics and kinetics, as well as modeling of the entire system. The thermodynamics was evaluated using the STANJAN chemical equilibrium solver program, which utilizes the principle of minimization of free Gibbs energy. The process kinetics was evaluated using published data on full‐scale Claus systems. The performance of various process modifications was evaluated by the use of numerical analyses. The major findings of the work are that oxygen‐based Claus systems can be used for: (1) Upgrading existing systems and for use in new systems; and (2) reducing the equipment size requirement in new systems. In addition the Claus process may be used for a variety of conditions met in various industries.

Computation of MHD equilibrium of tokamak plasma

Abstract Computation of the MHD equilibrium of a tokamak plasma is reviewed as comprehensively as possible. The basic equation of this problem is the Grad-Shafranov equation. General remarks on this equation and related issues are, first, summarized with historical survey of the MHD equilibrium solution, where some mathematical discussions on the numerical analysis of the problem are also presented. Distinguishing features of this problem are seen in treatment of the boundary condition and constraining conditions and we describe them in a rather detailed manner. In the main part of this review paper we present a concrete description on the numerical procedures of the MHD equilibrium solvers which are classified into two groups, that is, the real space solvers and the inverse equilibrium solvers. We also describe topics on more general equilibrium models, that is, the equilibrium with steady flow, anisotropic equilibria, equilibria with specified current sources, and equilibrium evolution. Brief comments on three-dimensional equilibrium solvers are also presented. As for application of the MHD equilibrium solvers we present only a small part, that is, beta limit optimization, design of external coils, analysis of positional instability, and analysis of experimentally obtained data from electromagnetic measurement. It is concluded that among various kinds of numerical solution methods we can usually find most adequate ones for the present problem.

Modelling the behaviour of the magnetic field diagnostic coils on the W VII-AS stellarator using a three-dimensional equilibrium code

Three-dimensional equilibrium solvers have been extensively used to investigate stellarator equilibria, but it was only recently that they have been adapted to model the magnetic diagnostics on a functioning experiment. The paper presents the new computer code NEMEC-DIAGNO which was written to model the diagnostic coils on the W VII-AS stellarator. This code extends the free-boundary equilibrium solver NEMEC to calculate the magnetic field and the vector potential due to the plasma currents alone and incorporates detailed models of the coils themselves. In its application to W VII-AS over a range of its operation, the code yields a set of simple curves which provide a prescription for finding the average plasma beta under varying conditions of pressure profile and net longitudinal ‘bootstrap’ current.