Magnetic PDE¶
Maxwell's Equations¶
The set of Maxwell's equations is based on the experiments of Ampére, Gauss and Faraday. The novelty of Maxwell was to introduce the displacement current density in Ampére's law, which enabled him to describe the propagation of electromagnetic waves. The full set of this equations can be written in local form as \begin{equation} \label{eq:Ampere} \textbf{curl}( \mathbf{H}) = \mathbf{J} + \frac{\partial \mathbf{D}}{\partial t}, \end{equation}
For the global form, we consider the arbitrary volume, depicted in the figure below and apply Gauss- and Stokes-theorem to (\ref{eq:Ampere})-(\ref{eq:solenoidalB}).
Volume domain for the global form of Maxwell's equations
This leads to \begin{equation} \label{eq:Ampereglobal} \int_{\partial \Gamma} \mathbf{H}\cdot \mathbf \tau_{\Gamma} = \int_{\Gamma} \mathbf{J} \cdot \mathbf n_{\Gamma} + \int_{\Gamma} \frac{\partial \mathbf{D}}{\partial t} \cdot \mathbf n_{\Gamma} , \end{equation}
\begin{equation} \int_{\partial \Omega} \mathbf{B} \cdot \mathbf n_{\Omega} = 0, \label{eq:solenoidalBglobal} \end{equation} which can be identified as the law of Ampére with Maxwell's contribution of the displacement current density (\ref{eq:Ampereglobal}), the law of Faraday (\ref{eq:Faradayglobal}), the law of Gauss (\ref{eq:Gaussglobal}) and the global solenoidal \mathbf B-field (\ref{eq:solenoidalBglobal}).
The equations above are additionally equipped with three constitutive laws, where isotropy of the electric conducticity \gamma, electric permittivity \epsilon and magnetic permeability \mu, is assumed \begin{equation} \mathbf{J} = \gamma (\mathbf{E} + \mathbf{v} \times \mathbf{B}), \label{eq:J} \end{equation} \begin{equation} \mathbf{D} = \epsilon \mathbf{E}, \end{equation} \begin{equation} \mathbf{B} = \mu \mathbf{H}. \label{eq:BH} \end{equation}
In the general case, the magnetic permeability \mu respectively reluctivity \nu=\frac{1}{\mu} is not constant for all values of \mathbf{B}, as it depends on the flux density \nu=\nu(|\mathbf{B}|).
The electric field intesity can then be split into a purely solenoidal part \mathbf E_s and an impressed one \mathbf E_{\textrm i}, which leads to the following splitting of the total current density (\ref{eq:J})
in which \gamma \mathbf{E}_{\textrm i} can be considered as an impressed (prescribed) current density \mathbf J_{\textrm i}.
An important fact, which follows from the law of Ampére and (\ref{eq:J}), is that the total current is divergence-free \textbf{div} (\mathbf J) = 0.
Currently, only the low frequency regime of the electromagnetic spectrum is implemented in openCFS, which means that the time derivative of the displacement current density is neglected and therefore wave propagation is not part of the consideration.
PDE Types¶
Depending on if the weak form is discretized with nodal (classical Lagrange) elements or with Nédélec edge elements, we have different PDE's implemented:
-
Discretization with nodal elements:
-
Discretization with Nédélec edge elements:
However for most purposes, the edge element version is way more versatile but only 3D problems can be considered. If a 2D example has to be considered, an extraction of the 2D surface into the depth direction with one element is recommended, resulting in a quasi-2D problem.
Boundary conditions and loads¶
- Essential (Dirichlet) boundary conditions on \Gamma_{\textrm D}, where \mathbf A \times \mathbf n and \mathbf A' \times \mathbf n is set to a zero. If these values are set to zero, we are in accordance with the considered function space V because it has zero tangential trace on the boundary. The fact that \mathbf A \times \mathbf n = 0 and $\mathbf A' \times \mathbf n = $ also results in
Boundary conditions¶
In openCFS this type of BC can be specified by
- Natural (Neumann) boundary conditions on \Gamma_{\textrm N}, where \nu \textbf{curl} \mathbf A \times \mathbf n is set to a certain value, corresponding to a surface current density.
not implemented yet
- Surface impedance boundary conditions (SIBC), where a linear relation between \mathbf E \times \mathbf n and \mathbf H \times \mathbf n is prescribed at the boundary, which can be used for surface excitation. This method presents an elegant way to drastically decrease the number of degrees of freedom in coil and inductor regions but it is only valid if the skin depth tends to zero, respectively if the frequency tends towards infinity.
partially implemented
Loads¶
In openCFS you can define these loads (list not completed):
- Magnetic excitation
Coils¶
Coils are used in openCFS to model current-carrying conductors which produce a magnetic field. Here is a simple example how to use it as a xml-input, be aware that they must be specified after <bcsAndLoads>
in <coilList>
.
It is recommended to prescibe the current density as current value (and use a wireCrossSection of 1). E.g. for a coil with a coilarea of "A", "n" windings and a wire current of "I" specify the current as shown:
<coilList>
<coil id = " coil1 "> <! -- reference id for coil results -- >
<source type = " current " value = "I*n/A" />
<part id = " 1 "> <! -- may have several parts -- >
<regionList>
<region name = "coilwire" />
</regionList>
<direction>
<analytic>
<comp dof = "z" value = " 1.0 " />
</analytic>
</direction>
<wireCrossSection area = "1" />
</part>
</coil>
</coilList>
Potentials and Formulations¶
According to the solenoidal (divergence-free/no magnetic monopoles) \mathbf B-field (\ref{eq:solenoidalB}), there exists a vector potential \mathbf A, called the magnetic vector potential, defined by
Plugging that into (\ref{eq:Faraday}), we obtain
Since the kernel of the curl operator consists of all gradients of scalar functions, the gradient of a potential, called the electric scalar potential V can be added to the electric field intensity without changing (\ref{eq:curlEpA})
The PDE-formulation containing both \mathbf A and V can be accessed via
Building on that, there are several ways how to fix the gauges, based on additional (constraint) equations. There is a large manifold of possible gauges, whereas only the most important ones for computational electromagnetics are mentioned here.
Coulomb gauge¶
This gauge states that
The huge benefit of this gauge is, that if we are only interested in the magnetic flux denstity \mathbf B, there is no need to additionally gauge the formulation because of (\ref{eq:bequalcurlA}). If however the magnetic vector potential itself is observed, gauging techniques must be considered. For the magnetostatic case, we add a small ficticous conductivity (relaxation parameter), which does not effect the solution dramatically but we can ensure uniqueness of \mathbf A. If there is no definition of the relaxation parameter, the default value is 1\cdot 10^{-6}.
Temporal-/Weyl-gauge -- Classical A Formulation¶
The more important gauge, regarding computational electromagnetics is the temporal- or Weyl-gauge, e.g. see p.48 in 1, which sets the electric scalar potential to zero V=0. Doing this, results in the classical \mathbf A formulation, defined by the following PDE-formulation
Using the Weyl-gauge, results in the following description of the electric field intensity
Combining (\ref{eq:intermediate1}), (\ref{eq:Ampere}), (\ref{eq:Faraday}) and the constitutive law (\ref{eq:J}) results in the well known curl-curl equation
Weak Form¶
The considered function space for the classical A formulation is the curl-conforming space, defined by
By multiplying the strong form (\ref{eq:curl1}) with all test functions \mathbf A'\in V and integrating over the domain \Omega, we obtain
and using the integration by parts for the first term, results in the final weak form
where the first term on the right hand side defines the boundary conditions. By identifying the expression in the brackets as
different boundary conditions can be stated.
Material¶
These material parameters are specified in the material-xml file, e.g. for copper
<material name="copper">
<magnetic>
<electricConductivity>
<linear>
<isotropic>5e7</isotropic>
</linear>
</electricConductivity>
<magneticPermeability>
<linear>
<isotropic> 1.256637061E-06 </isotropic>
</linear>
</magneticPermeability>
</magnetic>
</material>
An example of a field dependent magnetic permeanbility is defined as follows
<material name="iron">
<magnetic>
<electricConductivity>
<linear>
<isotropic>
<real>1.0E5</real>
</isotropic>
</linear>
</electricConductivity>
<permeability>
<linear>
<isotropic>
<real> 9.17290e-4 </real>
</isotropic>
</linear>
<nonlinear>
<isotropic>
<dependency></dependency>
<approxType>smoothSplines</approxType>
<measAccuracy> 0.05 </measAccuracy>
<maxApproxVal> 3.0 </maxApproxVal>
<dataName>bh1.fnc</dataName>
</isotropic>
</nonlinear>
</permeability>
</magnetic>
</material>
Analysis Types¶
Depending on the temporal setting, the formulation (\ref{eq:curl1}) can be used for
-
Static: \partial / \partial t (\cdot)= 0
-
Transient: \partial / \partial t (\cdot) \neq 0
- Harmonic: \partial / \partial t (\cdot) = j \omega (\cdot)
Postprocessing Results¶
Based on the solution A, different postprocessing results can be defined in openCFS, described in the following.
As in every other PDE in openCFS, the results are diveded into element-/ node-/ surface-element-/ surface-region and region-results, whereas we have to distinguish between the nodal- and edge-version.
Results for magEdgePDE (edge element version)¶
Element Results¶
-
Magnetic Vector Potential A (primary solution)
-
Time Derivative of Magnetic Vector Potential
- transient case: \begin{equation} \frac{\partial \mathbf A}{\partial t} \end{equation}
- harmonic case: \begin{equation} j \omega \mathbf A \end{equation}
-
Magnetic Flux Density B \begin{equation} \mathbf B = \textbf{curl} (\mathbf A) \end{equation}
-
Magnetic Field Strength H \begin{equation} \mathbf H = \nu \mathbf B = \nu \textbf{curl} (\mathbf A) \end{equation}
-
Electric Field Intensity
- transient case: \begin{equation} \mathbf E = \frac{\partial \mathbf A}{\partial t} \end{equation}
- harmonic case: \begin{equation} \mathbf E = j \omega \mathbf A \end{equation}
-
Magnetic Eddy Current Density
- transient case: \begin{equation} \mathbf J_{\textrm{eddy}}= \gamma \mathbf E = \gamma \frac{\partial \mathbf A}{\partial t} \end{equation}
- harmonic case: \begin{equation} \mathbf J_{\textrm{eddy}}= \gamma \mathbf E = \gamma j \omega \mathbf A \end{equation}
-
Magnetic Eddy Power Density (Joule Loss Power Density)
- transient case: \begin{equation} P_{\textrm{Joule}} = \mathbf J_{\textrm{eddy}} \cdot \mathbf E = \mathbf J_{\textrm{eddy}} \cdot \frac{\partial \mathbf A}{\partial t} =\gamma \frac{\partial \mathbf A}{\partial t} \cdot \frac{\partial \mathbf A}{\partial t} \end{equation}
- harmonic case: Averaged over one period of the base excitation frequency \begin{equation} \tilde{P}_{\textrm{Joule}} = \frac{1}{2} \gamma \omega^2 \mathbf A \cdot \mathbf A^* = \frac{1}{2} \gamma \omega^2 |\mathbf A| \end{equation}
-
Total Current Density \begin{equation} \mathbf J_{\textrm{total}} = \mathbf J_{\textrm{impressed}} + \mathbf J_{\textrm{eddy}} \end{equation}
-
Lorentz Force Density \begin{equation} \mathbf F_{\textrm{Lorentz}} = \mathbf J_{\textrm{total}} \times \mathbf B = \mathbf J_{\textrm{total}} \times \mathbf B \end{equation}
Surface Element Results¶
- Magnetic Force Maxwell Density on Surface
Surface Region Results¶
-
Magnetic Eddy Current on Surface \begin{equation} J = \int_{A} \mathbf J_{\textrm{eddy}} \cdot \mathbf n\, dS \end{equation}
-
Magnetic Flux on Surface \begin{equation} \Phi = \int_{A} \mathbf B \cdot \mathbf n\, dS \end{equation}
-
Magnetic Force (Virtual Work Principle) on Surface
-
Magnetic Force Maxwell on Surface
Region Results¶
-
Magnetic Energy Implemented as \begin{equation} E_{\textrm{mag}} = \frac{1}{2} \sum_{e=1}^{n_e} \int_{\Omega^e} \mathbf A^{\textrm{T}} \cdot \mathbf K^e \cdot \mathbf A, \end{equation} where n_e is the number of elements in the considered region, \mathbf A the element solution vector and \mathbf K^e the element stiffness matrix. This corresponds to \begin{equation} E_{\textrm{mag}} = \frac{1}{2} \int_{\Omega} \mathbf B \cdot \mathbf B, \end{equation}
-
Magnetic Eddy Power \begin{equation} E_{\textrm{eddy}} = \int_{\Omega} \mathbf J_{\textrm{eddy}} \cdot \mathbf E \end{equation}
-
Joule Loss Power
- transient case: The same as Magnetic Eddy Power \begin{equation} E_{\textrm{Joule}} = \int_{\Omega} P_{\textrm{Joule}} = \int_{\Omega} \mathbf J_{\textrm{eddy}} \cdot \mathbf E \end{equation}
- harmonic case: Averaged over one period of the base excitation frequency
\begin{equation}
E_{\textrm{Joule}} = \int_{\Omega} \tilde{P}_{\textrm{Joule}}
\end{equation}
-
Lorentz Force \begin{equation} F_{\textrm{Lorentz}} = \int_{\Omega} \mathbf F_{\textrm{Lorentz}} \end{equation}
Results for magneticPDE (nodal element version)¶
Tutorials¶
References¶
-
Miklos Kuczmann. Potential Formulations in Magnetics: Applying the Finite Element Method. Laboratory of Electromagnetic Fields, Szechenyi Istvan University, Gyor, Hungary, 2009. URL: http://maxwell.sze.hu/docs/C4.pdf. ↩