Solves static vector boundary value problem. More...

#include <static_vector_solver_ii.hpp>

Public Member Functions
	Solver2 (unsigned int p, unsigned int mapping_degree, const Triangulation< dim > &triangulation_T, const DoFHandler< dim > &dof_handler_T, const Vector< double > &solution_T, unsigned int type_of_pde_rhs=0, double eta_squared=0.0, std::string fname="data", const Function< dim > *exact_solution=nullptr, bool print_time_tables=false, bool project_exact_solution=false, bool write_higher_order_cells=false)
	The only constructor. More...

virtual void	fill_dirichlet_stack ()=0
	Initializes the data member StaticVectorSolver::Solver2::dirichlet_stack. More...

virtual void	solve ()=0
	Solves the system of linear equations.

void	setup ()
	Initializes system matrix and the right-hand side vector. More...

void	assemble ()
	Assembles the system matrix and the right-hand side vector.

void	compute_error_norms ()
	Computes error norms.

void	project_exact_solution_fcn ()
	Projects exact solution. More...

void	save () const
	Saves simulation results into a vtk file. More...

void	save_matrix_and_rhs_to_csv (std::string fname) const
	Saves the system matrix and the right-hand side into a csv file. More...

void	clear ()
	Releases computer memory associated with system matrix and right-hand side.

const DoFHandler< dim > &	get_dof_handler () const
	Returns a reference to a dof handler.

const Vector< double > &	get_solution () const
	Returns a reference to solution.

unsigned int	get_n_cells () const
	Returns the number of active cells in the mesh.

unsigned int	get_n_dofs () const
	Returns the total amount of the degrees of freedom.

double	get_L2_norm () const
	Returns \(L^2\) error norm.

double	get_Linfty_norm () const
	Returns \(L^{\infty}\) error norm.

unsigned int	get_mapping_degree () const
	Returns degree of the interpolating Lagrange polynomials used for mapping from the reference cell to the real mesh cell and back.

void	run ()
	Runs the simulation. More...

Protected Attributes
const Triangulation< dim > &	triangulation_T
	Reference to the mesh.

const DoFHandler< dim > &	dof_handler_T
	Reference to the dof handler that describes the current vector potential, \(\vec{T}\).

const Vector< double > &	solution_T
	Reference to the degrees of freedom that describe the current vector potential, \(\vec{T}\).

std::map< types::boundary_id, const Function< dim > * >	dirichlet_stack
	A map that contains pairs of boundary IDs and the corresponding Dirichlet boundary conditions. All boundary IDs must be odd.

const FE_Nedelec< dim >	fe
	The finite elements.

DoFHandler< dim >	dof_handler
	The finite elements.

Vector< double >	solution
	The solution vector, that is, degrees of freedom yielded by the simulation.

Vector< double >	projected_exact_solution
	The projected exact solution vector.

AffineConstraints< double >	constraints
	The constraints associated with the Dirichlet boundary conditions.

SparsityPattern	sparsity_pattern
	The sparsity pattern of the system matrix.

SparseMatrix< double >	system_matrix
	The system matrix.

Vector< double >	system_rhs
	The system right-hand side vector.

double	L2_norm
	The \(L^2\) norm.

double	Linfty_norm
	The \(L^{\infty}\) norm.

Detailed Description

template<int dim, int stage = 1>
class StaticVectorSolver::Solver2< dim, stage >

Solves static vector boundary value problem.

Implements the following recipes:

(1) Recipe for static vector solver in 3D
(2) Recipe for static vector solver in 2D

This class template is intended to be a general solver for problems in magnetostatics that can be formulated in terms of the magnetic vector potential, \(\vec{A}\). The Bossavit's diagrams below illustrate the partial differential equations that can be solved with a help of this class template.

This class template is very similar to StaticVectorSolver::Solver1. The main difference between StaticVectorSolver::Solver1 and StaticVectorSolver::Solver2 can be described as the following. All inputs to StaticVectorSolver::Solver1 must be given in a form of an analytical expressions. The same holds for StaticVectorSolver::Solver2 with an exception of the input that describes the right-hand side of the partial differential equation. The input on the right-hand side of the partial differential equation fed into the StaticVectorSolver::Solver2 class template must be a current vector potential, \(\vec{T}\), expressed in a form of a field function. That is, \(\vec{T}\) fed to StaticVectorSolver::Solver2 is a result of another deal.II simulation. This allows to solve the curl-curl equation for \(\vec{A}\) observing the compatibility condition. The diagram below illustrates how this can be done in two- and three- dimensions.

Recall that we do not gauge the magnetic vector potential, \(\vec{A}\), explicitly. Consequently, the conservative part of simulated \(\vec{A}\) is unknown. For this reason, it does not make mush sense to evaluate or visualize \(\vec{A}\). Instead, we evaluate the simulated \(\vec{A}\) indirectly by computing and evaluating the corresponding magnetic field, \(\vec{B}\), in the third stage, see the diagram above.

The magnetic vector potential is modeled by the FE_Nedelec finite elements.

A user of this class is supposed to do the following.

Derive a class from StaticVectorSolver::Solver2.
Override the virtual destructor.
Call constructor StaticVectorSolver::Solver2.
Override member functions
virtual void fill_dirichlet_stack() = 0;

virtual void solve() = 0;

StaticVectorSolver::Solver2::solve
virtual void solve()=0
Solves the system of linear equations.

StaticVectorSolver::Solver2::fill_dirichlet_stack
virtual void fill_dirichlet_stack()=0
Initializes the data member StaticVectorSolver::Solver2::dirichlet_stack.
Implement the member functions
void value_list(...);

of the classes StaticVectorSolver::TheCoefficient, StaticVectorSolver::Gamma, StaticVectorSolver::RobinRhs, and StaticVectorSolver::FreeSurfaceCurrent.
Implement member function
double value(...) const;

of the class StaticVectorSolver::Weight.
Explicitly call
StaticVectorSolver::run();

Alternatively, the user may call individual member functions (such as fill_dirichlet_stack(), setup(), etc.) in a proper order.

Note that there is no make_mesh() function this time. The StaticVectorSolver::Solver2 class template reuses the meshes created by other solvers. To this end, the reference to the mesh is passed as the third input parameter to the constructor.

It is assumed that the object that has been used to calculate \(\vec{T}\) (or \(T\) in 2D) is still in the computer memory such that the triangulation, the degrees of freedom, and the handler of the degrees of freedom are accessible while \(\vec{A}\) is computed. An object of the StaticVectorSolver::Solver2 class template reuses the triangulation by creating an additional dof handler associated with the FE_Nedelec finite elements. That is, \(\vec{T}\) and \(\vec{A}\) share the same triangulation. Two separate DoFHandler objects are used for \(\vec{T}\) and \(\vec{A}\).

The boundaries of the mesh must be labeled such that the boundary_id() member function of a face object returns the corresponding boundary ID. The boundary ID's must obey the following convention.

The Dirichlet boundary conditions are applied on the boundaries with odd boundary ID's.
The Robin boundary conditions are applied on the boundaries with even boundary ID's. The boundary ID's in this case must be greater than zero.
No boundary condition is applied on a boundary with zero ID. Applying no boundary condition is as good as applying the homogeneous Neumann boundary condition, \(\big(1/\mu\big)\hat{n}\times\big(\vec{\nabla}\times\vec{A}\big)=0\), as implicitly implied by the first term of the functional. This boundary condition can also be imposed by assigning to a boundary an even ID greater than zero, and setting \( \gamma \) and \(\vec{Q}\) to zero in the value_list methods of the classes StaticVectorSolver::Gamma and StaticVectorSolver::RobinRhs.

The constructor's argument type_of_pde_rhs switches the operation of the class template between following two modes:

The parameter type_of_pde_rhs equals any unsigned integer except for 2. In this mode both integrals of the functional associated with the right-hand side of the partial differential equation are computed. The integrals associated with the righ-hand side read
\[\iiint_{\Omega} \vec{T}\cdot\bigg(\vec{\nabla}\times\vec{A}\bigg) dV - \underbrace{ \iint_{\Gamma_{\Omega}}\vec{T}\cdot\bigg(\hat{n}\times\vec{A}\bigg)dS }_{\text{Boundary integral}} \]
in a three dimensional space and
\[ \iint_{\Omega}T\bigg(\vec{\nabla}\overset{S}{\times}\vec{A}\bigg)dS- \underbrace{ \int_{\Gamma_{\Omega}}T\bigg(\hat{n}\overset{S}{\times}\vec{A}\bigg)dl }_{\text{Boundary integral}} \]
in a two-dimensional space.
type_of_pde_rhs=2. In this mode the boundary integral in the expression above is not computed. This can be done if \(\vec{T}\) (or \(T\) in 2D) is forced to zero by the homogeneous Dirichlet boundary condition in the first stage of the simulation. This mode saves simulation time.

The algorithm walks synchronously through both DoFHandler objects and assembles the system matrix and the right-hand side. The algorithm utilizes the WorkStream technology of deal.II. The amount of threads used can be limited as the following.

#include <deal.II/base/multithread_info.h>
...
MultithreadInfo::set_thread_limit(nr_threads_max);

Note

Application examples:

mms-vt-i/, mms-vt-ii/, ssol-ii/, ssol-iii/.

Definition at line 223 of file static_vector_solver_ii.hpp.

Constructor & Destructor Documentation

◆ Solver2()

template<int dim, int stage = 1>

StaticVectorSolver::Solver2< dim, stage >::Solver2	(	unsigned int	p,
		unsigned int	mapping_degree,
		const Triangulation< dim > &	triangulation_T,
		const DoFHandler< dim > &	dof_handler_T,
		const Vector< double > &	solution_T,
		unsigned int	type_of_pde_rhs = `0`,
		double	eta_squared = `0.0`,
		std::string	fname = `"data"`,
		const Function< dim > *	exact_solution = `nullptr`,
		bool	print_time_tables = `false`,
		bool	project_exact_solution = `false`,
		bool	write_higher_order_cells = `false`
	)

inline

The only constructor.

Parameters

[in]	p	- Degree of the FE_Nedelec finite elements.
[in]	mapping_degree	- The degree of the interpolating Lagrange polynomials used for mapping. Setting it to 1 will do in the most of the cases. Note, that it makes sense to attach a meaningful manifold to the triangulation if this parameter is greater than 1.
[in]	triangulation_T	- A reference to the mesh on which the current vector potential, \(\vec{T}\), has been computed.
[in]	dof_handler_T	- A reference to the DoFHandler object that describes the current vector potential, \(\vec{T}\).
[in]	solution_T	- A reference to the degrees of freedom that describe the current vector potential, \(\vec{T}\).
[in]	type_of_pde_rhs	- If equals 2, the boundary integral \(I_{b3-2}\) is neglected.
[in]	eta_squared	- The gauging parameter, \(\eta^2\), in the partial differential equation.
[in]	fname	- The name of the output files without extension. Names of the output files will be generated by appending simulation conditions to this string.
[in]	exact_solution	- Points to an object that describes the exact solution to the problem. It is needed for calculating error norms. It is a responsibility of the user to make sure that the object exists at the time of the execution of run() or compute_error_norms().
[in]	print_time_tables	- If true, prints time tables on the screen.
[in]	project_exact_solution	- If true, projects the exact solution onto the space spanned by the Nedelec finite elements and saves the result into the output file next to the solution. This may be useful for debugging purposes as a comparison between the projected exact solution and the solution to the boundary value problem can yield a hint on where to search for bugs.
[in]	write_higher_order_cells	- Switches between the two modes of operation of the save() function, see the description of save().

Definition at line 265 of file static_vector_solver_ii.hpp.

Member Function Documentation

◆ fill_dirichlet_stack()

template<int dim, int stage = 1>

virtual void StaticVectorSolver::Solver2< dim, stage >::fill_dirichlet_stack ( )

pure virtual

Initializes the data member StaticVectorSolver::Solver2::dirichlet_stack.

This function must be overridden by the user. It must initialize the stack of the Dirichlet boundary conditions. For example:

using namespace dealii;
 
const types::boundary_id boundary_id_1 = 1;
const types::boundary_id boundary_id_2 = 3;
 
const DirichletFunction1<dim> dirichlet_bc_1;
const DirichletFunction2<dim> dirichlet_bc_2;
 
template<>
void SolverMyProblem<dim>::fill_dirichlet_stack()
{
  Solver<dim, stage>::dirichlet_stack =
    {{boundary_id_1, & dirichlet_bc_1},
     {boundary_id_2, & dirichlet_bc_2}};
}

The boundary IDs must be odd numbers, see above the convention on the boundary IDs.

◆ project_exact_solution_fcn()

template<int dim, int stage>

void StaticVectorSolver::Solver2< dim, stage >::project_exact_solution_fcn

Projects exact solution.

The mesh and the finite elements, are the same as are used for the numerical solution of the boundary value problem. The exact solution will be saved in the output file next to the numerical solution to the boundary value problem. This function works properly only if the exact solution is submitted to the constructor via the input parameter exact_solution and project_exact_solution=true.

Definition at line 1157 of file static_vector_solver_ii.hpp.

◆ run()

template<int dim, int stage = 1>

void StaticVectorSolver::Solver2< dim, stage >::run ( )

inline

Runs the simulation.

Executes the following member functions in a proper order: fill_dirichlet_stack(), setup(), assemble(), solve(), project_exact_solution_fcn(), compute_error_norms(), save().

Definition at line 478 of file static_vector_solver_ii.hpp.

◆ save()

template<int dim, int stage>

void StaticVectorSolver::Solver2< dim, stage >::save

Saves simulation results into a vtk file.

The following data are saved:

The calculated potential under the name "VectorField".
The \(L^2\) error norm associated with the calculated potential under the name "L2norm". One value per mesh cell is saved.
The \(L^{\infty}\) error norm associated with the calculated potential under the name "LinftyNorm". One value per mesh cell is saved.
The exact solution expressed as a linear combination of the shape functions of the FE_Nedelec finite elements is saved under the name "VectorFieldExact". The "VectorField" and "VectorFieldExact" are modeled by exactly the same finite elements.

The "L2norm", \(L^{\infty}\), and "VectorFieldExact" are saved only if an exact solution is submitted to the constructor. Moreover, "VectorFieldExact" is calculated and saved only if project_exact_solution = true.

If write_higher_order_cells = false, the name of the file is computed by appending ".vtk" to the string contained by the parameter fname passed to the constructor. The file can be inspected with a help of VisIt or Paraview. Higher-order cells are saved as regular quadrilaterals and hexahedra. If write_higher_order_cells = true, the name of the file is computed by appending ".vtu" to the string contained by the parameter fname. The data is saved preserving the higher-order cells. The file can be viewed with a help of Paraview version 5.5.0 or higher. VisIt cannot load higher-order cells.

Definition at line 1174 of file static_vector_solver_ii.hpp.

◆ save_matrix_and_rhs_to_csv()

template<int dim, int stage>

void StaticVectorSolver::Solver2< dim, stage >::save_matrix_and_rhs_to_csv ( std::string fname ) const

Saves the system matrix and the right-hand side into a csv file.

All the zeros are included into the csv files. This is a very dumb and inefficient way of saving sparse matrices. On the positive side - it is very easy and straightforward to read the csv files. This function may be useful for debugging. One can assemble the system on a coarse mesh (so there are a few mesh cells and the system matrix is small) and export the system matrix together with the right-hand side into another program such as Matlab or GNU Octave for an analysis.

Parameters

[in] fname - A stem of the names of the output files. The matrix will be saved into fname_matrix.csv file. The right-hand side will be save into fname_rhs.csv file.

Definition at line 1227 of file static_vector_solver_ii.hpp.

◆ setup()

template<int dim, int stage>

void StaticVectorSolver::Solver2< dim, stage >::setup

Initializes system matrix and the right-hand side vector.

Initialises StaticVectorSolver::Solver2::system_matrix, StaticVectorSolver::Solver2::system_rhs and other arrays. Applies the Dirichlet boundary conditions by constraining the system matrix. Distributes degrees of freedom.

Definition at line 691 of file static_vector_solver_ii.hpp.

The documentation for this class was generated from the following file:

static/shared/include/static_vector_solver_ii.hpp

Public Member Functions

Protected Attributes

Detailed Description

template<int dim, int stage = 1> class StaticVectorSolver::Solver2< dim, stage >

Constructor & Destructor Documentation

◆ Solver2()

Member Function Documentation

◆ fill_dirichlet_stack()

◆ project_exact_solution_fcn()

◆ run()

◆ save()

◆ save_matrix_and_rhs_to_csv()

◆ setup()

template<int dim, int stage = 1>
class StaticVectorSolver::Solver2< dim, stage >