Quadratic Flux Minimizing Surfaces
Quadratic-Flux-Minimizing (QFM) surfaces are a useful tool in vacuum stellarator optimization, where, given a coil magnetic field, one may extract a surface in space which minimizes the quadratic magnetic flux passing through it, that is,
Where \(B_n = \mathbf{B} \cdot \mathbf{n}\) is the normal magnetic field on the surface \(S\), and the optimization is understood to be done under the constraint of the surface enclosing some amount of net toroidal magnetic flux. The surface \(S\) itself is described typically with a double Fourier series in the poloidal and toroidal angles (See the docs on basis functions for more details) describing the cylindrical \(R\) and \(Z\) coordinates of the surface, so the optimization is done by varying the Fourier coefficients describing the surface (or potentially, the magnetic field’s degrees of freedom as well). This short notebook will give an example of using QFM surfaces to extrace a surface from a magnetic field, and then solving the interior equilibrum with DESC and comparing it to the original magnetic field.
[1]:
import sys
import os
sys.path.insert(0, os.path.abspath("."))
sys.path.append(os.path.abspath("../../../"))
If you have access to a GPU, uncomment the following two lines before any DESC or JAX related imports. You should see about an order of magnitude speed improvement with only these two lines of code!
[ ]:
# from desc import set_device
# set_device("gpu")
As mentioned in DESC Documentation on performance tips, one can use compilation cache directory to reduce the compilation overhead time. Note: One needs to create jax-caches folder manually.
[3]:
# import jax
# jax.config.update("jax_compilation_cache_dir", "../jax-caches")
# jax.config.update("jax_persistent_cache_min_entry_size_bytes", -1)
# jax.config.update("jax_persistent_cache_min_compile_time_secs", 0)
[4]:
from desc.objectives import (
SurfaceQuadraticFlux,
ToroidalFlux,
ObjectiveFunction,
)
from desc.magnetic_fields import SplineMagneticField
from desc.geometry import FourierRZToroidalSurface
from desc.optimize import Optimizer
from desc.grid import LinearGrid
from desc.plotting import plot_2d, poincare_plot, plot_boundary, plot_surfaces
import numpy as np
First, we will load a magnetic field from an mgrid file (the same as is used in the stellarator free boundary notebook)
[5]:
extcur = [4700.0, 1000.0]
field = SplineMagneticField.from_mgrid(
"../../../tests/inputs/mgrid_test.nc", extcur=extcur
)
Our initial surface is a circular surface, which we see is not well-approximating any flux surfaces in this vacuum field.
[6]:
surface0 = FourierRZToroidalSurface(
R_lmn=[0.7, 0.1],
modes_R=[[0, 0], [1, 0]],
Z_lmn=[-0.1],
modes_Z=[[-1, 0]],
sym=True,
NFP=field.NFP,
M=6,
N=6,
)
fig, ax = plot_surfaces(surface0, rho=1.0, theta=0)
poincare_plot(
field,
R0=np.linspace(0.7, 0.8, 10),
Z0=np.zeros(10),
ax=ax,
NFP=field.NFP,
ntransit=50,
);
Let’s setup the optimization problem to find the QFM surface.
[7]:
# SurfaceQuadraticFlux objective is the one to use when trying to minimize quadratic flux for a surface
qflux = SurfaceQuadraticFlux(
surface0,
field,
eval_grid=LinearGrid(M=2 * surface0.M, N=2 * surface0.M, NFP=surface0.NFP),
field_fixed=True, # field is not being optimized
bs_chunk_size=10,
)
# Must include a ToroidalFlux or a Volume target to ensure we don't get a trivial solution of the surface collapsing to a point
target_psi = -0.035
tflux = ToroidalFlux(
surface0,
field,
target=target_psi, # the same toroidal flux as used in the free bdry eq
eq_fixed=False, # we will allow the thing passed in (QFM surface in this case) to vary
field_fixed=True, # field is not being optimized, so we fix it
)
obj = ObjectiveFunction(qflux)
opt = Optimizer("lsq-auglag")
(qfm_surf,), _ = opt.optimize(
surface0,
objective=obj,
constraints=(tflux,),
ftol=1e-5,
maxiter=500,
verbose=3,
copy=True,
)
Building objective: Surface Quadratic Flux
Precomputing transforms
Timer: Precomputing transforms = 934 ms
Timer: Objective build = 978 ms
Building objective: toroidal-flux
Precomputing transforms
Timer: Precomputing transforms = 530 ms
Timer: Objective build = 1.22 sec
Number of parameters: 169
Number of objectives: 625
Number of equality constraints: 1
Number of inequality constraints: 0
Timer: Initializing the optimization = 2.37 sec
Starting optimization
Using method: lsq-auglag
Solver options:
------------------------------------------------------------
Maximum Function Evaluations : 2501
Maximum Allowed Total Δx Norm : inf
Scaled Termination : True
Trust Region Method : qr
Initial Trust Radius : 2.519e-01
Maximum Trust Radius : inf
Minimum Trust Radius : 2.220e-16
Trust Radius Increase Ratio : 4.000e+00
Trust Radius Decrease Ratio : 2.500e-01
Trust Radius Increase Threshold : 7.500e-01
Trust Radius Decrease Threshold : 5.000e-01
Alpha Omega : 1.000e+00
Beta Omega : 1.000e+00
Alpha Eta : 1.000e-01
Beta Eta : 9.000e-01
Omega : 1.000e-02
Eta : 1.000e-02
Tau : 9.000e-01
------------------------------------------------------------
Iteration Total nfev Cost Cost reduction Step norm Optimality Constr viol. Penalty param max(|mltplr|)
0 1 6.632e+01 9.991e+00 1.383e-02 1.000e+01 0.000e+00
1 2 6.099e+01 5.329e+00 2.309e-03 9.376e+00 1.441e-02 1.000e+01 0.000e+00
2 3 4.257e+01 1.843e+01 9.279e-03 7.151e+00 1.665e-02 1.000e+01 0.000e+00
3 4 5.031e+00 3.753e+01 3.739e-02 1.583e+00 2.437e-02 1.000e+01 0.000e+00
4 5 5.034e+00 -2.830e-03 1.384e-01 1.873e+00 2.534e-03 1.000e+01 0.000e+00
5 6 1.305e-01 4.904e+00 2.762e-02 2.340e-01 1.459e-04 1.000e+01 0.000e+00
6 8 3.707e-03 1.268e-01 8.902e-03 3.072e-02 9.054e-04 1.000e+01 0.000e+00
7 10 6.144e-04 3.092e-03 3.706e-03 1.288e-02 1.066e-03 1.000e+01 0.000e+00
8 12 2.904e-04 3.240e-04 3.565e-03 9.099e-03 7.678e-04 1.000e+01 0.000e+00
9 13 1.215e-04 1.689e-04 3.095e-03 4.758e-03 3.251e-04 1.000e+01 0.000e+00
10 15 2.765e-05 9.384e-05 6.619e-04 3.506e-04 2.391e-04 1.000e+01 2.391e-03
11 17 2.617e-05 1.487e-06 8.471e-04 3.374e-04 8.699e-06 1.000e+01 2.391e-03
12 19 2.494e-05 1.225e-06 5.958e-04 3.091e-04 5.704e-05 1.000e+01 2.391e-03
13 21 2.389e-05 1.046e-06 6.454e-04 2.430e-04 7.572e-05 1.000e+01 2.391e-03
14 23 2.307e-05 8.274e-07 6.645e-04 3.017e-04 8.494e-05 1.000e+01 2.391e-03
15 25 2.229e-05 7.765e-07 6.630e-04 3.207e-04 9.258e-05 1.000e+01 2.391e-03
16 27 2.165e-05 6.393e-07 7.073e-04 3.783e-04 9.970e-05 1.000e+01 2.391e-03
17 28 2.114e-05 5.065e-07 7.595e-04 4.403e-04 1.063e-04 1.000e+01 2.391e-03
18 29 2.065e-05 4.970e-07 7.594e-04 4.348e-04 1.121e-04 1.000e+01 2.391e-03
19 30 2.016e-05 4.917e-07 7.323e-04 3.903e-04 1.171e-04 1.000e+01 2.391e-03
20 31 1.966e-05 4.920e-07 6.975e-04 3.288e-04 1.213e-04 1.000e+01 2.391e-03
21 32 1.917e-05 4.923e-07 6.650e-04 2.634e-04 1.249e-04 1.000e+01 2.391e-03
22 33 1.871e-05 4.593e-07 6.489e-04 2.590e-04 1.283e-04 1.000e+01 2.391e-03
23 34 1.832e-05 3.949e-07 6.559e-04 2.786e-04 1.317e-04 1.000e+01 2.391e-03
24 35 1.798e-05 3.370e-07 6.694e-04 2.926e-04 1.352e-04 1.000e+01 2.391e-03
25 36 1.768e-05 2.969e-07 6.753e-04 2.903e-04 1.387e-04 1.000e+01 2.391e-03
26 37 1.742e-05 2.630e-07 6.762e-04 2.783e-04 1.421e-04 1.000e+01 2.391e-03
27 38 1.718e-05 2.355e-07 6.771e-04 2.643e-04 1.454e-04 1.000e+01 2.391e-03
28 39 1.697e-05 2.189e-07 6.793e-04 2.504e-04 1.484e-04 1.000e+01 2.391e-03
29 40 1.676e-05 2.100e-07 6.822e-04 2.419e-04 1.512e-04 1.000e+01 2.391e-03
30 41 1.648e-05 2.733e-07 1.695e-04 1.496e-05 1.523e-04 1.000e+01 8.679e-04
31 42 1.640e-05 8.553e-08 7.492e-04 2.359e-04 8.296e-06 1.000e+01 8.679e-04
32 43 1.611e-05 2.846e-07 1.698e-04 1.431e-05 5.094e-06 1.000e+01 8.679e-04
33 45 1.606e-05 4.795e-08 1.708e-04 1.412e-05 5.143e-06 1.000e+01 8.679e-04
34 47 1.602e-05 4.680e-08 1.713e-04 1.381e-05 5.616e-06 1.000e+01 8.679e-04
35 49 1.597e-05 4.607e-08 1.717e-04 1.360e-05 6.150e-06 1.000e+01 8.679e-04
36 51 1.593e-05 4.542e-08 1.721e-04 1.341e-05 6.690e-06 1.000e+01 8.679e-04
37 53 1.588e-05 4.478e-08 1.725e-04 1.325e-05 7.238e-06 1.000e+01 8.679e-04
38 55 1.584e-05 4.414e-08 1.729e-04 1.310e-05 7.782e-06 1.000e+01 8.679e-04
39 57 1.579e-05 4.350e-08 1.733e-04 1.294e-05 8.301e-06 1.000e+01 8.679e-04
40 59 1.575e-05 4.285e-08 1.737e-04 1.280e-05 8.813e-06 1.000e+01 8.679e-04
41 61 1.571e-05 4.219e-08 1.742e-04 1.266e-05 9.320e-06 1.000e+01 8.679e-04
42 63 1.567e-05 4.151e-08 1.747e-04 1.252e-05 9.823e-06 1.000e+01 8.679e-04
43 65 1.563e-05 4.083e-08 1.751e-04 1.238e-05 1.032e-05 1.000e+01 8.679e-04
44 67 1.559e-05 4.016e-08 1.757e-04 1.224e-05 1.081e-05 1.000e+01 8.679e-04
45 69 1.555e-05 3.948e-08 1.762e-04 1.213e-05 1.130e-05 1.000e+01 8.679e-04
46 71 1.551e-05 3.879e-08 1.767e-04 1.203e-05 1.179e-05 1.000e+01 8.679e-04
47 73 1.547e-05 3.810e-08 1.772e-04 1.193e-05 1.227e-05 1.000e+01 8.679e-04
48 75 1.543e-05 3.741e-08 1.777e-04 1.183e-05 1.274e-05 1.000e+01 8.679e-04
49 77 1.540e-05 3.672e-08 1.782e-04 1.173e-05 1.320e-05 1.000e+01 8.679e-04
50 79 1.536e-05 3.606e-08 1.787e-04 1.163e-05 1.367e-05 1.000e+01 8.679e-04
51 81 1.532e-05 3.541e-08 1.793e-04 1.154e-05 1.413e-05 1.000e+01 8.679e-04
52 83 1.529e-05 3.478e-08 1.798e-04 1.145e-05 1.459e-05 1.000e+01 8.679e-04
53 85 1.526e-05 3.416e-08 1.805e-04 1.135e-05 1.504e-05 1.000e+01 8.679e-04
54 87 1.522e-05 3.357e-08 1.811e-04 1.122e-05 1.546e-05 1.000e+01 8.679e-04
55 89 1.519e-05 3.300e-08 1.817e-04 1.106e-05 1.587e-05 1.000e+01 8.679e-04
56 91 1.516e-05 3.246e-08 1.823e-04 1.087e-05 1.625e-05 1.000e+01 8.679e-04
57 93 1.512e-05 3.196e-08 1.830e-04 1.066e-05 1.661e-05 1.000e+01 8.679e-04
58 95 1.509e-05 3.148e-08 1.837e-04 1.043e-05 1.696e-05 1.000e+01 8.679e-04
59 97 1.506e-05 3.104e-08 1.845e-04 1.042e-05 1.729e-05 1.000e+01 8.679e-04
60 99 1.503e-05 3.064e-08 1.853e-04 1.082e-05 1.763e-05 1.000e+01 8.679e-04
61 101 1.500e-05 3.028e-08 1.861e-04 1.119e-05 1.795e-05 1.000e+01 8.679e-04
62 103 1.497e-05 2.996e-08 1.870e-04 1.156e-05 1.829e-05 1.000e+01 8.679e-04
63 105 1.494e-05 2.967e-08 1.880e-04 1.200e-05 1.866e-05 1.000e+01 8.679e-04
64 107 1.491e-05 2.942e-08 1.890e-04 1.240e-05 1.903e-05 1.000e+01 8.679e-04
65 109 1.488e-05 2.914e-08 1.897e-04 1.268e-05 1.942e-05 1.000e+01 8.679e-04
66 111 1.485e-05 2.883e-08 1.900e-04 1.280e-05 1.983e-05 1.000e+01 8.679e-04
67 113 1.483e-05 2.849e-08 1.900e-04 1.283e-05 2.024e-05 1.000e+01 8.679e-04
68 115 1.480e-05 2.810e-08 1.895e-04 1.293e-05 2.067e-05 1.000e+01 8.679e-04
69 117 1.477e-05 2.767e-08 1.887e-04 1.378e-05 2.109e-05 1.000e+01 8.679e-04
70 119 1.474e-05 2.719e-08 1.876e-04 1.455e-05 2.154e-05 1.000e+01 8.679e-04
71 121 1.472e-05 2.666e-08 1.862e-04 1.519e-05 2.199e-05 1.000e+01 8.679e-04
72 123 1.469e-05 2.607e-08 1.843e-04 1.569e-05 2.245e-05 1.000e+01 8.679e-04
73 125 1.466e-05 2.545e-08 1.821e-04 1.599e-05 2.291e-05 1.000e+01 8.679e-04
74 127 1.464e-05 2.482e-08 1.798e-04 1.607e-05 2.338e-05 1.000e+01 8.679e-04
75 129 1.462e-05 2.418e-08 1.776e-04 1.596e-05 2.384e-05 1.000e+01 8.679e-04
76 131 1.459e-05 2.354e-08 1.753e-04 1.559e-05 2.429e-05 1.000e+01 8.679e-04
77 133 1.457e-05 2.294e-08 1.729e-04 1.495e-05 2.473e-05 1.000e+01 8.679e-04
78 135 1.455e-05 2.237e-08 1.710e-04 1.402e-05 2.517e-05 1.000e+01 8.679e-04
79 137 1.452e-05 2.185e-08 1.693e-04 1.292e-05 2.559e-05 1.000e+01 8.679e-04
80 139 1.450e-05 2.140e-08 1.677e-04 1.243e-05 2.598e-05 1.000e+01 8.679e-04
81 141 1.448e-05 2.100e-08 1.663e-04 1.269e-05 2.635e-05 1.000e+01 8.679e-04
82 143 1.446e-05 2.066e-08 1.652e-04 1.285e-05 2.671e-05 1.000e+01 8.679e-04
83 145 1.444e-05 2.038e-08 1.643e-04 1.288e-05 2.705e-05 1.000e+01 8.679e-04
84 147 1.442e-05 2.015e-08 1.637e-04 1.273e-05 2.736e-05 1.000e+01 8.679e-04
85 149 1.440e-05 1.994e-08 1.633e-04 1.235e-05 2.766e-05 1.000e+01 8.679e-04
86 151 1.438e-05 1.976e-08 1.630e-04 1.174e-05 2.793e-05 1.000e+01 8.679e-04
87 153 1.436e-05 1.957e-08 1.631e-04 1.098e-05 2.819e-05 1.000e+01 8.679e-04
88 155 1.434e-05 1.940e-08 1.635e-04 1.037e-05 2.844e-05 1.000e+01 8.679e-04
89 157 1.432e-05 1.923e-08 1.643e-04 1.080e-05 2.868e-05 1.000e+01 8.679e-04
90 159 1.430e-05 1.906e-08 1.652e-04 1.172e-05 2.891e-05 1.000e+01 8.679e-04
91 161 1.429e-05 1.890e-08 1.660e-04 1.264e-05 2.913e-05 1.000e+01 8.679e-04
92 163 1.427e-05 1.876e-08 1.667e-04 1.348e-05 2.934e-05 1.000e+01 8.679e-04
93 165 1.425e-05 1.863e-08 1.670e-04 1.421e-05 2.954e-05 1.000e+01 8.679e-04
94 167 1.423e-05 1.854e-08 1.671e-04 1.478e-05 2.973e-05 1.000e+01 8.679e-04
95 169 1.421e-05 1.849e-08 1.666e-04 1.512e-05 2.991e-05 1.000e+01 8.679e-04
96 171 1.419e-05 1.849e-08 1.658e-04 1.523e-05 3.006e-05 1.000e+01 8.679e-04
97 173 1.417e-05 1.854e-08 1.645e-04 1.514e-05 3.020e-05 1.000e+01 8.679e-04
98 175 1.416e-05 1.862e-08 1.631e-04 1.494e-05 3.037e-05 1.000e+01 8.679e-04
99 177 1.414e-05 1.875e-08 1.614e-04 1.455e-05 3.053e-05 1.000e+01 8.679e-04
100 179 1.412e-05 1.890e-08 1.595e-04 1.405e-05 3.067e-05 1.000e+01 8.679e-04
101 181 1.410e-05 1.906e-08 1.576e-04 1.342e-05 3.081e-05 1.000e+01 8.679e-04
102 183 1.408e-05 1.922e-08 1.558e-04 1.277e-05 3.093e-05 1.000e+01 8.679e-04
103 185 1.406e-05 1.935e-08 1.541e-04 1.211e-05 3.105e-05 1.000e+01 8.679e-04
104 187 1.404e-05 1.944e-08 1.526e-04 1.147e-05 3.114e-05 1.000e+01 8.679e-04
105 189 1.402e-05 1.947e-08 1.513e-04 1.090e-05 3.124e-05 1.000e+01 8.679e-04
106 191 1.400e-05 1.944e-08 1.503e-04 1.040e-05 3.134e-05 1.000e+01 8.679e-04
107 193 1.398e-05 1.934e-08 1.497e-04 9.995e-06 3.145e-05 1.000e+01 5.534e-04
108 195 1.395e-05 3.518e-08 1.653e-04 1.080e-05 1.911e-06 1.000e+01 5.534e-04
109 197 1.393e-05 1.974e-08 1.489e-04 9.317e-06 5.394e-07 1.000e+01 5.534e-04
110 199 1.391e-05 1.853e-08 1.488e-04 9.286e-06 5.649e-07 1.000e+01 5.534e-04
111 201 1.389e-05 1.806e-08 1.487e-04 9.594e-06 6.653e-07 1.000e+01 5.534e-04
112 203 1.387e-05 1.758e-08 1.489e-04 1.028e-05 7.743e-07 1.000e+01 5.534e-04
113 205 1.386e-05 1.705e-08 1.491e-04 1.097e-05 8.872e-07 1.000e+01 5.534e-04
114 207 1.384e-05 1.650e-08 1.494e-04 1.166e-05 1.003e-06 1.000e+01 5.534e-04
115 209 1.382e-05 1.592e-08 1.498e-04 1.235e-05 1.118e-06 1.000e+01 5.534e-04
116 211 1.381e-05 1.533e-08 1.506e-04 1.301e-05 1.237e-06 1.000e+01 5.534e-04
117 213 1.379e-05 1.473e-08 1.519e-04 1.359e-05 1.357e-06 1.000e+01 5.534e-04
118 215 1.378e-05 1.412e-08 1.537e-04 1.408e-05 1.479e-06 1.000e+01 5.534e-04
119 217 1.377e-05 1.353e-08 1.560e-04 1.451e-05 1.599e-06 1.000e+01 5.534e-04
120 219 1.375e-05 1.295e-08 1.588e-04 1.481e-05 1.712e-06 1.000e+01 5.534e-04
121 221 1.374e-05 1.240e-08 1.623e-04 1.497e-05 1.819e-06 1.000e+01 5.534e-04
122 223 1.373e-05 1.187e-08 1.660e-04 1.502e-05 1.919e-06 1.000e+01 5.534e-04
123 225 1.372e-05 1.137e-08 1.697e-04 1.495e-05 2.020e-06 1.000e+01 5.534e-04
124 227 1.371e-05 1.092e-08 1.737e-04 1.473e-05 2.114e-06 1.000e+01 5.534e-04
125 229 1.370e-05 1.050e-08 1.778e-04 1.647e-05 2.190e-06 1.000e+01 5.534e-04
126 231 1.369e-05 1.013e-08 1.817e-04 1.815e-05 2.261e-06 1.000e+01 5.534e-04
127 233 1.368e-05 9.794e-09 1.854e-04 1.978e-05 2.331e-06 1.000e+01 5.534e-04
128 235 1.367e-05 9.498e-09 1.889e-04 2.124e-05 2.395e-06 1.000e+01 5.534e-04
129 237 1.366e-05 9.246e-09 1.921e-04 2.262e-05 2.458e-06 1.000e+01 5.534e-04
130 239 1.365e-05 9.040e-09 1.948e-04 2.393e-05 2.521e-06 1.000e+01 5.534e-04
131 241 1.364e-05 8.862e-09 1.972e-04 2.511e-05 2.580e-06 1.000e+01 5.534e-04
132 243 1.363e-05 8.710e-09 1.992e-04 2.616e-05 2.641e-06 1.000e+01 5.534e-04
133 245 1.362e-05 8.569e-09 2.009e-04 2.706e-05 2.701e-06 1.000e+01 5.534e-04
134 247 1.361e-05 8.431e-09 2.022e-04 2.767e-05 2.749e-06 1.000e+01 5.534e-04
135 249 1.361e-05 8.320e-09 2.031e-04 2.814e-05 2.795e-06 1.000e+01 5.534e-04
136 251 1.360e-05 8.211e-09 2.039e-04 2.852e-05 2.843e-06 1.000e+01 5.534e-04
137 253 1.359e-05 8.114e-09 2.044e-04 2.874e-05 2.904e-06 1.000e+01 5.534e-04
138 255 1.358e-05 8.044e-09 2.047e-04 2.886e-05 2.977e-06 1.000e+01 5.534e-04
139 257 1.357e-05 7.964e-09 2.049e-04 2.892e-05 3.048e-06 1.000e+01 5.534e-04
140 259 1.357e-05 7.890e-09 2.050e-04 2.893e-05 3.117e-06 1.000e+01 5.534e-04
141 261 1.356e-05 7.814e-09 2.049e-04 2.883e-05 3.200e-06 1.000e+01 5.534e-04
142 263 1.355e-05 7.758e-09 2.049e-04 2.870e-05 3.288e-06 1.000e+01 5.534e-04
143 265 1.354e-05 7.697e-09 2.049e-04 2.858e-05 3.381e-06 1.000e+01 5.534e-04
144 267 1.353e-05 7.634e-09 2.047e-04 2.837e-05 3.479e-06 1.000e+01 5.534e-04
145 269 1.353e-05 7.592e-09 2.045e-04 2.805e-05 3.580e-06 1.000e+01 5.534e-04
146 271 1.352e-05 7.528e-09 2.042e-04 2.768e-05 3.700e-06 1.000e+01 5.534e-04
147 273 1.351e-05 7.455e-09 2.041e-04 2.743e-05 3.823e-06 1.000e+01 5.534e-04
148 275 1.350e-05 7.395e-09 2.038e-04 2.711e-05 3.954e-06 1.000e+01 5.534e-04
149 277 1.350e-05 7.332e-09 2.036e-04 2.675e-05 4.087e-06 1.000e+01 5.534e-04
150 279 1.349e-05 7.275e-09 2.034e-04 2.635e-05 4.216e-06 1.000e+01 5.534e-04
151 281 1.348e-05 7.211e-09 2.032e-04 2.591e-05 4.342e-06 1.000e+01 5.534e-04
152 283 1.348e-05 7.143e-09 2.032e-04 2.546e-05 4.466e-06 1.000e+01 5.534e-04
153 285 1.347e-05 7.073e-09 2.032e-04 2.502e-05 4.591e-06 1.000e+01 5.534e-04
154 287 1.346e-05 6.999e-09 2.032e-04 2.457e-05 4.719e-06 1.000e+01 5.534e-04
155 289 1.345e-05 6.926e-09 2.034e-04 2.415e-05 4.849e-06 1.000e+01 5.534e-04
156 291 1.345e-05 6.842e-09 2.037e-04 2.375e-05 4.984e-06 1.000e+01 5.534e-04
157 293 1.344e-05 6.747e-09 2.040e-04 2.336e-05 5.122e-06 1.000e+01 5.534e-04
158 295 1.343e-05 6.622e-09 2.044e-04 2.302e-05 5.266e-06 1.000e+01 5.534e-04
159 297 1.343e-05 6.479e-09 2.038e-04 2.327e-05 5.397e-06 1.000e+01 5.534e-04
160 299 1.342e-05 6.340e-09 2.036e-04 2.342e-05 5.512e-06 1.000e+01 5.534e-04
161 301 1.342e-05 6.168e-09 2.033e-04 2.353e-05 5.625e-06 1.000e+01 5.534e-04
162 303 1.341e-05 5.996e-09 2.030e-04 2.369e-05 5.731e-06 1.000e+01 5.534e-04
163 305 1.340e-05 5.847e-09 2.034e-04 2.368e-05 5.821e-06 1.000e+01 5.534e-04
164 307 1.340e-05 5.666e-09 2.042e-04 2.363e-05 5.916e-06 1.000e+01 5.534e-04
165 309 1.339e-05 5.481e-09 2.042e-04 2.404e-05 5.991e-06 1.000e+01 5.534e-04
166 311 1.339e-05 5.267e-09 2.043e-04 2.458e-05 6.074e-06 1.000e+01 5.534e-04
167 313 1.338e-05 5.034e-09 2.033e-04 2.545e-05 6.138e-06 1.000e+01 5.534e-04
168 315 1.338e-05 4.790e-09 2.023e-04 2.637e-05 6.195e-06 1.000e+01 5.534e-04
169 317 1.337e-05 4.556e-09 2.017e-04 2.739e-05 6.255e-06 1.000e+01 5.534e-04
170 319 1.337e-05 4.345e-09 2.014e-04 2.834e-05 6.311e-06 1.000e+01 5.534e-04
171 321 1.336e-05 4.164e-09 2.012e-04 2.900e-05 6.370e-06 1.000e+01 5.534e-04
172 323 1.336e-05 3.975e-09 2.010e-04 2.959e-05 6.415e-06 1.000e+01 5.534e-04
173 325 1.336e-05 3.825e-09 2.008e-04 3.007e-05 6.459e-06 1.000e+01 5.534e-04
174 327 1.335e-05 3.694e-09 2.009e-04 3.041e-05 6.507e-06 1.000e+01 5.534e-04
175 329 1.335e-05 3.565e-09 2.008e-04 3.056e-05 6.549e-06 1.000e+01 5.534e-04
176 330 1.335e-05 3.449e-09 2.006e-04 3.074e-05 6.570e-06 1.000e+01 5.534e-04
177 331 1.334e-05 3.332e-09 2.003e-04 3.095e-05 6.580e-06 1.000e+01 5.534e-04
178 332 1.334e-05 3.262e-09 2.002e-04 3.125e-05 6.578e-06 1.000e+01 5.534e-04
179 333 1.334e-05 3.216e-09 2.005e-04 3.147e-05 6.584e-06 1.000e+01 5.534e-04
180 334 1.333e-05 3.130e-09 2.010e-04 3.165e-05 6.596e-06 1.000e+01 5.534e-04
181 335 1.333e-05 3.017e-09 2.011e-04 3.194e-05 6.607e-06 1.000e+01 5.534e-04
182 336 1.333e-05 2.939e-09 2.011e-04 3.232e-05 6.610e-06 1.000e+01 5.534e-04
183 337 1.332e-05 2.857e-09 2.008e-04 3.270e-05 6.620e-06 1.000e+01 5.534e-04
184 338 1.332e-05 2.767e-09 2.007e-04 3.330e-05 6.630e-06 1.000e+01 5.534e-04
185 339 1.332e-05 2.714e-09 2.005e-04 3.381e-05 6.641e-06 1.000e+01 5.534e-04
186 340 1.332e-05 2.661e-09 2.003e-04 3.429e-05 6.649e-06 1.000e+01 5.534e-04
187 341 1.331e-05 2.633e-09 2.003e-04 3.469e-05 6.651e-06 1.000e+01 5.534e-04
188 342 1.331e-05 2.571e-09 2.004e-04 3.517e-05 6.655e-06 1.000e+01 5.534e-04
189 343 1.331e-05 2.515e-09 2.006e-04 3.567e-05 6.663e-06 1.000e+01 5.534e-04
190 344 1.331e-05 2.546e-09 2.008e-04 3.578e-05 6.670e-06 1.000e+01 5.534e-04
191 345 1.330e-05 2.511e-09 2.010e-04 3.618e-05 6.677e-06 1.000e+01 5.534e-04
192 346 1.330e-05 2.505e-09 2.012e-04 3.640e-05 6.686e-06 1.000e+01 5.534e-04
193 347 1.330e-05 2.553e-09 2.011e-04 3.626e-05 6.706e-06 1.000e+01 5.534e-04
194 348 1.330e-05 2.501e-09 2.011e-04 3.640e-05 6.722e-06 1.000e+01 5.534e-04
195 349 1.329e-05 2.439e-09 2.014e-04 3.670e-05 6.736e-06 1.000e+01 5.534e-04
196 350 1.329e-05 2.543e-09 2.011e-04 3.615e-05 6.751e-06 1.000e+01 5.534e-04
197 351 1.329e-05 2.502e-09 2.008e-04 3.579e-05 6.768e-06 1.000e+01 5.534e-04
198 352 1.329e-05 2.474e-09 2.007e-04 3.558e-05 6.785e-06 1.000e+01 5.534e-04
199 353 1.328e-05 2.468e-09 2.008e-04 3.566e-05 6.803e-06 1.000e+01 5.534e-04
200 354 1.328e-05 2.522e-09 2.007e-04 3.553e-05 6.823e-06 1.000e+01 5.534e-04
201 355 1.328e-05 2.496e-09 2.010e-04 3.576e-05 6.842e-06 1.000e+01 5.534e-04
202 356 1.328e-05 2.608e-09 2.012e-04 3.527e-05 6.865e-06 1.000e+01 5.534e-04
203 357 1.327e-05 2.615e-09 2.012e-04 3.503e-05 6.884e-06 1.000e+01 5.534e-04
204 358 1.327e-05 2.623e-09 2.010e-04 3.508e-05 6.897e-06 1.000e+01 5.534e-04
205 359 1.327e-05 2.716e-09 2.006e-04 3.455e-05 6.914e-06 1.000e+01 5.534e-04
206 360 1.326e-05 2.719e-09 2.007e-04 3.423e-05 6.930e-06 1.000e+01 5.534e-04
207 361 1.326e-05 2.726e-09 2.004e-04 3.427e-05 6.953e-06 1.000e+01 5.534e-04
208 362 1.326e-05 2.783e-09 2.000e-04 3.388e-05 6.983e-06 1.000e+01 5.534e-04
209 363 1.326e-05 2.732e-09 2.000e-04 3.395e-05 7.009e-06 1.000e+01 5.534e-04
210 364 1.325e-05 2.779e-09 2.000e-04 3.379e-05 7.034e-06 1.000e+01 5.534e-04
211 365 1.325e-05 2.763e-09 2.001e-04 3.379e-05 7.056e-06 1.000e+01 5.534e-04
212 366 1.325e-05 2.781e-09 2.002e-04 3.368e-05 7.080e-06 1.000e+01 5.534e-04
213 367 1.325e-05 2.774e-09 2.002e-04 3.353e-05 7.099e-06 1.000e+01 5.534e-04
214 368 1.324e-05 2.752e-09 2.001e-04 3.367e-05 7.106e-06 1.000e+01 5.534e-04
215 369 1.324e-05 2.743e-09 2.002e-04 3.402e-05 7.100e-06 1.000e+01 5.534e-04
216 370 1.324e-05 2.800e-09 2.000e-04 3.379e-05 7.094e-06 1.000e+01 5.534e-04
217 371 1.323e-05 2.768e-09 1.992e-04 3.361e-05 7.093e-06 1.000e+01 5.534e-04
218 372 1.323e-05 2.738e-09 1.985e-04 3.332e-05 7.102e-06 1.000e+01 5.534e-04
219 373 1.323e-05 2.697e-09 1.980e-04 3.308e-05 7.110e-06 1.000e+01 5.534e-04
220 374 1.323e-05 2.643e-09 1.979e-04 3.328e-05 7.121e-06 1.000e+01 5.534e-04
221 375 1.322e-05 2.635e-09 1.981e-04 3.348e-05 7.130e-06 1.000e+01 5.534e-04
222 376 1.322e-05 2.605e-09 1.983e-04 3.399e-05 7.133e-06 1.000e+01 5.534e-04
223 377 1.322e-05 2.566e-09 1.989e-04 3.484e-05 7.128e-06 1.000e+01 5.534e-04
224 378 1.322e-05 2.604e-09 1.994e-04 3.555e-05 7.119e-06 1.000e+01 5.534e-04
225 379 1.321e-05 2.692e-09 1.995e-04 3.603e-05 7.101e-06 1.000e+01 5.534e-04
226 380 1.321e-05 2.781e-09 1.990e-04 3.626e-05 7.068e-06 1.000e+01 5.534e-04
227 381 1.321e-05 2.844e-09 1.983e-04 3.628e-05 7.038e-06 1.000e+01 5.534e-04
228 382 1.320e-05 2.945e-09 1.978e-04 3.603e-05 7.002e-06 1.000e+01 5.534e-04
229 383 1.320e-05 3.091e-09 1.971e-04 3.593e-05 6.957e-06 1.000e+01 5.534e-04
230 384 1.320e-05 3.338e-09 1.957e-04 3.549e-05 6.905e-06 1.000e+01 5.534e-04
231 385 1.319e-05 3.657e-09 1.937e-04 3.459e-05 6.838e-06 1.000e+01 5.534e-04
232 386 1.319e-05 4.098e-09 1.909e-04 3.316e-05 6.753e-06 1.000e+01 5.534e-04
233 388 1.319e-05 4.743e-09 1.866e-04 3.079e-05 6.654e-06 1.000e+01 5.534e-04
234 390 1.318e-05 5.558e-09 1.807e-04 2.747e-05 6.533e-06 1.000e+01 5.534e-04
235 392 1.317e-05 6.548e-09 1.721e-04 2.317e-05 6.403e-06 1.000e+01 5.534e-04
236 394 1.317e-05 7.518e-09 1.627e-04 1.877e-05 6.268e-06 1.000e+01 5.534e-04
237 396 1.316e-05 8.370e-09 1.545e-04 1.487e-05 6.134e-06 1.000e+01 5.534e-04
238 398 1.315e-05 8.942e-09 1.484e-04 1.178e-05 6.006e-06 1.000e+01 5.534e-04
239 400 1.314e-05 8.938e-09 1.461e-04 1.119e-05 5.887e-06 1.000e+01 5.534e-04
240 402 1.313e-05 8.394e-09 1.470e-04 1.055e-05 5.787e-06 1.000e+01 5.534e-04
241 404 1.312e-05 7.450e-09 1.499e-04 9.603e-06 5.710e-06 1.000e+01 5.534e-04
242 406 1.312e-05 6.334e-09 1.544e-04 8.654e-06 5.657e-06 1.000e+01 5.534e-04
243 408 1.311e-05 5.210e-09 1.614e-04 8.583e-06 5.625e-06 1.000e+01 5.534e-04
244 410 1.311e-05 4.126e-09 1.728e-04 1.138e-05 5.621e-06 1.000e+01 5.534e-04
245 412 1.311e-05 3.079e-09 1.918e-04 1.483e-05 5.636e-06 1.000e+01 5.534e-04
246 414 1.310e-05 2.265e-09 2.128e-04 1.845e-05 5.640e-06 1.000e+01 5.534e-04
247 415 1.310e-05 1.612e-09 2.449e-04 2.412e-05 5.666e-06 1.000e+01 5.534e-04
248 416 1.310e-05 2.507e-09 6.512e-05 1.648e-06 5.660e-06 1.000e+01 5.534e-04
249 418 1.310e-05 4.816e-10 6.733e-05 1.766e-06 5.662e-06 1.000e+01 5.534e-04
250 420 1.310e-05 4.580e-10 6.928e-05 1.841e-06 5.665e-06 1.000e+01 5.534e-04
251 422 1.310e-05 4.409e-10 7.097e-05 1.900e-06 5.670e-06 1.000e+01 5.534e-04
252 424 1.310e-05 4.275e-10 7.243e-05 1.960e-06 5.677e-06 1.000e+01 5.534e-04
253 426 1.310e-05 4.155e-10 7.365e-05 2.021e-06 5.685e-06 1.000e+01 5.534e-04
254 428 1.310e-05 4.051e-10 7.455e-05 2.081e-06 5.691e-06 1.000e+01 5.534e-04
255 430 1.310e-05 3.960e-10 7.520e-05 2.144e-06 5.694e-06 1.000e+01 5.534e-04
256 432 1.310e-05 3.859e-10 7.588e-05 2.203e-06 5.706e-06 1.000e+01 5.534e-04
257 434 1.309e-05 3.793e-10 7.646e-05 2.255e-06 5.723e-06 1.000e+01 5.534e-04
258 436 1.309e-05 3.767e-10 7.692e-05 2.300e-06 5.740e-06 1.000e+01 5.534e-04
259 438 1.309e-05 3.750e-10 7.729e-05 2.319e-06 5.755e-06 1.000e+01 5.534e-04
260 440 1.309e-05 3.735e-10 7.773e-05 2.337e-06 5.774e-06 1.000e+01 5.534e-04
261 442 1.309e-05 3.713e-10 7.794e-05 2.372e-06 5.797e-06 1.000e+01 5.534e-04
262 444 1.309e-05 3.734e-10 7.802e-05 2.391e-06 5.815e-06 1.000e+01 5.534e-04
263 446 1.309e-05 3.745e-10 7.776e-05 2.391e-06 5.825e-06 1.000e+01 5.534e-04
264 448 1.309e-05 3.734e-10 7.754e-05 2.390e-06 5.835e-06 1.000e+01 5.534e-04
265 450 1.309e-05 3.730e-10 7.733e-05 2.387e-06 5.843e-06 1.000e+01 5.534e-04
266 452 1.309e-05 3.771e-10 7.743e-05 2.391e-06 5.845e-06 1.000e+01 5.534e-04
267 454 1.309e-05 3.792e-10 7.759e-05 2.398e-06 5.846e-06 1.000e+01 5.534e-04
268 456 1.309e-05 3.794e-10 7.767e-05 2.410e-06 5.847e-06 1.000e+01 5.534e-04
269 458 1.309e-05 3.757e-10 7.769e-05 2.435e-06 5.851e-06 1.000e+01 5.534e-04
270 460 1.309e-05 3.745e-10 7.777e-05 2.464e-06 5.858e-06 1.000e+01 5.534e-04
271 462 1.309e-05 3.754e-10 7.794e-05 2.495e-06 5.864e-06 1.000e+01 5.534e-04
272 464 1.309e-05 3.787e-10 7.821e-05 2.541e-06 5.874e-06 1.000e+01 5.534e-04
273 466 1.309e-05 3.828e-10 7.862e-05 2.598e-06 5.888e-06 1.000e+01 5.534e-04
274 468 1.309e-05 3.867e-10 7.888e-05 2.637e-06 5.908e-06 1.000e+01 5.534e-04
275 470 1.309e-05 3.919e-10 7.914e-05 2.680e-06 5.934e-06 1.000e+01 5.534e-04
276 472 1.309e-05 4.006e-10 7.949e-05 2.723e-06 5.964e-06 1.000e+01 5.534e-04
277 474 1.309e-05 4.117e-10 7.983e-05 2.768e-06 5.996e-06 1.000e+01 5.534e-04
278 476 1.309e-05 4.244e-10 8.009e-05 2.816e-06 6.026e-06 1.000e+01 5.534e-04
279 478 1.309e-05 4.382e-10 8.034e-05 2.847e-06 6.055e-06 1.000e+01 5.534e-04
280 480 1.309e-05 4.546e-10 8.065e-05 2.858e-06 6.081e-06 1.000e+01 5.534e-04
281 482 1.309e-05 4.671e-10 8.090e-05 2.873e-06 6.106e-06 1.000e+01 5.534e-04
282 484 1.309e-05 4.817e-10 8.118e-05 2.900e-06 6.133e-06 1.000e+01 5.534e-04
283 486 1.308e-05 5.077e-10 8.147e-05 2.931e-06 6.158e-06 1.000e+01 5.534e-04
284 488 1.308e-05 5.372e-10 8.205e-05 2.935e-06 6.191e-06 1.000e+01 5.534e-04
285 490 1.308e-05 5.710e-10 8.245e-05 2.943e-06 6.227e-06 1.000e+01 5.534e-04
286 492 1.308e-05 6.049e-10 8.278e-05 2.951e-06 6.268e-06 1.000e+01 5.534e-04
287 494 1.308e-05 6.398e-10 8.305e-05 2.951e-06 6.313e-06 1.000e+01 5.534e-04
288 496 1.308e-05 6.768e-10 8.325e-05 2.946e-06 6.359e-06 1.000e+01 5.534e-04
289 498 1.308e-05 7.136e-10 8.341e-05 2.937e-06 6.406e-06 1.000e+01 5.534e-04
290 500 1.308e-05 7.471e-10 8.353e-05 2.927e-06 6.453e-06 1.000e+01 5.534e-04
291 502 1.308e-05 7.779e-10 8.342e-05 2.920e-06 6.493e-06 1.000e+01 5.534e-04
292 504 1.308e-05 8.030e-10 8.331e-05 2.907e-06 6.532e-06 1.000e+01 5.534e-04
293 506 1.308e-05 8.234e-10 8.318e-05 2.897e-06 6.572e-06 1.000e+01 5.534e-04
294 508 1.308e-05 8.441e-10 8.307e-05 2.884e-06 6.613e-06 1.000e+01 5.534e-04
295 510 1.308e-05 8.657e-10 8.296e-05 2.871e-06 6.649e-06 1.000e+01 5.534e-04
296 512 1.308e-05 8.831e-10 8.286e-05 2.860e-06 6.687e-06 1.000e+01 5.534e-04
297 514 1.307e-05 8.963e-10 8.274e-05 2.850e-06 6.731e-06 1.000e+01 5.534e-04
298 516 1.307e-05 9.111e-10 8.264e-05 2.838e-06 6.777e-06 1.000e+01 5.534e-04
299 518 1.307e-05 9.254e-10 8.257e-05 2.823e-06 6.824e-06 1.000e+01 5.534e-04
300 520 1.307e-05 9.393e-10 8.257e-05 2.801e-06 6.874e-06 1.000e+01 5.534e-04
301 522 1.307e-05 9.539e-10 8.257e-05 2.779e-06 6.925e-06 1.000e+01 5.534e-04
302 524 1.307e-05 9.702e-10 8.254e-05 2.755e-06 6.973e-06 1.000e+01 5.534e-04
303 526 1.307e-05 9.817e-10 8.256e-05 2.729e-06 7.020e-06 1.000e+01 5.534e-04
304 528 1.307e-05 9.936e-10 8.239e-05 2.708e-06 7.054e-06 1.000e+01 5.534e-04
305 530 1.307e-05 9.994e-10 8.228e-05 2.683e-06 7.087e-06 1.000e+01 5.534e-04
306 532 1.307e-05 1.009e-09 8.223e-05 2.652e-06 7.128e-06 1.000e+01 5.534e-04
307 533 1.307e-05 -7.987e-11 3.244e-04 3.863e-05 7.312e-06 1.000e+01 5.534e-04
308 534 1.306e-05 5.242e-09 8.301e-05 2.532e-06 7.238e-06 1.000e+01 5.534e-04
309 535 1.306e-05 2.055e-10 3.250e-04 3.662e-05 7.368e-06 1.000e+01 5.534e-04
310 536 1.306e-05 5.047e-09 8.306e-05 2.348e-06 7.241e-06 1.000e+01 5.534e-04
311 537 1.305e-05 3.240e-10 3.238e-04 3.578e-05 7.372e-06 1.000e+01 5.534e-04
312 538 1.305e-05 4.788e-09 8.258e-05 2.139e-06 7.262e-06 1.000e+01 5.534e-04
313 539 1.305e-05 2.853e-10 3.207e-04 3.811e-05 7.396e-06 1.000e+01 5.534e-04
314 540 1.305e-05 4.621e-09 8.220e-05 2.231e-06 7.264e-06 1.000e+01 5.534e-04
315 541 1.304e-05 9.738e-11 3.171e-04 3.988e-05 7.392e-06 1.000e+01 5.534e-04
316 542 1.304e-05 4.465e-09 8.148e-05 2.451e-06 7.250e-06 1.000e+01 5.534e-04
317 544 1.304e-05 8.692e-10 8.065e-05 2.546e-06 7.235e-06 1.000e+01 5.534e-04
318 546 1.304e-05 8.477e-10 8.010e-05 2.612e-06 7.227e-06 1.000e+01 5.534e-04
319 548 1.304e-05 8.291e-10 7.975e-05 2.651e-06 7.225e-06 1.000e+01 5.534e-04
320 550 1.304e-05 8.143e-10 7.949e-05 2.678e-06 7.224e-06 1.000e+01 5.534e-04
321 552 1.304e-05 7.990e-10 7.930e-05 2.698e-06 7.224e-06 1.000e+01 5.534e-04
322 554 1.304e-05 7.848e-10 7.908e-05 2.708e-06 7.220e-06 1.000e+01 5.534e-04
323 556 1.303e-05 7.696e-10 7.889e-05 2.714e-06 7.217e-06 1.000e+01 5.534e-04
324 558 1.303e-05 7.538e-10 7.878e-05 2.721e-06 7.214e-06 1.000e+01 5.534e-04
325 560 1.303e-05 7.374e-10 7.867e-05 2.727e-06 7.218e-06 1.000e+01 5.534e-04
326 562 1.303e-05 7.256e-10 7.865e-05 2.722e-06 7.220e-06 1.000e+01 5.534e-04
327 564 1.303e-05 7.107e-10 7.854e-05 2.696e-06 7.225e-06 1.000e+01 5.534e-04
328 566 1.303e-05 6.852e-10 7.826e-05 2.670e-06 7.246e-06 1.000e+01 5.534e-04
329 568 1.303e-05 6.671e-10 7.785e-05 2.642e-06 7.267e-06 1.000e+01 5.534e-04
330 570 1.303e-05 6.504e-10 7.753e-05 2.618e-06 7.287e-06 1.000e+01 5.534e-04
331 572 1.303e-05 6.317e-10 7.716e-05 2.591e-06 7.308e-06 1.000e+01 5.534e-04
332 574 1.303e-05 6.133e-10 7.665e-05 2.616e-06 7.333e-06 1.000e+01 5.534e-04
333 576 1.303e-05 6.028e-10 7.641e-05 2.637e-06 7.356e-06 1.000e+01 5.534e-04
334 578 1.303e-05 5.897e-10 7.617e-05 2.655e-06 7.382e-06 1.000e+01 5.534e-04
335 580 1.303e-05 5.767e-10 7.581e-05 2.663e-06 7.408e-06 1.000e+01 5.534e-04
336 582 1.303e-05 5.604e-10 7.544e-05 2.655e-06 7.439e-06 1.000e+01 5.534e-04
337 584 1.303e-05 5.428e-10 7.488e-05 2.649e-06 7.473e-06 1.000e+01 5.534e-04
338 586 1.303e-05 5.260e-10 7.427e-05 2.644e-06 7.507e-06 1.000e+01 5.534e-04
339 588 1.302e-05 5.086e-10 7.357e-05 2.643e-06 7.540e-06 1.000e+01 5.534e-04
340 590 1.302e-05 4.844e-10 7.271e-05 2.640e-06 7.583e-06 1.000e+01 5.534e-04
341 592 1.302e-05 4.660e-10 7.180e-05 2.654e-06 7.624e-06 1.000e+01 5.534e-04
342 594 1.302e-05 4.473e-10 7.089e-05 2.656e-06 7.665e-06 1.000e+01 5.534e-04
343 596 1.302e-05 4.263e-10 6.983e-05 2.668e-06 7.706e-06 1.000e+01 5.534e-04
344 598 1.302e-05 4.054e-10 6.863e-05 2.693e-06 7.748e-06 1.000e+01 5.534e-04
345 600 1.302e-05 3.859e-10 6.751e-05 2.704e-06 7.789e-06 1.000e+01 5.534e-04
346 602 1.302e-05 3.678e-10 6.636e-05 2.698e-06 7.826e-06 1.000e+01 5.534e-04
347 604 1.302e-05 3.490e-10 6.521e-05 2.679e-06 7.863e-06 1.000e+01 5.534e-04
348 606 1.302e-05 3.304e-10 6.404e-05 2.646e-06 7.899e-06 1.000e+01 5.534e-04
349 608 1.302e-05 3.127e-10 6.286e-05 2.599e-06 7.936e-06 1.000e+01 5.534e-04
350 610 1.302e-05 2.960e-10 6.168e-05 2.538e-06 7.971e-06 1.000e+01 5.534e-04
351 612 1.302e-05 2.804e-10 6.054e-05 2.471e-06 8.007e-06 1.000e+01 5.534e-04
352 614 1.302e-05 2.661e-10 5.949e-05 2.399e-06 8.041e-06 1.000e+01 5.534e-04
353 616 1.302e-05 2.530e-10 5.844e-05 2.353e-06 8.074e-06 1.000e+01 5.534e-04
354 618 1.302e-05 2.433e-10 5.748e-05 2.311e-06 8.102e-06 1.000e+01 5.534e-04
355 620 1.302e-05 2.318e-10 5.672e-05 2.309e-06 8.129e-06 1.000e+01 5.534e-04
356 622 1.302e-05 2.224e-10 5.550e-05 2.226e-06 8.149e-06 1.000e+01 5.534e-04
357 624 1.302e-05 2.112e-10 5.422e-05 2.117e-06 8.167e-06 1.000e+01 5.534e-04
358 626 1.302e-05 2.013e-10 5.311e-05 2.021e-06 8.184e-06 1.000e+01 5.534e-04
359 628 1.302e-05 1.940e-10 5.239e-05 1.994e-06 8.196e-06 1.000e+01 5.534e-04
360 630 1.302e-05 1.864e-10 5.171e-05 1.964e-06 8.207e-06 1.000e+01 5.534e-04
361 632 1.302e-05 1.782e-10 5.101e-05 1.912e-06 8.218e-06 1.000e+01 5.534e-04
362 634 1.302e-05 1.708e-10 5.037e-05 1.853e-06 8.228e-06 1.000e+01 5.534e-04
363 636 1.302e-05 1.636e-10 4.989e-05 1.794e-06 8.239e-06 1.000e+01 5.534e-04
364 638 1.302e-05 1.554e-10 4.989e-05 1.837e-06 8.258e-06 1.000e+01 5.534e-04
365 640 1.302e-05 1.524e-10 5.004e-05 1.930e-06 8.274e-06 1.000e+01 5.534e-04
366 642 1.302e-05 1.492e-10 5.021e-05 2.023e-06 8.290e-06 1.000e+01 5.534e-04
367 644 1.302e-05 1.464e-10 5.036e-05 2.103e-06 8.305e-06 1.000e+01 5.534e-04
368 646 1.302e-05 1.441e-10 5.074e-05 2.213e-06 8.321e-06 1.000e+01 5.534e-04
369 648 1.302e-05 1.421e-10 5.118e-05 2.332e-06 8.337e-06 1.000e+01 5.534e-04
370 650 1.302e-05 1.406e-10 5.155e-05 2.432e-06 8.354e-06 1.000e+01 5.534e-04
371 652 1.302e-05 1.392e-10 5.188e-05 2.516e-06 8.370e-06 1.000e+01 5.534e-04
372 654 1.302e-05 1.388e-10 5.222e-05 2.601e-06 8.387e-06 1.000e+01 5.534e-04
373 656 1.302e-05 1.375e-10 5.273e-05 2.719e-06 8.405e-06 1.000e+01 5.534e-04
374 658 1.302e-05 1.382e-10 5.243e-05 2.690e-06 8.418e-06 1.000e+01 5.534e-04
375 660 1.302e-05 1.370e-10 5.213e-05 2.636e-06 8.428e-06 1.000e+01 5.534e-04
376 662 1.302e-05 1.357e-10 5.200e-05 2.602e-06 8.436e-06 1.000e+01 5.534e-04
377 664 1.302e-05 1.338e-10 5.199e-05 2.592e-06 8.443e-06 1.000e+01 5.534e-04
378 666 1.302e-05 1.323e-10 5.205e-05 2.594e-06 8.451e-06 1.000e+01 5.534e-04
379 668 1.302e-05 1.309e-10 5.216e-05 2.605e-06 8.457e-06 1.000e+01 5.534e-04
380 670 1.302e-05 1.294e-10 5.230e-05 2.621e-06 8.465e-06 1.000e+01 5.534e-04
Optimization terminated successfully.
`ftol` condition satisfied. (ftol=1.00e-05)
Current function value: 1.302e-05
Constraint violation: 8.465e-06
Total delta_x: 9.170e-02
Iterations: 380
Function evaluations: 670
Jacobian evaluations: 383
Timer: Solution time = 1.28 min
Timer: Avg time per step = 202 ms
==============================================================================================================
Start --> End
Total (sum of squares): 6.632e+01 --> 1.302e-05,
Maximum absolute QFM surface normal field error: 3.970e-02 --> 5.304e-05 (T m^2)
Minimum absolute QFM surface normal field error: 9.255e-18 --> 1.113e-17 (T m^2)
Average absolute QFM surface normal field error: 1.676e-02 --> 4.210e-06 (T m^2)
Maximum absolute QFM surface normal field error: 5.672e-01 --> 7.577e-04 (normalized)
Minimum absolute QFM surface normal field error: 1.322e-16 --> 1.590e-16 (normalized)
Average absolute QFM surface normal field error: 2.394e-01 --> 6.014e-05 (normalized)
Maximum Toroidal Flux: -2.117e-02 --> -3.501e-02 (Wb)
Minimum Toroidal Flux: -2.117e-02 --> -3.501e-02 (Wb)
Average Toroidal Flux: -2.117e-02 --> -3.501e-02 (Wb)
Maximum Toroidal Flux: -2.117e-02 --> -3.501e-02 (normalized)
Minimum Toroidal Flux: -2.117e-02 --> -3.501e-02 (normalized)
Average Toroidal Flux: -2.117e-02 --> -3.501e-02 (normalized)
==============================================================================================================
We can see that after the optimization, the normal field error is quite small
[8]:
plot_2d(qfm_surf, "B*n", field=field, cmap="viridis", log=True);
We also see that the optimized surface is now a good approximation of a flux surface in this vacuum field
[9]:
data = qfm_surf.compute(["R", "Z"], grid=LinearGrid(rho=1.0, theta=0, zeta=0))
fig, ax = plot_surfaces(qfm_surf, rho=1.0, theta=0)
poincare_plot(
field, R0=np.append(data["R"], np.linspace(0.7, 0.8, 9)), Z0=np.zeros(10), ax=ax
);
We can use this QFM surface as the surface for a fixed-boundary DESC equilibrium solve, to find the equilibrium field which matches the vacuum field inside this surface. (Alternatively, one could also perform a free-boundary solve, though we know the surface should not change much from this QFM surface)
[10]:
from desc.equilibrium import Equilibrium
eq = Equilibrium(surface=qfm_surf, Psi=-0.035, L=8, M=6, N=6)
eq.solve(verbose=3, ftol=1e-8);
Building objective: force
Precomputing transforms
Timer: Precomputing transforms = 805 ms
Timer: Objective build = 1.14 sec
Building objective: lcfs R
Building objective: lcfs Z
Building objective: fixed Psi
Building objective: fixed pressure
Building objective: fixed current
Building objective: fixed sheet current
Building objective: self_consistency R
Building objective: self_consistency Z
Building objective: lambda gauge
Building objective: axis R self consistency
Building objective: axis Z self consistency
Timer: Objective build = 613 ms
Timer: LinearConstraintProjection build = 5.00 sec
Number of parameters: 628
Number of objectives: 3250
Timer: Initializing the optimization = 6.86 sec
Starting optimization
Using method: lsq-exact
Solver options:
------------------------------------------------------------
Maximum Function Evaluations : 501
Maximum Allowed Total Δx Norm : inf
Scaled Termination : True
Trust Region Method : qr
Initial Trust Radius : 4.264e+01
Maximum Trust Radius : inf
Minimum Trust Radius : 2.220e-16
Trust Radius Increase Ratio : 2.000e+00
Trust Radius Decrease Ratio : 2.500e-01
Trust Radius Increase Threshold : 7.500e-01
Trust Radius Decrease Threshold : 2.500e-01
------------------------------------------------------------
Iteration Total nfev Cost Cost reduction Step norm Optimality
0 1 3.949e-01 6.805e-01
1 2 2.375e-02 3.711e-01 4.824e-01 5.959e-02
2 3 6.834e-04 2.307e-02 1.043e-01 9.521e-03
3 4 2.552e-05 6.579e-04 7.154e-02 2.345e-03
4 5 1.664e-05 8.882e-06 5.916e-02 1.543e-03
5 6 5.739e-06 1.090e-05 4.538e-02 1.006e-03
6 7 2.633e-07 5.476e-06 1.325e-02 1.555e-04
7 8 1.160e-07 1.473e-07 1.331e-02 9.814e-05
8 9 4.818e-08 6.786e-08 1.989e-03 7.832e-06
9 10 4.781e-08 3.649e-10 4.063e-03 3.204e-06
10 11 4.773e-08 7.782e-11 5.620e-04 4.549e-07
11 12 4.773e-08 2.613e-12 1.619e-03 3.917e-07
12 13 4.773e-08 1.685e-12 3.555e-04 9.319e-08
13 14 4.773e-08 2.860e-13 6.708e-04 5.214e-08
14 15 4.773e-08 1.322e-13 2.244e-04 2.553e-08
15 16 4.773e-08 5.043e-14 2.973e-04 1.175e-08
16 17 4.773e-08 2.478e-14 1.352e-04 6.736e-09
`gtol` condition satisfied. (gtol=1.00e-08)
Current function value: 4.773e-08
Total delta_x: 4.076e-01
Iterations: 16
Function evaluations: 17
Jacobian evaluations: 17
Timer: Solution time = 20.6 sec
Timer: Avg time per step = 1.21 sec
==============================================================================================================
Start --> End
Total (sum of squares): 3.949e-01 --> 4.773e-08,
Maximum absolute Force error: 9.665e+04 --> 5.936e+01 (N)
Minimum absolute Force error: 1.263e-01 --> 1.017e-04 (N)
Average absolute Force error: 4.629e+03 --> 1.594e+00 (N)
Maximum absolute Force error: 2.254e-01 --> 1.385e-04 (normalized)
Minimum absolute Force error: 2.946e-07 --> 2.372e-10 (normalized)
Average absolute Force error: 1.080e-02 --> 3.718e-06 (normalized)
R boundary error: 0.000e+00 --> 1.395e-17 (m)
Z boundary error: 0.000e+00 --> 1.813e-18 (m)
Fixed Psi error: 0.000e+00 --> 0.000e+00 (Wb)
Fixed pressure profile error: 0.000e+00 --> 0.000e+00 (Pa)
Fixed current profile error: 0.000e+00 --> 0.000e+00 (A)
Fixed sheet current error: 0.000e+00 --> 0.000e+00 (~)
==============================================================================================================
Finally, we compare the solve equilibrium’s flux surfaces with the vacuum field Poincare plot and confirm that indeed, the equilibrium matches the vacuum field as expected.
[11]:
grid = LinearGrid(L=5)
data = eq.compute(["R", "Z"], grid=grid)
fig, ax = plot_surfaces(eq, rho=grid.nodes[:, 0], rho_lw=2)
poincare_plot(field, R0=data["R"], Z0=data["Z"], ax=ax, NFP=eq.NFP, size=10);