Test Case : Feelpp CFPDEs - Magnetostatic HTS with erf

1. Introduction

This is the test case of High-Temperature Superconductors with the A-V Formulation and Gauge Condition on a bulk cylinder geometry surrounded by air in 2D coordinates.

2. Run the Calculation

The command line to run this case is :

    mpirun -np 16 feelpp_toolbox_coefficientformpdes --config-file=mqs_erf.cfg --cfpdes.solver=Picard --cfpdes.verbose_solvertimer=1

This simulation uses the last version 109 of Feelpp.

3. Data Files

The case data files are available in Github here :

4. Equation

A-V Formulation in 2D coordinates
\[- \begin{matrix} \Delta A_z = -\mu J & \text{on } \Omega & \text{(AV 2D)} \\ \end{matrix}\]

With :

  • \(A_{z}\) : \(z\) component of potential magnetic field

  • \(J=J_c \text{erf}\left(\frac{-A_z}{A_r}\right)\) : current density \(A/m^2\)

  • \(\mu\) : electric permeability \(kg/A^2/S^2\)

5. Geometry

The geometry is a circle in 2D coordinates \((x,y)\) representing a bulk cylinder, surrounded by air.

circle
Geometry

The geometrical domains are :

  • Conductor : the cylinder

  • Air : the air surrounding Conductor

    • Infty : the Air 's boundary

Symbol

Description

value

unit

\(R\)

radius of cylinder

\(0.001\)

m

\(R_{inf}\)

radius of infty border

\(0.01\)

m

6. Parameters

The parameters of the problem are :

  • On Conductor :

Symbol

Description

Value

Unit

\(\mu=\mu_0\)

magnetic permeability of vacuum

\(4\pi.10^{-7}\)

\(kg \, m / A^2 / S^2\)

\(A_r\)

parameter resulting from the combination of \(E_0\) and the time it takes to reach the peak of AC excitation

\(1.10^{-7}\)

\(Wb/m\)

\(b_{ext}\)

external applied field

\(0.02\)

\(T\)

\(J_c\)

critical current density

\(1.10^8\)

\(A/m^2\)

The error function erf is defined by :

\[erf(x)=\frac{2}{\sqrt{\pi}}\int_0^x \exp(-t^2)dt\]
Erf plot

This function is not implemented on feelpp, so we use a fit on a .csv with a large panel of values calculated with the function.

  • On Air :

Symbol

Description

Value

Unit

\(\mu=\mu_0\)

magnetic permeability of vacuum

\(4\pi.10^{-7}\)

\(kg \, m / A^2 / S^2\)

On JSON file, the parameters are written :

Parameters on JSON file
    "Parameters":
    {
        "mu":"4*pi*1e-7",
        "A_r":1e-7,
	    "Bext":0.02,
        "jc":1e8,
        "erf": {
            "type":"fit",
            "filename":"$cfgdir/erf.csv",
            "abscissa":"X",
            "ordinate":"Erf",
            "interpolation":"P1",
            "expr":"-magnetic_A/A_r:magnetic_A:A_r"
        }
    },

7. Boundary Conditions

For the Dirichlet boundary conditions, we want to impose the applied magnetic field :

We have \(B_{ext}=0.02 T\), and in 2D Cartesian coordinates, \(B=\nabla\times A\) becomes \(B_x=\partial A/\partial y\) and \(B_y=-\partial A/\partial x\). So in order to obtain a magnetic field \(B_{ext}\) along \(y\), it is sufficient to impose \(A_z=-xB_{ext}\).

Finally we have :

  • On Infty : \(A = -x B_{ext}\)

On JSON file, the boundary conditions are written :

Boundary conditions on JSON file
    "BoundaryConditions":
    {
        "magnetic":
        {
            "Dirichlet":
            {
                "magdir":
                {
                    "markers":["Infty"],
                    "expr":"-x*Bext:x:Bext"
                }
            }
        }
    },

8. Weak Formulation

Weak formulation of A Formulation in two dimensions approximation
\[ \int_{\Omega}{ - \frac{1}{\mu} \nabla A_z \cdot \nabla \phi \ dxdy } = \int_{\Omega_c}{ J_c erf\left(\frac{-A_z}{A_r}\right) \cdot \phi \ dxdy } + \int_{\Gamma_D}{ -B_0 x } \hspace{1cm} \text{(Weak A 2D)}\]

9. Coefficient Form PDEs

The Feelpp toolboxe Coefficient Form PDEs is used here. The coefficients associated to the Weak Formulation are :

  • On Conductor :

Coefficient

Description

Expression

\(c\)

diffusion coefficient

\(1\)

\(f\)

source term

\(\mu J_c \text{erf}\left(\frac{-A}{A_r}\right)\)

The bilinear form was formulated as a non-linear problem, which in CFPDES requires the source term to be multiplied by the unknown A. Hence, for the sake of consistency with the model, the source term is written as a reaction coefficient and multiplied by the term \((1/A)\). If \(A=0\), the source term is multiplied by 1, that’s why a notzero parameter is introduced.

Coefficient

Description

Expression

\(c\)

diffusion coefficient

\(1\)

\(a\)

absorption or reaction coefficient

\(-\frac{1}{A}\mu J_c \text{erf}\left(\frac{-A}{A_r}\right)\)

  • On Air :

Coefficient

Description

Expression

\(c\)

diffusion coefficient

\(1\)

On JSON file, the coefficients are written :

CFPDEs coefficients on JSON file
    "Materials":
    {
        "Conductor":
        {
            "markers":["Conductor"],


            "notzero":"(1/magnetic_A)^(1-(magnetic_A>-1E-100)*(magnetic_A<1E-100)):magnetic_A",

	        "magnetic_c":"1",
            "magnetic_a":"-mu*jc*erf*notzero:mu:jc:erf:notzero"
        },
        "Air":
        {
	        "markers":["Air"],
            "magnetic_c":"1"
        }
    },

10. Mesh & Solver

  • Mesh :

circlemesh
Mesh of Geometry
  • Non-linear Solver : Picard

    • Maximum number of iterations : \(600\)

11. Results

The results that we obtain with this formulation with Feelpp are compared to the results of the article A numerical model to introduce student to AC loss calculation in superconductors where the finite element software FreeFEM is used.

11.1. Electric current density

The electric current density \(J\) is defined by :

\[ J= J_c \text{erf}\left(\frac{-A}{A_r}\right)\]
2D erf j
Figure 1. "Electric current density \(J (A/m^2)\)

We compare the current density profiles with Feelpp and FreeFEM on the \(O_r\) axis, on the diameter of the cylinder, for a maximum applied field of 0.02 T.

L2 Relative Error Norm : \(3.4 \%\)

11.2. Magnetic flux density

The magnetic flux density \(B\) is defined by:

\[ B=\nabla\times A =\begin{pmatrix}\partial_y A\\ -\partial_x A\\ 0\end{pmatrix}\]

Therefore, \(B_y\), the y-component of the magnetic flux density is defined as \(-\partial_x A\) :

2D erf bz
Figure 2. "y-component of the magnetic flux density \(B_y (T)\)

12. References

  • A numerical model to introduce students to AC loss calculation in superconductors. Francesco Grilli and Enrico Rizzo 2020 Eur. J. Phys. 41 045203