MFEM  v4.6.0
Finite element discretization library
tesla_solver.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_TESLA_SOLVER
13 #define MFEM_TESLA_SOLVER
14 
15 #include "../common/pfem_extras.hpp"
16 #include "../common/mesh_extras.hpp"
17 #include "electromagnetics.hpp"
18 
19 #ifdef MFEM_USE_MPI
20 
21 #include <string>
22 #include <map>
23 
24 namespace mfem
25 {
26 
27 using common::H1_ParFESpace;
28 using common::ND_ParFESpace;
29 using common::RT_ParFESpace;
30 using common::ParDiscreteGradOperator;
31 using common::ParDiscreteCurlOperator;
32 using common::DivergenceFreeProjector;
33 
34 namespace electromagnetics
35 {
36 
37 class SurfaceCurrent;
39 {
40 public:
41  TeslaSolver(ParMesh & pmesh, int order, Array<int> & kbcs,
42  Array<int> & vbcs, Vector & vbcv,
43  Coefficient & muInvCoef,
44  void (*a_bc )(const Vector&, Vector&),
45  void (*j_src)(const Vector&, Vector&),
46  void (*m_src)(const Vector&, Vector&));
47  ~TeslaSolver();
48 
50 
51  void PrintSizes();
52 
53  void Assemble();
54 
55  void Update();
56 
57  void Solve();
58 
59  void GetErrorEstimates(Vector & errors);
60 
62 
63  void WriteVisItFields(int it = 0);
64 
65  void InitializeGLVis();
66 
67  void DisplayToGLVis();
68 
69  const ParGridFunction & GetVectorPotential() { return *a_; }
70 
71 private:
72 
73  int myid_;
74  int num_procs_;
75  int order_;
76 
77  ParMesh * pmesh_;
78 
79  VisItDataCollection * visit_dc_;
80 
81  H1_ParFESpace * H1FESpace_;
82  ND_ParFESpace * HCurlFESpace_;
83  RT_ParFESpace * HDivFESpace_;
84 
85  ParBilinearForm * curlMuInvCurl_;
86  ParBilinearForm * hCurlMass_;
87  ParMixedBilinearForm * hDivHCurlMuInv_;
88  ParMixedBilinearForm * weakCurlMuInv_;
89 
92 
93  ParGridFunction * a_; // Vector Potential (HCurl)
94  ParGridFunction * b_; // Magnetic Flux (HDiv)
95  ParGridFunction * h_; // Magnetic Field (HCurl)
96  ParGridFunction * jr_; // Raw Volumetric Current Density (HCurl)
97  ParGridFunction * j_; // Volumetric Current Density (HCurl)
98  ParGridFunction * k_; // Surface Current Density (HCurl)
99  ParGridFunction * m_; // Magnetization (HDiv)
100  ParGridFunction * bd_; // Dual of B (HCurl)
101  ParGridFunction * jd_; // Dual of J, the rhs vector (HCurl)
102 
103  DivergenceFreeProjector * DivFreeProj_;
104  SurfaceCurrent * SurfCur_;
105 
106  Coefficient * muInvCoef_; // Dia/Paramagnetic Material Coefficient
107  VectorCoefficient * aBCCoef_; // Vector Potential BC Function
108  VectorCoefficient * jCoef_; // Volume Current Density Function
109  VectorCoefficient * mCoef_; // Magnetization Vector Function
110 
111  void (*a_bc_ )(const Vector&, Vector&);
112  void (*j_src_)(const Vector&, Vector&);
113  void (*m_src_)(const Vector&, Vector&);
114 
115  Array<int> ess_bdr_;
116  Array<int> ess_bdr_tdofs_;
117  Array<int> non_k_bdr_;
118 
119  std::map<std::string,socketstream*> socks_;
120 };
121 
123 {
124 public:
127  Array<int> & kbcs, Array<int> & vbcs, Vector & vbcv);
128  ~SurfaceCurrent();
129 
130  void InitSolver() const;
131 
133 
134  void Update();
135 
136  ParGridFunction * GetPsi() { return psi_; }
137 
138 private:
139  int myid_;
140 
141  ParFiniteElementSpace * H1FESpace_;
142  ParDiscreteGradOperator * grad_;
143  Array<int> * kbcs_;
144  Array<int> * vbcs_;
145  Vector * vbcv_;
146 
147  ParBilinearForm * s0_;
148  ParGridFunction * psi_;
149  ParGridFunction * rhs_;
150 
151  HypreParMatrix * S0_;
152  mutable Vector Psi_;
153  mutable Vector RHS_;
154 
155  mutable HypreBoomerAMG * amg_;
156  mutable HyprePCG * pcg_;
157 
158  Array<int> ess_bdr_, ess_bdr_tdofs_;
159  Array<int> non_k_bdr_;
160 };
161 
162 } // namespace electromagnetics
163 
164 } // namespace mfem
165 
166 #endif // MFEM_USE_MPI
167 
168 #endif // MFEM_TESLA_SOLVER
Base class for vector Coefficients that optionally depend on time and space.
Abstract parallel finite element space.
Definition: pfespace.hpp:28
void GetErrorEstimates(Vector &errors)
The BoomerAMG solver in hypre.
Definition: hypre.hpp:1590
SurfaceCurrent(ParFiniteElementSpace &H1FESpace, ParDiscreteGradOperator &Grad, Array< int > &kbcs, Array< int > &vbcs, Vector &vbcv)
Data collection with VisIt I/O routines.
void RegisterVisItFields(VisItDataCollection &visit_dc)
PCG solver in hypre.
Definition: hypre.hpp:1215
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:41
HYPRE_Int HYPRE_BigInt
Class for parallel bilinear form using different test and trial FE spaces.
void j_src(const Vector &x, double t, Vector &j)
Definition: maxwell.cpp:108
const ParGridFunction & GetVectorPotential()
Class for parallel bilinear form.
void ComputeSurfaceCurrent(ParGridFunction &k)
Vector data type.
Definition: vector.hpp:58
Class for parallel grid function.
Definition: pgridfunc.hpp:32
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:343
Class for parallel meshes.
Definition: pmesh.hpp:32
TeslaSolver(ParMesh &pmesh, int order, Array< int > &kbcs, Array< int > &vbcs, Vector &vbcv, Coefficient &muInvCoef, void(*a_bc)(const Vector &, Vector &), void(*j_src)(const Vector &, Vector &), void(*m_src)(const Vector &, Vector &))