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,

\[\begin{equation} S = \arg\min_S \int_S |B_n|^2 dS \end{equation}\]

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,
);
../../_images/notebooks_tutorials_QFM_surface_10_0.png

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);
../../_images/notebooks_tutorials_QFM_surface_14_0.png

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
);
../../_images/notebooks_tutorials_QFM_surface_16_0.png

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);
../../_images/notebooks_tutorials_QFM_surface_20_0.png