"local" force

Functions/Subroutines
subroutine	lforce (lvol, iocons, ideriv, Ntz, dBB, XX, YY, length, DDl, MMl, iflag)
	Computes \(B^2\), and the spectral condensation constraints if required, on the interfaces, \({\cal I}_i\). More...

Detailed Description

Function/Subroutine Documentation

lforce()

subroutine lforce	(	integer, intent(in)	lvol,
		integer, intent(in)	iocons,
		integer, intent(in)	ideriv,
		integer, intent(in)	Ntz,
		real, dimension(1:ntz, -1:2)	dBB,
		real, dimension(1:ntz)	XX,
		real, dimension(1:ntz)	YY,
		real, dimension(1:ntz)	length,
		real	DDl,
		real	MMl,
		integer, intent(in)	iflag
	)

Computes \(B^2\), and the spectral condensation constraints if required, on the interfaces, \({\cal I}_i\).

field strength

• The field strength is given by \(B^2 = B^s B_s + B^\theta B_\theta + B^\zeta B_\zeta\), and on the interfaces \(B^s=0\) by construction.
• The magnetic field is \(\sqrt g \; {\bf B} = (\partial_\theta A_\zeta - \partial_\zeta A_\theta ) {\bf e}_s - \partial_s A_\zeta {\bf e}_\theta + \partial_s A_\theta {\bf e}_\zeta\).
• The covariant components of the field are computed via \(B_\theta = B^\theta g_{\theta\theta} + B^\zeta g_{\theta\zeta}\) and \(B_\zeta = B^\theta g_{\theta\zeta} + B^\zeta g_{\zeta\zeta}\).
• The expression for \(B^2\) is

  \begin{eqnarray} (\sqrt g)^2 B^2 = A_\zeta^\prime \; A_\zeta^\prime \; g_{\theta\theta} - 2 \; A_\zeta^\prime \; A_\theta^\prime \; g_{\theta\zeta} + A_\theta^\prime \; A_\theta^\prime \; g_{\zeta\zeta}, \end{eqnarray}

  where the \("\prime"\) denotes derivative with respect to \(s\).
• The quantity returned is

  \begin{eqnarray} F \equiv \texttt{pscale} \times \frac{P}{V^\gamma} + \frac{B^2}{2}, \end{eqnarray}

  where \(P\equiv\) adiabatic and \(V\equiv\) volume.

spectral constraints

• In addition to the physical-force-balance constraints, namely that \([[p+B^2/2]]=0\) across the interfaces, additional angle constraints are required to obtain a unique Fourier representation of the interface geometry.
• Introducing the angle functional: a weighted combination of the "polar" constraint; the normalized, poloidal, spectral width (Hirshman & Meier (1985) [3], Hirshman & Breslau (1998) [2]) the poloidal-angle origin constraint; and the "length" of the angle curves

  \begin{eqnarray} F \!\equiv\! \sum_{i=1}^{N-1} \alpha_i \underbrace{\oint\!\!\!\oint \!d\theta d\zeta \,\,\, \!\! \frac{1}{\Theta_{i,\theta}}}_{polar-angle} + \sum_{i=1}^{N-1} \beta_i \, \underbrace{\;\; M_i \;\;}_{spectral-width} + \sum_{i=1}^{N-1} \gamma_i \int_{0}^{2\pi} \! \frac{1}{2}\left[Z_i(0,\zeta)-Z_{i,0}\right]^2 d\zeta + {\color{red}{\oint\!\!\!\oint \!d\theta d\zeta \,\,\, \!\! \sum_{i=1}^{N} \delta_i \, \underbrace{\;\;\;\; L_i \;\;\;\;}_{poloidal-length}}} \end{eqnarray}

  where \(i\) labels the interfaces, and

  \begin{eqnarray} \Theta_{i,\theta} & \equiv & \frac{x \, y_\theta - x_\theta \, y}{x^2+y^2}, \\ M_i & \equiv & \frac{\sum_j m_j^p ( R_{j,i}^2 + Z_{j,i}^2 )}{\sum_j ( R_{j,i}^2 + Z_{j,i}^2 )} ,\\ {\color{red}{L_i}} & \equiv & {\color{red}{\sqrt{ [R_{i}(\theta,\zeta)-R_{i-1}(\theta,\zeta)]^2+[Z_{i}(\theta,\zeta)-Z_{i-1}(\theta,\zeta)]^2}}}, \end{eqnarray}

  and where \(j\) labels the Fourier harmonics. The \(\alpha_i\), \(\beta_i\), \(\gamma_i\) and \(\delta_i \equiv \) sweight are user-supplied weight factors.
• The polar constraint is derived from defining \(\tan \Theta \equiv y/x\), where

  \begin{eqnarray} x(\theta,\zeta) & \equiv & R_{i}(\theta,\zeta)-R_{i,0}(\zeta), \\ y(\theta,\zeta) & \equiv & Z_{i}(\theta,\zeta)-Z_{i,0}(\zeta), \end{eqnarray}

  and where the geometric center of each interface is given by the arc-length weighted integrals, see rzaxis(),

  \begin{eqnarray} R_{i,0} & \equiv & \int_0^{2\pi} \!\!\!\! d\theta \; R_i(\theta,\zeta) \sqrt{R_{i,\theta}(\theta,\zeta)^2+Z_{i,\theta}(\theta,\zeta)^2}, \\ Z_{i,0} & \equiv & \int_0^{2\pi} \!\!\!\! d\theta \; Z_i(\theta,\zeta) \sqrt{R_{i,\theta}(\theta,\zeta)^2+Z_{i,\theta}(\theta,\zeta)^2}, \end{eqnarray}

  and \(\cos \Theta = x / \sqrt{x^2+y^2}\) has been used to simplify the expressions and to avoid divide-by-zero.
• Only "poloidal tangential" variations will be allowed to find the extremum of \(F\), which are described by

  \begin{eqnarray} \delta R_i(\theta,\zeta) & \equiv & R_{i,\theta}(\theta,\zeta) \, \delta u_i(\theta,\zeta),\\ \delta Z_i(\theta,\zeta) & \equiv & Z_{i,\theta}(\theta,\zeta) \, \delta u_i(\theta,\zeta), \end{eqnarray}

  from which it follows that the variation in each Fourier harmonic is

  \begin{eqnarray} \delta R_{j,i} = \displaystyle \oint\!\!\!\oint \!d\theta d\zeta \,\,\, R_{i,\theta}(\theta,\zeta) \, \delta u_i(\theta,\zeta) \, \cos(m_j\theta-n_j\zeta), \\ \delta Z_{j,i} = \displaystyle \oint\!\!\!\oint \!d\theta d\zeta \,\,\, Z_{i,\theta}(\theta,\zeta) \, \delta u_i(\theta,\zeta) \, \sin(m_j\theta-n_j\zeta), \end{eqnarray}

  and

  \begin{eqnarray} \delta R_{i,\theta}(\theta,\zeta) & \equiv & R_{i,\theta\theta}(\theta,\zeta) \, \delta u_i(\theta,\zeta) + R_{i,\theta}(\theta,\zeta) \, \delta u_{i,\theta} (\theta,\zeta) \\ \delta Z_{i,\theta}(\theta,\zeta) & \equiv & Z_{i,\theta\theta}(\theta,\zeta) \, \delta u_i(\theta,\zeta) + Z_{i,\theta}(\theta,\zeta) \, \delta u_{i,\theta} (\theta,\zeta) \end{eqnarray}

• The variation in \(F\) is

  \begin{eqnarray} \delta F & = & \sum_{i=1}^{N-1} \alpha_i \;\;\;\oint\!\!\!\oint \!d\theta d\zeta \,\,\, \!\! \left( \frac{-2\Theta_{i,\theta\theta}}{\Theta_{i,\theta}^2} \right) \delta u_i \nonumber \\ & + & \sum_{i=1}^{N-1} \beta_i \;\;\;\oint\!\!\!\oint \!d\theta d\zeta \,\,\, \!\! \left( R_{i,\theta} X_i + Z_{i,\theta} Y_i \right) \delta u_i \nonumber \\ & + & \sum_{i=1}^{N-1} \gamma_i \;\;\; \int \! d\zeta \; \left( Z_{i}(0,\zeta)-Z_{i,0} \right) Z_{i,\theta} \; \delta u_i \nonumber \\ & + & {\color{red}{\sum_{i=1}^{N-1} \delta_i \;\;\;\oint\!\!\!\oint \!d\theta d\zeta \,\,\, \!\! \left( \frac{\Delta R_{i } R_{i,\theta} + \Delta Z_{i } Z_{i,\theta}}{L_{i }} \right) \delta u_i}} \nonumber \\ & - & {\color{red}{\sum_{i=1}^{N-1} \delta_{i+1} \oint\!\!\!\oint \!d\theta d\zeta \,\,\, \!\! \left( \frac{\Delta R_{i+1} R_{i,\theta} + \Delta Z_{i+1} Z_{i,\theta}}{L_{i+1}} \right) \delta u_i}} \label{eq:firstvariation_lforce} \end{eqnarray}

  where, for the stellarator symmetric case,

  \begin{eqnarray} X_i & \equiv & \displaystyle \sum\nolimits_j ( m_j^p - M_i ) \, R_{j,i} \cos(m_j\theta-n_j\zeta), \\ Y_i & \equiv & \displaystyle \sum\nolimits_j ( m_j^p - M_i ) \, Z_{j,i} \sin(m_j\theta-n_j\zeta), \end{eqnarray}

  and

  \begin{eqnarray} {\color{red}{\Delta R_{i}}} & \equiv & {\color{red}{R_{i}(\theta,\zeta)-R_{i-1}(\theta,\zeta)}},\\ {\color{red}{\Delta Z_{i}}} & \equiv & {\color{red}{Z_{i}(\theta,\zeta)-Z_{i-1}(\theta,\zeta)}}, \end{eqnarray}

• The spectral constraints derived from Eqn. \((\ref{eq:firstvariation_lforce})\) are

  \begin{eqnarray} I_i(\theta,\zeta) & \equiv & - 2 \alpha_i \frac{\Theta_{i,\theta\theta}}{\Theta_{i,\theta}^{2}} + \beta_i \left( R_{i,\theta} X_i + Z_{i,\theta} Y_i \right) + \gamma_i \left( Z_{i}(0,\zeta) - Z_{i,0} \right) Z_{i,\theta}(0,\zeta) \nonumber \\ & + & {\color{red}{\delta_{i } \frac{\Delta R_{i } R_{i,\theta} + \Delta Z_{i } Z_{i,\theta}}{L_i}}} - {\color{red}{\delta_{i+1} \frac{\Delta R_{i+1} R_{i,\theta} + \Delta Z_{i+1} Z_{i,\theta}}{L_{i+1} } } } \label{eq:spectralconstraints_lforce} \end{eqnarray}

• Note that choosing \(p=2\) gives \(X=-R_{\theta\theta}\) and \(Y=-Z_{\theta\theta}\), and the spectrally condensed angle constraint, \(R_\theta X + Z_\theta Y=0\), becomes \(\partial_\theta (R_\theta^2+Z_\theta^2)=0\), which defines the equal arc length angle.
• The poloidal-angle origin term, namely \(\gamma_i \left( Z_{i}(0,\zeta) - Z_{i,0} \right) Z_{i,\theta}(0,\zeta)\) is only used to constrain the \(m_j=0\) harmonics.
• The construction of the angle functional was influenced by the following considerations:
  • The minimal spectral width constraint is very desirable as it reduces the required Fourier resolution, but it does not constrain the \(m=0\) harmonics and the minimizing spectral-width poloidal-angle may not be consistent with the poloidal angle used on adjacent interfaces.
  • The regularization of the vector potential and the coordinate interpolation near the coordinate origin (see elsewhere) assumes that the poloidal angle is the polar angle.
  • The user will provide the Fourier harmonics of the boundary, and thus the user will implicitly define the poloidal angle used on the boundary.
  • Minimizing the length term will ensure that the poloidal angle used on each interface is smoothly connected to the poloidal angle used on adjacent interfaces.
• A suitable choice of the weight factors, \(\alpha_i\), \(\beta_i\), \(\gamma_i\) and \(\delta_i\), will ensure that the polar constraint dominates for the innermost surfaces and that this constraint rapidly becomes insignificant away from the origin; that the minimal spectral constraint dominates in the "middle"; and that the minimizing length constraint will be significant near the origin and dominant near the edge, so that the minimizing spectral width angle will be continuously connected to the polar angle on the innermost surfaces and the user-implied angle at the plasma boundary. The length constraint should not be insignificant where the spectral constraint is dominant (so that the \(m=0\) harmonics are constrained).
• The polar constraint does not need normalization. The spectral width constraint has already been normalized. The length constraint is not yet normalized, but perhaps it should be.
• The spectral constraints given in Eqn. \((\ref{eq:spectralconstraints_lforce})\) need to be differentiated with respect to the interface Fourier harmonics, \(R_{j,i}\) and \(Z_{j,i}\). The first and second terms lead to a block diagonal hessian, and the length term leads to a block tri-diagonal hessian.
• Including the poloidal-angle origin constraint means that the polar angle constraint can probably be ignored, i.e. \(\alpha_i=0\).

Parameters

[in]	lvol
[in]	iocons
[in]	ideriv
[in]	Ntz
	dBB
	XX
	YY
	length
	DDl
	MMl
[in]	iflag

References inputlist::adiabatic, allglobal::ate, allglobal::ato, allglobal::aze, allglobal::azo, allglobal::bemn, allglobal::bomn, allglobal::cfmn, allglobal::comn, coords(), allglobal::cpus, allglobal::drij, allglobal::dzij, allglobal::efmn, allglobal::evmn, inputlist::gamma, allglobal::guvij, constants::half, allglobal::iemn, inputlist::igeometry, allglobal::ijimag, allglobal::ijreal, allglobal::im, allglobal::in, allglobal::iomn, allglobal::irbc, allglobal::irbs, allglobal::irij, allglobal::izbc, allglobal::izbs, allglobal::izij, allglobal::jiimag, allglobal::jireal, inputlist::lcheck, allglobal::lcoordinatesingularity, lforce(), inputlist::lrad, allglobal::mmpp, allglobal::mn, allglobal::mpi_comm_spec, allglobal::myid, allglobal::ncpu, allglobal::notstellsym, allglobal::nt, inputlist::nvol, allglobal::nz, allglobal::odmn, allglobal::ofmn, constants::one, fileunits::ounit, allglobal::pemn, allglobal::pomn, inputlist::pscale, allglobal::regumm, allglobal::rtt, allglobal::semn, allglobal::sfmn, allglobal::sg, allglobal::simn, allglobal::somn, allglobal::trij, allglobal::tt, constants::two, allglobal::tzij, allglobal::vvolume, allglobal::yesstellsym, and constants::zero.

Referenced by dfp200(), evaluate_dbb(), and lforce().

Here is the call graph for this function:

Here is the caller graph for this function: