Functions/Subroutines
subroutine	pc00aa (NGdof, position, Nvol, mn, ie04dgf)
	Use preconditioned conjugate gradient method to find minimum of energy functional. More...

subroutine	pc00ab (mode, NGdof, Position, Energy, Gradient, nstate, iuser, ruser)
	Returns the energy functional and it's derivatives with respect to geometry. More...

Detailed Description

Function/Subroutine Documentation

◆ pc00aa()

subroutine pc00aa	(	integer, intent(in)	NGdof,
		real, dimension(0:ngdof), intent(inout)	position,
		integer, intent(in)	Nvol,
		integer, intent(in)	mn,
		integer	ie04dgf
	)

Use preconditioned conjugate gradient method to find minimum of energy functional.

energy functional

The energy functional is described in pc00ab() .

relevant input variables

The following input variables control the operation of E04DGF :
- epsilon : weighting of "spectral energy"; see pc00ab()
- maxstep : this is given to E04DGF for the \(\texttt{Maximum Step Length}\)
- maxiter : upper limit on derivative calculations used in the conjugate gradient iterations
- verify : if verify=1, then E04DGF will confirm user supplied gradients (provided by pc00ab() ) are correct;
Todo:
Unfortunately, E04DGF seems to require approximately \(3 N\) function evaluations before proceeding to minimize the energy functional, where there are \(N\) degrees of freedom. I don't know how to turn this off!

Parameters

[in]	NGdof
[in,out]	position
[in]	Nvol
[in]	mn
	ie04dgf

References allglobal::cpus, allglobal::energy, allglobal::forceerr, inputlist::forcetol, allglobal::myid, allglobal::ncpu, fileunits::ounit, pc00aa(), pc00ab(), constants::ten, and constants::zero.

Referenced by pc00aa().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ pc00ab()

subroutine pc00ab	(	integer	mode,
		integer	NGdof,
		real, dimension(1:ngdof)	Position,
		real	Energy,
		real, dimension(1:ngdof)	Gradient,
		integer	nstate,
		integer, dimension(1:2)	iuser,
		real, dimension(1:1)	ruser
	)

Returns the energy functional and it's derivatives with respect to geometry.

Energy functional

The energy functional is
\begin{eqnarray} F \equiv \sum_{l=1}^{N} \int_{\cal V} \left( \frac{p}{\gamma-1} + \frac{B^2}{2} \right) dv, \label{eq:energyfunctional_pc00ab} \end{eqnarray}
where \(N \equiv\,\)Nvol is the number of interfaces.
Assuming that the toroidal and poloidal fluxes, \(\psi_t\) and \(\psi_p\), the helicity, \({\cal K}\), the helicity multiplier, \(\mu\), and/or the interface rotational-transforms, \({{\,\iota\!\!\!}-}\), are appropriately constrained, the Beltrami fields in each volume depend only the geometry of the adjacent interfaces. So, the energy functional is assumed to be a function of "position", i.e. \(F = F(R_{l,j},Z_{l,j})\).
Introducing a ficitious time, \(t\), the position may be advanced according to
\begin{eqnarray} \begin{array}{cccccccccccccccccccccccc} \displaystyle \frac{\partial R_j}{\partial t} & \equiv & \displaystyle - \frac{\partial }{\partial R_j} \sum_{l=1}^{N} \int \left( \frac{p}{\gamma-1} + \frac{B^2}{2} \right) dv,\\ \displaystyle \frac{\partial Z_j}{\partial t} & \equiv & \displaystyle - \frac{\partial }{\partial Z_j} \sum_{l=1}^{N} \int \left( \frac{p}{\gamma-1} + \frac{B^2}{2} \right) dv. \end{array} \label{eq:descent_pc00ab} \end{eqnarray}
There remain degrees of freedom in the angle representation of the interfaces.

Spectral energy minimization

Consider variations which do not affect the geometry of the surfaces,
\begin{eqnarray} \delta R &=& R_\theta \; u,\\ \delta Z &=& Z_\theta \; u, \end{eqnarray}
where \(u\) is a angle variation.
The corresponding variation in each of the Fourier harmonics is
\begin{eqnarray} \delta R_j &\equiv& \oint\!\!\!\oint \!d\theta d\zeta \,\,\, R_\theta \; u \; \cos \alpha_j,\\ \delta Z_j &\equiv& \oint\!\!\!\oint \!d\theta d\zeta \,\,\, Z_\theta \; u \; \sin \alpha_j, \end{eqnarray}
Following Hirshman et al., introducing the normalized spectral width
\begin{eqnarray} M \equiv \frac{\sum_j ( m_j^p + n_j^q ) ( R_{l,j}^2+Z_{l,j}^2 )}{\sum_j ( R_{l,j}^2+Z_{l,j}^2 )}, \end{eqnarray}
Using the notation
\begin{eqnarray} N &\equiv& \sum_j \lambda_j ( R_{l,j}^2+Z_{l,j}^2 ), \\ D &\equiv& \sum_j ( R_{l,j}^2+Z_{l,j}^2 ), \end{eqnarray}
where \(\lambda_j \equiv m_j^p + n_j^q\), the variation in the normalized spectral width is
\begin{eqnarray} \delta M = (\delta N - M \delta D)/D. \end{eqnarray}
For tangential variations,
\begin{eqnarray} \delta N &=& 2 \oint\!\!\!\oint \!d\theta d\zeta \,\,\, \; u \left(R_\theta \sum_j \lambda_j R_j \cos \alpha_j + Z_\theta \sum_j \lambda_j Z_j \sin \alpha_j \right),\\ \delta D &=& 2 \oint\!\!\!\oint \!d\theta d\zeta \,\,\, \; u \left(R_\theta \sum_j R_j \cos \alpha_j + Z_\theta \sum_j Z_j \sin \alpha_j \right). \end{eqnarray}
The "tangential spectral-width descent direction" is thus
\begin{eqnarray} \frac{\partial u}{\partial t} &=& -\left[ R_\theta \sum_j(\lambda_j - M) R_j \cos \alpha_j / D + Z_\theta \sum_j(\lambda_j - M)Z_j \sin \alpha_j / D \right]. \end{eqnarray}
This suggests that position should be advanced according to
\begin{eqnarray} \frac{\partial R_j}{\partial t} & \equiv & - \frac{\partial }{\partial R_j} \sum_{l=1}^{N} \int \left( \frac{p}{\gamma-1} + \frac{B^2}{2} \right) dv - [R_\theta (R_\theta X + Z_\theta Y)]_j,\\ \frac{\partial Z_j}{\partial t} & \equiv & - \frac{\partial }{\partial Z_j} \sum_{l=1}^{N} \int \left( \frac{p}{\gamma-1} + \frac{B^2}{2} \right) dv - [Z_\theta (R_\theta X + Z_\theta Y)]_j, \end{eqnarray}
where \(X \equiv \sum_j (\lambda_j - M)R_j \cos\alpha_j / D\) and \(Y \equiv \sum_j (\lambda_j - M)Z_j \sin\alpha_j / D\).

numerical implementation

The spectral condensation terms,
\begin{eqnarray} R_\theta (R_\theta X + Z_\theta Y) &=& \sum_{j,k,l} m_j m_k (\lambda_l-M) R_j ( + R_k R_l \sin\alpha_j\sin\alpha_k\cos\alpha_l - Z_k Z_l \sin\alpha_j\cos\alpha_k\sin\alpha_l)/D,\\ Z_\theta (R_\theta X + Z_\theta Y) &=& \sum_{j,k,l} m_j m_k (\lambda_l-M) Z_j ( - R_k R_l \cos\alpha_j\sin\alpha_k\cos\alpha_l + Z_k Z_l \cos\alpha_j\cos\alpha_k\sin\alpha_l)/D, \end{eqnarray}

are calculated using triple angle expressions...
Todo:
IT IS VERY LIKELY THAT FFTs WOULD BE FASTER!!!

References allglobal::cpus, allglobal::dbbdrz, dforce(), allglobal::diidrz, inputlist::epsilon, allglobal::forceerr, inputlist::forcetol, constants::half, inputlist::igeometry, allglobal::lbbintegral, allglobal::mn, allglobal::myid, inputlist::nvol, constants::one, fileunits::ounit, pc00ab(), allglobal::yesstellsym, and constants::zero.