 |
Logbook
(07-04-2025)
Static problems
|
|
|
Logbook
(07-04-2025)
Static problems
|
|
Classes:
1) BatchSCH
2) SolverSCH
3) ExactSolutionSCH_PHI
4) SettingsSCH
Files:
1) sch/src/main.cpp
2) sch/include/solver.hpp
3) sch/src/solver.cpp
4) sch/include/exact_solution.hpp
5) sch/src/exact_solution.cpp
6) sch/src/static_scalar_input.cpp
7) sch/include/settings.hpp
In the first four numerical experiments the surface free-charge density \(\kappa_f\) is set to zero. As discussed in the Method of manufactured solutions (mms/) numerical experiment, the reason for this that it is quite difficult to write down a manufactured or exact solution that captures all the aspects of the static scalar boundary value problem. Consequently, the test of the solver implemented by the class template StaticScalarSolver::Solver is split in parts. Most of its functionality is tested in the Method of manufactured solutions (mms/) numerical experiment. The interface between dissimilar dielectric materials is tested in Interface between dielectrics (int/) numerical experiment. This numerical experiment is made to test the ability of the solver to handle problems with the surface free-charge density, \(\kappa_f\).
We presume that the surface free-charge density is built up on a certain surface within the problem domain. We will call this surface an interface. We presume that all parts of the interface are made up of cell faces. That is, no segment of the interface runs though a cell.
The following textbook problem is solved in the two-dimensional version of the experiment. (Section 6.1.9. Surface charge on a cylindrical interface)
Free charge is uniformly distributed on a surface of an infinity long cylinder of a radius \(a\). This cylinder is enclosed by an infinitely long conducting tube with an inner radius of \(b\). The cylinder and the tube are coaxial. The tube is grounded. The density of the surface free charge, \(\kappa_f\), is high enough to induce a potential of \(\Phi_0\) Volts on the surface of the inner cylinder. All materials inside the conducting tube have permittivity of \(\epsilon_0\). The figure below illustrates this configuration.
We would like to calculate the electric potential.
The solution to this problem reads
\begin{equation} \kappa_f = \frac{\epsilon_0 \Phi_0}{a\ln(b/a)}, \end{equation}
\begin{equation} \Phi = \left\{ \begin{aligned} &\Phi_0 &\text{if } \text{ } & r \le a,\\ &\frac{\ln(b/r)}{\ln(b/a)} \Phi_0 &\text{if } \text{ } & r \ge a,\\ \end{aligned} \right. \end{equation}
\begin{equation} \vec{\nabla} \Phi = \left\{ \begin{aligned} & 0 &\text{if } \text{ } & r \le a,\\ & -\frac{\Phi_0}{\ln(b/a)} \frac{1}{r}\hat{r}&\text{if } \text{ } & r \ge a.\\ \end{aligned} \right. \end{equation}
As soon as this problem exhibits translation symmetry, we solve it numerically on a two-dimensional planar domain. We remove the conducting tube from the problem domain and model it by setting up a Dirichlet boundary condition \(\Phi = 0\) on the outer boundary of the problem domain.
The following textbook problem is solved in the three-dimensional version of the experiment. (Section 6.1.10. Surface charge on a spherical interface)
Free charge is uniformly distributed on a surface of a sphere of a radius \(a\). The sphere is enclosed by a conducting spherical shell with an inner radius of \(b\). The sphere and the shell are concentric. The shell is grounded. The density of the surface free charge, \(\kappa_f\), is high enough to induce a potential of \(\Phi_0\) Volts on the surface of the inner sphere. All materials inside the conducting spherical shell have permittivity of \(\epsilon_0\). The figure below illustrates this configuration.
We would like to calculate the electric potential.
The solution to this problem reads
\begin{equation} \kappa_f = \frac{b}{a(b - a)} \epsilon_0 \Phi_0, \end{equation}
\begin{equation} \Phi = \left\{ \begin{aligned} &\Phi_0 &\text{if } \text{ } & r \le a,\\ &\frac{ab}{b-a}\Bigg(\frac{1}{r}-\frac{1}{b}\Bigg) \Phi_0 &\text{if } \text{ } & r \ge a,\\ \end{aligned} \right. \end{equation}
\begin{equation} \vec{\nabla} \Phi = \left\{ \begin{aligned} &0 &\text{if } \text{ } & r \le a,\\ &-\frac{ab}{b-a} \frac{\Phi_0}{r^2}\hat{r} &\text{if } \text{ } & r \ge a.\\ \end{aligned} \right. \end{equation}
We solve this problem numerically in a three-dimensional domain. We remove the conducting shell from the problem domain and model it by setting up a Dirichlet boundary condition \(\Phi = 0\) on the outer boundary of the problem domain.
Both versions of this numerical experiment, the two- and three- dimensional, solve the following boundary value problem by the finite element method:
\begin{equation} \begin{array}{lrcll} \text{ }&- \vec{\nabla} \cdot \big( \epsilon_0 \vec{\nabla} \Phi \big)= 0 & \text{in} & \Omega & \text{(i)},\\ \text{(e)} &\Phi = 0 & \text{on} & \Gamma_{D1} & \text{(ii)},\\ \text{(e)} &\Phi_{+} = \Phi_{-} & \text{on} & \Gamma_{I1} & \text{(iii)},\\ \text{(n)}&\epsilon_{0}\hat{n}\cdot\vec{\nabla}\Phi_{+}-\epsilon_{0}\hat{n}\cdot\vec{\nabla}\Phi_{-}=-k_f&\text{on}&\Gamma_{I1}&\text{(iv)}. \end{array} \end{equation}
This boundary value problem is a special case of the general static scalar boundary value problem.
This experiment is implemented in accordance with the base code structure. This time we use second-order mapping and attach a spherical manifold to all faces of the mesh. That is, we declare the spherical manifold in the SolverSCH template,
and attach it to the mesh as the following:
That is, the spherical manifold is attached to a face if all the vertices of the face are equidistant from the origin and the center of the cell the face belongs to is outside a circle (ball) of a radius of \(d1\) (the half-side of the square (cube) in the middle of the mesh). Then we pass 2 as the second parameter to the constructor of SolverSCH class template to enable second-order mapping,
The build process generates two executable files: sch-circle and sch-sphere. That is, one executable file for the two-dimensional experiment, sch-circle, and one executable file for the three-dimensional experiment, sch-sphere. To rebuild them change into sch/build/Release directory and execute the following:
Then all executable files must be executed again. This can be done by changing into the sch/bin/Release directory and executing run-all script there,
This will generate various files in the sch/bin/Release/Data directory. Among the generated files there are vtu files that can be viewed with a help of ParaView software package of the Kitware, Inc. The data files in tex and txt format contain the convergence tables. The sch/bin/Release/Data directory also contains gpi files. They can be used to visualize the calculated potential along the \(x\) axis.
Note that executable files require a set of meshes to be present in the sch/gmsh/data directory. If they are missing, they can be generated anew. This can be done by changing into sch/gmsh directory and executing the following:
This will generate a set of globally refined meshes in sch/gmsh/data.
The SettingsSCH 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. The last can be useful when debugging.
To work properly, the StaticScalarSolver::Solver needs some way of discriminating the cell faces that belong to the interface \(\Gamma_{I1}\), so it can distribute the surface free-charge density, \(\kappa_f\), on it. This discrimination is accomplished by invoking the mechanism of cell and face user ID which is available in deal.II. The StaticScalarSolver::Solver class template calculates integrals related to \(\kappa_f\) only if both, the user cell ID and the user face ID do not equal zero. Furthermore, StaticScalarSolver::Solver passes the values of the user cell and face IDs to StaticScalarSolver::FreeSurfaceCharge::value_list, so the algorithm implemented by this function knows which cell face is currently processed. Therefore, from the user's standpoint \(\kappa_f\) is handled in two places in the code: (i) in the private function SolverSCH::mark_materials (it is called just after a mesh is uploaded) where cells and faces on the interface are marked; and (ii) StaticScalarSolver::TheFreeSurfaceCharge, where proper values of \(\kappa_f\) are calculated based on cell and face user IDs.
Note that each face on an interface always belongs to two cells. Therefore, a care must be taken that each face on the interface is handled only ones. This is done by marking a face on the interface only if it belongs to a cell center of which is located within the circle (ball) of radius \(a\). The code snippet below explains this in more detail.
Then the StaticScalarSolver::FreeSurfaceCharge<2>::value_list and StaticScalarSolver::FreeSurfaceCharge<3>::value_list calculate the values of \(\kappa_f\) in the two- and three-dimensional versions of the experiment, respectively:
All calculation have been done with the following settings: \(\Phi_0 = 1.0[V]\), \(a=0.5[m]\), and \(b=1.0[m]\), see SettingsSCH.
The plot below illustrates the numerically calculated electric potential.
The table below illustrates the error norms and the convergence orders.
| p | r | cells | dofs | \(\|e\|_{L^2}\) | \(\alpha_{L^2}\) | \(\|e\|_{H^1}\) | \(\alpha_{H^1}\) |
|---|---|---|---|---|---|---|---|
| 1 | 5 | 320 | 337 | 5.51e-03 | - | 1.59e-01 | - |
| 1 | 6 | 500 | 521 | 3.54e-03 | 1.98 | 1.27e-01 | 0.99 |
| 1 | 7 | 720 | 745 | 2.46e-03 | 1.99 | 1.06e-01 | 0.99 |
| 1 | 8 | 980 | 1009 | 1.81e-03 | 1.99 | 9.11e-02 | 0.99 |
| 2 | 5 | 320 | 1313 | 1.54e-04 | - | 8.00e-03 | - |
| 2 | 6 | 500 | 2041 | 7.93e-05 | 2.96 | 5.15e-03 | 1.97 |
| 2 | 7 | 720 | 2929 | 4.61e-05 | 2.98 | 3.59e-03 | 1.98 |
| 2 | 8 | 980 | 3977 | 2.91e-05 | 2.98 | 2.65e-03 | 1.98 |
| 3 | 5 | 320 | 2929 | 7.81e-06 | - | 4.67e-04 | - |
| 3 | 6 | 500 | 4561 | 3.23e-06 | 3.96 | 2.42e-04 | 2.95 |
| 3 | 7 | 720 | 6553 | 1.56e-06 | 3.97 | 1.41e-04 | 2.96 |
| 3 | 8 | 980 | 8905 | 8.47e-07 | 3.98 | 8.92e-05 | 2.97 |
The following notations were used in the table header:
The plot below illustrates a numerically calculated electric potential.
The table below illustrates the error norms and the convergence orders generated by the three-dimensional version of the numerical experiment.
| p | r | cells | dofs | \(\|e\|_{L^2}\) | \(\alpha_{L^2}\) | \(\|e\|_{H^1}\) | \(\alpha_{H^1}\) |
|---|---|---|---|---|---|---|---|
| 1 | 5 | 832 | 909 | 1.29e-02 | - | 3.86e-01 | - |
| 1 | 6 | 1625 | 1736 | 8.36e-03 | 1.96 | 3.10e-01 | 0.98 |
| 1 | 7 | 2808 | 2959 | 5.83e-03 | 1.97 | 2.59e-01 | 0.99 |
| 1 | 8 | 4459 | 4656 | 4.30e-03 | 1.98 | 2.22e-01 | 0.99 |
| 2 | 5 | 832 | 6905 | 5.85e-04 | - | 3.00e-02 | - |
| 2 | 6 | 1625 | 13371 | 3.01e-04 | 2.97 | 1.94e-02 | 1.96 |
| 2 | 7 | 2808 | 22981 | 1.75e-04 | 2.98 | 1.35e-02 | 1.97 |
| 2 | 8 | 4459 | 36359 | 1.10e-04 | 2.99 | 9.98e-03 | 1.98 |
| 3 | 5 | 832 | 22981 | 1.57e-04 | - | 3.38e-03 | - |
| 3 | 6 | 1625 | 44656 | 6.51e-05 | 3.95 | 1.75e-03 | 2.95 |
| 3 | 7 | 2808 | 76915 | 3.16e-05 | 3.96 | 1.02e-03 | 2.97 |
| 3 | 8 | 4459 | 121864 | 1.71e-05 | 3.97 | 6.43e-04 | 2.98 |
The following notations were used in the table header:
The results of this numerical experiment suggest that the solver implemented by the class template StaticScalarSolver::Solver handles the problems with surface charge density well.
1.9.1