v0.6.20
Classes | Public Member Functions | Public Attributes | List of all members
MagneticElement Struct Reference

Implementation of magneto-static problem (basic Implementation)Look for theory and details here: More...

#include <users_modules/basic_finite_elements/magnetostatic/src/MagneticElement.hpp>

Collaboration diagram for MagneticElement:
[legend]

Classes

struct  BlockData
 data structure storing material constants, model parameters, matrices, etc. More...
 
struct  OpCurlCurl
 calculate and assemble CurlCurl matrix

\[ \mathbf{A} = \int_\Omega \mu^{-1} \left( \nabla \times \mathbf{u} \cdot \nabla \times \mathbf{v} \right) \textrm{d}\Omega \]

where

\[ \mathbf{u} = \nabla \times \mathbf{B} \]

where \(\mathbf{B}\) is magnetic flux and \(\mu\) is magnetic permeability. More...

 
struct  OpNaturalBC
 calculate essential boundary conditions

\[ \mathbf{A} = \int_{\partial\Omega} \ \mathbf{u} \cdot \mathbf{j}_i \textrm{d}{\partial\Omega} \]

where \(\mathbf{j}_i\) is current density function. More...

 
struct  OpPostProcessCurl
 calculate and assemble CurlCurl matrix More...
 
struct  OpStab
 calculate and assemble stabilization matrix

\[ \mathbf{A} = \int_\Omega \epsilon \mathbf{u} \cdot \mathbf{v} \textrm{d}\Omega \]

where \(\epsilon\) is regularization parameter. More...

 
struct  TriFE
 define surface element More...
 
struct  VolumeFE
 definition of volume element More...
 

Public Member Functions

 MagneticElement (MoFEM::Interface &m_field)
 
virtual ~MagneticElement ()
 
MoFEMErrorCode getNaturalBc ()
 get natural boundary conditions More...
 
MoFEMErrorCode getEssentialBc ()
 get essential boundary conditions (only homogenous case is considered) More...
 
MoFEMErrorCode createFields ()
 build problem data structures More...
 
MoFEMErrorCode createElements ()
 create finite elementsCreate volume and surface element. Surface element is used to integrate natural boundary conditions. More...
 
MoFEMErrorCode createProblem ()
 create problemProblem is collection of finite elements. With the information on which fields finite elements operates the matrix and left and right hand side vector could be created. More...
 
MoFEMErrorCode destroyProblem ()
 destroy Distributed mesh manager More...
 
MoFEMErrorCode solveProblem (const bool regression_test=false)
 solve problemCreate matrices; integrate over elements; solve linear system of equations More...
 
MoFEMErrorCode postProcessResults ()
 post-process results, i.e. save solution on the mesh More...
 

Public Attributes

MoFEM::InterfacemField
 
BlockData blockData
 

Detailed Description

Implementation of magneto-static problem (basic Implementation)

Look for theory and details here:

[26] <www.hpfem.jku.at/publications/szthesis.pdf>

Election file and all other problem related file are here maxwell_element.

Todo:

Extension for mix formulation

Use appropriate pre-conditioner for large problems

Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 42 of file MagneticElement.hpp.

Constructor & Destructor Documentation

◆ MagneticElement()

MagneticElement::MagneticElement ( MoFEM::Interface m_field)
Examples:
MagneticElement.hpp.

Definition at line 70 of file MagneticElement.hpp.

70 : mField(m_field) {}
MoFEM::Interface & mField

◆ ~MagneticElement()

virtual MagneticElement::~MagneticElement ( )
virtual
Examples:
MagneticElement.hpp.

Definition at line 71 of file MagneticElement.hpp.

71 {}

Member Function Documentation

◆ createElements()

MoFEMErrorCode MagneticElement::createElements ( )

create finite elementsCreate volume and surface element. Surface element is used to integrate natural boundary conditions.

Returns
error code
Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 242 of file MagneticElement.hpp.

242  {
244  // //Elements
253  "MESH_NODE_POSITIONS");
262  blockData.feNaturalBCName, "MESH_NODE_POSITIONS");
264  blockData.feName);
267  // build finite elemnts
269  // build adjacencies
272  }
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string &name_row)=0
set field col which finite element use
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:491
virtual MoFEMErrorCode build_finite_elements(int verb=-1)=0
Build finite elementsBuild finite element data structures. Have to be run before problem and adjacenc...
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Common.hpp:147
MoFEM::Interface & mField
#define CHKERR
Inline error check.
Definition: definitions.h:610
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL)=0
add finite element
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string &name_row)=0
set field row which finite element use
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:436
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=-1)=0
build adjacencies

◆ createFields()

MoFEMErrorCode MagneticElement::createFields ( )

build problem data structures

Returns
error code
Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 166 of file MagneticElement.hpp.

166  {
168 
169  // Set entities bit level. each entity has bit level depending for example
170  // on refinement level. In this case we do not refine mesh or not do
171  // topological changes, simply set refinement level to zero on all entities.
172 
173  CHKERR mField.getInterface<BitRefManager>()->setBitRefLevelByDim(
174  0, 3,BitRefLevel().set(0));
175 
176  // add fields
178  1);
179  CHKERR mField.add_field("MESH_NODE_POSITIONS", H1, AINSWORTH_LEGENDRE_BASE,
180  3);
181  // meshset consisting all entities in mesh
182  EntityHandle root_set = mField.get_moab().get_root_set();
183  // add entities to field
184  CHKERR mField.add_ents_to_field_by_type(root_set, MBTET,
186 
187  // // The higher-order gradients can be gauged by locally skipping the
188  // // corresponding degrees of freedom and basis functions in the
189  // higher-order
190  // // edge-based, face-based and cell-based finite element subspaces.
191  //
192  // Range tris,edges;
193  // CHKERR mField.get_moab().get_entities_by_type(root_set,MBTRI,tris,true);
194  // mField.get_moab().get_entities_by_type(root_set,MBEDGE,edges,true);
195  //
196  // // Set order in volume
197  // Range bc_ents = unite(blockData.naturalBc,blockData.essentialBc);
198  // Range vol_ents = subtract(unite(tris,edges),bc_ents);
199  // CHKERR
200  // mField.set_field_order(vol_ents,blockData.fieldName,blockData.oRder); int
201  // gauged_order = 1; CHKERR
202  // mField.set_field_order(bc_ents,blockData.fieldName,gauged_order);
203 
204  // Set order on tets
206  blockData.oRder);
208  blockData.oRder);
209  CHKERR mField.set_field_order(root_set, MBEDGE, blockData.fieldName,
210  blockData.oRder);
211 
212  // Set geometry approximation ordered
213  CHKERR mField.add_ents_to_field_by_type(root_set, MBTET,
214  "MESH_NODE_POSITIONS");
215  CHKERR mField.set_field_order(root_set, MBTET, "MESH_NODE_POSITIONS", 2);
216  CHKERR mField.set_field_order(root_set, MBTRI, "MESH_NODE_POSITIONS", 2);
217  CHKERR mField.set_field_order(root_set, MBEDGE, "MESH_NODE_POSITIONS", 2);
218  CHKERR mField.set_field_order(root_set, MBVERTEX, "MESH_NODE_POSITIONS", 1);
219 
220  // build field
222 
223  // get HO geometry for 10 node tets
224  // This method takes coordinates form edges mid nodes in 10 node tet and
225  // project values on 2nd order hierarchical basis used to approx. geometry.
226  Projection10NodeCoordsOnField ent_method_material(mField,
227  "MESH_NODE_POSITIONS");
228  CHKERR mField.loop_dofs("MESH_NODE_POSITIONS", ent_method_material);
229 
231  }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:491
virtual MoFEMErrorCode add_field(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=-1)=0
Add field.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Common.hpp:147
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
int oRder
approximation order
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=-1)=0
Add entities to field meshsetThe lower dimension entities are added depending on the space type...
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=-1)=0
Set order approximation of the entities in the field.
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:141
field with continuous tangents
Definition: definitions.h:177
MoFEM::Interface & mField
#define CHKERR
Inline error check.
Definition: definitions.h:610
virtual MoFEMErrorCode build_fields(int verb=-1)=0
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, EntMethod &method, int lower_rank, int upper_rank, int verb=-1)=0
Make a loop over entities.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:436
continuous field
Definition: definitions.h:176
Construction of base is by Demkowicz .
Definition: definitions.h:144

◆ createProblem()

MoFEMErrorCode MagneticElement::createProblem ( )

create problemProblem is collection of finite elements. With the information on which fields finite elements operates the matrix and left and right hand side vector could be created.

Here we use Distributed mesh manager from PETSc as a inteface.

Returns
error code
Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 286 of file MagneticElement.hpp.

286  {
288  // set up DM
289  DMType dm_name = "MAGNETIC_PROBLEM";
290  CHKERR DMRegister_MoFEM(dm_name);
291  CHKERR DMCreate(PETSC_COMM_WORLD, &blockData.dM);
292  CHKERR DMSetType(blockData.dM, dm_name);
294  BitRefLevel().set(0));
295  CHKERR DMSetFromOptions(blockData.dM);
296  // add elements to blockData.dM
299  CHKERR DMSetUp(blockData.dM);
300  // create matrices and vectors
301  CHKERR DMCreateGlobalVector(blockData.dM, &blockData.D);
303  }
PetscErrorCode DMMoFEMAddElement(DM dm, const char fe_name[])
add element to dm
Definition: DMMMoFEM.cpp:449
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:491
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Common.hpp:147
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:91
MoFEM::Interface & mField
PetscErrorCode DMMoFEMCreateMoFEM(DM dm, MoFEM::Interface *m_field_ptr, const char problem_name[], const MoFEM::BitRefLevel bit_level, const MoFEM::BitRefLevel bit_mask=MoFEM::BitRefLevel().set())
Must be called by user to set MoFEM data structures.
Definition: DMMMoFEM.cpp:150
#define CHKERR
Inline error check.
Definition: definitions.h:610
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:436

◆ destroyProblem()

MoFEMErrorCode MagneticElement::destroyProblem ( )

destroy Distributed mesh manager

Returns
[description]
Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 310 of file MagneticElement.hpp.

310  {
311 
313  ierr = DMDestroy(&blockData.dM);
314  CHKERRG(ierr);
316  }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:515
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:558
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:521

◆ getEssentialBc()

MoFEMErrorCode MagneticElement::getEssentialBc ( )

get essential boundary conditions (only homogenous case is considered)

Returns
error code
Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 135 of file MagneticElement.hpp.

135  {
138  if (bit->getName().compare(0, 10, "ESSENTIALBC") == 0) {
139  Range faces;
140  CHKERR mField.get_moab().get_entities_by_type(bit->meshset, MBTRI,
141  faces, true);
142  CHKERR mField.get_moab().get_adjacencies(
143  faces, 1, true, blockData.essentialBc, moab::Interface::UNION);
144  blockData.essentialBc.merge(faces);
145  }
146  }
147  if (blockData.essentialBc.empty()) {
148  Range tets;
149  CHKERR mField.get_moab().get_entities_by_type(0, MBTET, tets);
150  Skinner skin(&mField.get_moab());
151  Range skin_faces; // skin faces from 3d ents
152  CHKERR skin.find_skin(0, tets, false, skin_faces);
153  skin_faces = subtract(skin_faces, blockData.naturalBc);
154  CHKERR mField.get_moab().get_adjacencies(
155  skin_faces, 1, true, blockData.essentialBc, moab::Interface::UNION);
156  blockData.essentialBc.merge(skin_faces);
157  }
159  }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:491
MoFEM::Interface & mField
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field...
#define CHKERR
Inline error check.
Definition: definitions.h:610
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:436

◆ getNaturalBc()

MoFEMErrorCode MagneticElement::getNaturalBc ( )

get natural boundary conditions

Returns
error code
Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 115 of file MagneticElement.hpp.

115  {
118  if (bit->getName().compare(0, 9, "NATURALBC") == 0) {
119  Range faces;
120  CHKERR mField.get_moab().get_entities_by_type(bit->meshset, MBTRI,
121  faces, true);
122  CHKERR mField.get_moab().get_adjacencies(
123  faces, 1, true, blockData.naturalBc, moab::Interface::UNION);
124  blockData.naturalBc.merge(faces);
125  }
126  }
128  }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:491
MoFEM::Interface & mField
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field...
#define CHKERR
Inline error check.
Definition: definitions.h:610
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:436

◆ postProcessResults()

MoFEMErrorCode MagneticElement::postProcessResults ( )

post-process results, i.e. save solution on the mesh

Returns
[description]
Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 414 of file MagneticElement.hpp.

414  {
415 
418  CHKERR post_proc.generateReferenceElementMesh();
419  CHKERR post_proc.addFieldValuesPostProc("MESH_NODE_POSITIONS");
420  CHKERR post_proc.addFieldValuesPostProc(blockData.fieldName);
421  post_proc.getOpPtrVector().push_back(new OpPostProcessCurl(
422  blockData, post_proc.postProcMesh, post_proc.mapGaussPts));
424  &post_proc);
425  CHKERR post_proc.writeFile("out_values.h5m");
427  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:491
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:531
Post processing.
MoFEM::Interface & mField
#define CHKERR
Inline error check.
Definition: definitions.h:610
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:436

◆ solveProblem()

MoFEMErrorCode MagneticElement::solveProblem ( const bool  regression_test = false)

solve problemCreate matrices; integrate over elements; solve linear system of equations

Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 324 of file MagneticElement.hpp.

324  {
326  CHKERR DMCreateMatrix(blockData.dM, &blockData.A);
327  CHKERR DMCreateGlobalVector(blockData.dM, &blockData.F);
328  CHKERR VecDuplicate(blockData.F, &blockData.D);
329 
330  VolumeFE vol_fe(mField);
331  vol_fe.getOpPtrVector().push_back(new OpCurlCurl(blockData));
332  vol_fe.getOpPtrVector().push_back(new OpStab(blockData));
333  TriFE tri_fe(mField);
334  tri_fe.getOpPtrVector().push_back(new OpNaturalBC(blockData));
335 
336  CHKERR MatZeroEntries(blockData.A);
337  CHKERR VecZeroEntries(blockData.F);
338  CHKERR VecGhostUpdateBegin(blockData.F, ADD_VALUES, SCATTER_FORWARD);
339  CHKERR VecGhostUpdateEnd(blockData.F, ADD_VALUES, SCATTER_FORWARD);
340 
342  &vol_fe);
344  blockData.feNaturalBCName.c_str(), &tri_fe);
345 
346  CHKERR MatAssemblyBegin(blockData.A, MAT_FINAL_ASSEMBLY);
347  CHKERR MatAssemblyEnd(blockData.A, MAT_FINAL_ASSEMBLY);
348  CHKERR VecGhostUpdateBegin(blockData.F, ADD_VALUES, SCATTER_REVERSE);
349  CHKERR VecGhostUpdateEnd(blockData.F, ADD_VALUES, SCATTER_REVERSE);
350  CHKERR VecAssemblyBegin(blockData.F);
351  CHKERR VecAssemblyEnd(blockData.F);
352 
353  // Boundary conditions
354  std::vector<int> dofs_bc_indices;
355  const MoFEM::Problem *problem_ptr;
356  ierr = DMMoFEMGetProblemPtr(blockData.dM, &problem_ptr);
357  for (Range::iterator eit = blockData.essentialBc.begin();
358  eit != blockData.essentialBc.end(); eit++) {
360  problem_ptr, blockData.fieldName, *eit, mField.get_comm_rank(),
361  dof_ptr)) {
362  std::bitset<8> pstatus(dof_ptr->get()->getPStatus());
363  if (pstatus.test(0))
364  continue; // only local
365  dofs_bc_indices.push_back(dof_ptr->get()->getPetscGlobalDofIdx());
366  }
367  }
368 
369  const double diag = 1;
370  CHKERR MatZeroRowsColumns(
371  blockData.A, dofs_bc_indices.size(),
372  dofs_bc_indices.empty() ? PETSC_NULL : &*dofs_bc_indices.begin(), diag,
373  PETSC_NULL, PETSC_NULL);
374 
375  KSP solver;
376  CHKERR KSPCreate(PETSC_COMM_WORLD, &solver);
377  // CHKERR KSPSetDM(solver,blockData.dM); CHKERRG(ierr);
378  CHKERR KSPSetOperators(solver, blockData.A, blockData.A);
379  CHKERR KSPSetFromOptions(solver);
380  CHKERR KSPSetUp(solver);
381  CHKERR KSPSolve(solver, blockData.F, blockData.D);
382  CHKERR VecGhostUpdateBegin(blockData.D, INSERT_VALUES, SCATTER_FORWARD);
383  CHKERR VecGhostUpdateEnd(blockData.D, INSERT_VALUES, SCATTER_FORWARD);
384  CHKERR KSPDestroy(&solver);
385  CHKERR VecDestroy(&blockData.F);
386  CHKERR MatDestroy(&blockData.A);
387 
388  CHKERR VecGhostUpdateBegin(blockData.D, INSERT_VALUES, SCATTER_FORWARD);
389  CHKERR VecGhostUpdateEnd(blockData.D, INSERT_VALUES, SCATTER_FORWARD);
391  SCATTER_REVERSE);
392 
393  if (regression_test) {
394  // This test is for order = 1 only
395  double nrm2;
396  CHKERR VecNorm(blockData.D, NORM_2, &nrm2);
397  const double expected_value = 4.6772e+01;
398  const double tol = 1e-4;
399  if ((nrm2 - expected_value) / expected_value > tol) {
400  SETERRQ2(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
401  "Regression test failed %6.4e != %6.4e", nrm2, expected_value);
402  }
403  }
404 
405  CHKERR VecDestroy(&blockData.D);
407  }
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:469
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:491
keeps basic data about problemThis is low level structure with information about problem, what elements compose problem and what DOFs are on rows and columns.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method)
Executes FEMethod for finite elements in DM.
Definition: DMMMoFEM.cpp:531
virtual int get_comm_rank() const =0
#define _IT_NUMEREDDOF_ROW_BY_NAME_ENT_PART_FOR_LOOP_(PROBLEMPTR, NAME, ENT, PART, IT)
use with loops to iterate row DOFs
MoFEM::Interface & mField
#define CHKERR
Inline error check.
Definition: definitions.h:610
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:436
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:372

Member Data Documentation

◆ blockData

BlockData MagneticElement::blockData
Examples:
MagneticElement.hpp, and magnetostatic.cpp.

Definition at line 108 of file MagneticElement.hpp.

◆ mField

MoFEM::Interface& MagneticElement::mField
Examples:
MagneticElement.hpp.

Definition at line 44 of file MagneticElement.hpp.


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