Introduction

In this numerical experiment we will convert a two-dimensional free-current density, $\vec{J}_f$, into a current vector potential, $T$. The following figure illustrates the idea:

The two-dimensional curl-curl equation for the magnetic vector potential makes sense only in the planar configuration in which the magnetic vector potential, $\vec{A}$, is an in-plane vector. Consequently, only one flavor of two-dimensional current vector potential makes any sense: the out-of-plane vector $ T $, i.e., a scalar, on a planar problem domain. In all other configurations (if they do exist), the curl-curl equation is replaced with an equivalent div-grad equation (see Section 3.4. Magnetic vector potential in two dimensions). The div-grad partial differential equation has no need in the current vector potential.

The compatibility condition for the two-dimensional curl-curl equation,

\[ \vec{\nabla}\overset{V}{\times}\bigg(\frac{1}{\mu}\vec{\nabla}\overset{S}{\times}\vec{A}\bigg) = \vec{J}_f, \]

requires the free-current density, $\vec{J}_f$, to be derived as a curl of a current vector potential,

\[ \vec{J}_f = \vec{\nabla}\overset{V}{\times}T. \]

(see Section 5.5. Compatibility conditions). Then the initial curl-curl equation becomes

\[ \vec{\nabla}\overset{V}{\times}\bigg(\frac{1}{\mu}\vec{\nabla}\overset{S}{\times}\vec{A}\bigg) = \vec{\nabla}\overset{V}{\times} T. \]

In order to solve the last equation, one needs to convert the free-current density, $\vec{J}_f$, given by the formulation of the problem into a current vector potential, $T$. The most straightforward way to do so is to solve the following partial differential equation:

\[ \vec{\nabla}\overset{S}{\times}\bigg(\vec{\nabla}\overset{V}{\times} T\bigg) = \vec{\nabla}\overset{S}{\times}\vec{J}_f. \]

The last equation can be replaced with the following equivalent equation:

\[ -\vec{\nabla}\cdot\bigg( \vec{\nabla} T\bigg) = \vec{\nabla}\overset{S}{\times}\vec{J}_f. \]

This equation can be solved with a help of StaticScalarSolver::Solver class template. This numerical experiment illustrates how to do this.

First we will solve the last partial differential equation for $ T $. Second, we will convert the numerically calculated $T$ back into the free-current density $\vec{J}_{f}$. The Bossavit's diagrams below illustrate these two calculations.

Third, we will compare the numerically calculated $T$ and $\vec{J}_{f}$ to the corresponding closed-form analytical expressions by computing $L^2$ error norms. As per usual, we will repeat the calculations on a family of progressively refined meshes and observe the behavior of the $L^2$ error norms.

In the current configuration the current vector potential $T$ is implicitly gauged by the translation symmetry of the problem. For this reason we can write down an closed-form analytical expression for $T$ and compare it to its numerically computed counterpart.

Note that the current vector potential computed in this numerical experiment can only be used in two-dimensional problems. How to solve for a three-dimensional current vector potential is illustrated in the Current vector potential (cvp-i/) numerical experiment.

Implementation

The problem domain and the exact solution

The problem domain is assumed to be a disk concentric with the origin. The radius of the disk equals $ d_2 $. The boundary of the disk represents infinity. We apply the homogeneous Dirichlet boundary condition,

\[ T = 0. \]

The surface free-current density, $ \vec{J}_f $, is assumed to be confined to a circular ring concentric with the origin. The inner and outer radii of the ring equal $a$ and $b$, respectively. The figure below illustrates this configuration.

The free-current density is assumed to have the following form:

\[ \vec{J}_f = \left\{ \begin{aligned} & 0 &&\text{if }&& r < a\\ & -y\hat{i} + x\hat{j} && \text{if } && a < r < b\\ & 0 &&\text{if }&& r > b, \end{aligned} \right. \]

where

\[ r = \sqrt{x^2+y^2}. \]

In this configuration the current vector potential is implicitly gauged by the Coulomb gauge. For this reason we can write down a closed form analytical expression for $T$,

\[ T = \left\{ \begin{aligned} & \frac{1}{2} (b^2-a^2) &&\text{if }&& r \le a \\ & - \frac{1}{2} (x^2+y^2-b^2) && \text{if } && a \le r \le b \\ & 0 &&\text{if }&& r \ge b. \end{aligned} \right. \]

The boundary value problem

The following boundary value problem is solved in this numerical experiment:

\begin{equation} \begin{array}{rrcll} -\vec{\nabla}\cdot\bigg( \vec{\nabla} T\bigg) = \vec{\nabla}\overset{S}{\times}\vec{J}_f & \text{in} & \Omega & \text{(i)},\\ T = 0 & \text{on} & \Gamma_{D1} & \text{(ii)}. \\ \end{array} \end{equation}

To verify the solution, the numerically computed current vector potential is converted back to the free-current density as

\[ \vec{J}_{f} = \vec{\nabla}\overset{V}{\times}T. \]

The program

This experiment is implemented in accordance with the base code structure. The build process generates one executable file: cvp-ii-linear. To rebuild it, change into cvp-ii/build/Release directory and execute the following:

user@computer .../cvp-ii/build/Release$ ./clean

user@computer .../cvp-ii/build/Release$ ./build

Then the experiment needs to be executed again. To do so, change into the cvp-ii/bin/Release directory and execute the executable file,

user@computer .../cvp-ii/build/Release$ cd ../../bin/Release

user@computer .../cvp-ii/bin/Release$ ./cvp-ii-linear

This will generate vtu files in the cvp-ii/bin/Release/Data/circle-linear/ directory. They can be viewed with a help of a fresh version of ParaView. The VisIt visualization software will not do as this numerical experiment uses second-degree mapping.

Note that executable files require a set of meshes to be present in the cvp-ii/gmsh/data directory. If they are missing, they can be generated anew. This can be done by changing into cvp-ii/gmsh directory and executing the following:

user@computer .../cvp-ii/gmsh$ ./clean

user@computer .../cvp-ii/gmsh$ ./build

This will generate a set of globally refined meshes in cvp-ii/gmsh/data.

The SettingsCVPII class allows switching on three useful features: printing time tables on the computer screen, logging convergence data of the conjugate gradient solver, and saving the exact solution into the vtu files next to the numerical solution.

It worth mentioning that the free-current density, $\vec{J}_f$, equals zero on the boundary by the definition of the problem. For this reason we can neglect the boundary integral $I_{b3-2}$ in the Recipe for static scalar solver in 2D (current vect. potential). Consequently, we pass the input parameter type_of_pde_pde = 2 to the constructor of the StaticScalarSolver::Solver:

...
 
SolverCVPII(unsigned int p,
            unsigned int mapping_degree,
            unsigned int r,
            std::string fname)
  : Solver<2>(p,
              mapping_degree,
              2,  // That is, type_of_pde_rhs = 2.
              fname,
              &exact_solution,
              false,
              true,
              SettingsCVPII::print_time_tables,
              SettingsCVPII::project_exact_solution,
              true)
  , fname(fname)
  , r(r)
{
  StaticScalarSolver::Solver<2>::run();
}
 
  ...

Simulation results

The experiment was executed with the following settings: $a = 0.3[m]$, $b = 0.5[m]$, and $d_2=1.0[m]$.

Current vector potential

The figure below illustrates the calculated current vector potential, $T$.

The corresponding convergence table is presented below.

p	r	cells	dofs	$\\|e\\|_{L^2}$	$\alpha_{L^2}$	$\\|e\\|_{H^1}$	$\alpha_{H^1}$
1	12	3388	3433	1.94e-05	-	3.72e-03	-
1	13	4032	4081	1.63e-05	2.00	3.41e-03	1.00
1	14	4732	4785	1.39e-05	2.00	3.15e-03	1.00
1	15	5488	5545	1.20e-05	2.00	2.92e-03	1.00
2	12	3388	13641	6.21e-09	-	1.01e-06	-
2	13	4032	16225	4.38e-09	4.00	7.76e-07	3.00
2	14	4732	19033	3.18e-09	4.00	6.10e-07	3.00
2	15	5488	22065	2.37e-09	4.00	4.88e-07	3.00
3	12	3388	30625	6.21e-09	-	1.01e-06	-
3	13	4032	36433	4.39e-09	4.00	7.76e-07	3.00
3	14	4732	42745	3.18e-09	4.00	6.10e-07	3.00
3	15	5488	49561	2.37e-09	4.00	4.88e-07	3.00

The following notations were used in the table header:

p - the degree of the interpolating Lagrange polynomials that constitute the shape functions.
r - the number of nodes on transfinite lines, see discussion in here.
cells - the total amount of active cells.
dofs - the amount of degrees of freedom in the active cells.
$\|e\|_{L^2}$ - the $L^2$ error norm.
$\alpha_{L^2}$ - the order of convergence of the $L^2$ error norm.

Free-current density

The figure below illustrates the calculated free-current density, $\vec{J}_{f}$, and its magnitude, $|\vec{J}_{f}|$.