MFEM  v3.3.2
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
tesla_solver.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
11 
12 #ifndef MFEM_TESLA_SOLVER
13 #define MFEM_TESLA_SOLVER
14 
15 #include "../common/pfem_extras.hpp"
16 
17 #ifdef MFEM_USE_MPI
18 
19 #include <string>
20 #include <map>
21 
22 namespace mfem
23 {
24 
25 using miniapps::H1_ParFESpace;
26 using miniapps::ND_ParFESpace;
27 using miniapps::RT_ParFESpace;
28 using miniapps::ParDiscreteGradOperator;
29 using miniapps::ParDiscreteCurlOperator;
30 using miniapps::DivergenceFreeProjector;
31 
32 namespace electromagnetics
33 {
34 
35 // Physical Constants
36 // Permeability of Free Space (units H/m)
37 static double mu0_ = 4.0e-7*M_PI;
38 
39 class SurfaceCurrent;
41 {
42 public:
43  TeslaSolver(ParMesh & pmesh, int order, Array<int> & kbcs,
44  Array<int> & vbcs, Vector & vbcv,
45  Coefficient & muInvCoef,
46  void (*a_bc )(const Vector&, Vector&),
47  void (*j_src)(const Vector&, Vector&),
48  void (*m_src)(const Vector&, Vector&));
49  ~TeslaSolver();
50 
51  HYPRE_Int GetProblemSize();
52 
53  void PrintSizes();
54 
55  void Assemble();
56 
57  void Update();
58 
59  void Solve();
60 
61  void GetErrorEstimates(Vector & errors);
62 
64 
65  void WriteVisItFields(int it = 0);
66 
67  void InitializeGLVis();
68 
69  void DisplayToGLVis();
70 
71  const ParGridFunction & GetVectorPotential() { return *a_; }
72 
73 private:
74 
75  int myid_;
76  int num_procs_;
77  int order_;
78 
79  ParMesh * pmesh_;
80 
81  VisItDataCollection * visit_dc_;
82 
83  H1_ParFESpace * H1FESpace_;
84  ND_ParFESpace * HCurlFESpace_;
85  RT_ParFESpace * HDivFESpace_;
86 
87  ParBilinearForm * curlMuInvCurl_;
88  ParBilinearForm * hCurlMass_;
89  ParMixedBilinearForm * hDivHCurlMuInv_;
90  ParMixedBilinearForm * weakCurlMuInv_;
91 
94 
95  ParGridFunction * a_; // Vector Potential (HCurl)
96  ParGridFunction * b_; // Magnetic Flux (HDiv)
97  ParGridFunction * h_; // Magnetic Field (HCurl)
98  ParGridFunction * jr_; // Raw Volumetric Current Density (HCurl)
99  ParGridFunction * j_; // Volumetric Current Density (HCurl)
100  ParGridFunction * k_; // Surface Current Density (HCurl)
101  ParGridFunction * m_; // Magnetization (HDiv)
102  ParGridFunction * bd_; // Dual of B (HCurl)
103  ParGridFunction * jd_; // Dual of J, the rhs vector (HCurl)
104 
105  DivergenceFreeProjector * DivFreeProj_;
106  SurfaceCurrent * SurfCur_;
107 
108  Coefficient * muInvCoef_; // Dia/Paramagnetic Material Coefficient
109  VectorCoefficient * aBCCoef_; // Vector Potential BC Function
110  VectorCoefficient * jCoef_; // Volume Current Density Function
111  VectorCoefficient * mCoef_; // Magnetization Vector Function
112 
113  void (*a_bc_ )(const Vector&, Vector&);
114  void (*j_src_)(const Vector&, Vector&);
115  void (*m_src_)(const Vector&, Vector&);
116 
117  Array<int> ess_bdr_;
118  Array<int> ess_bdr_tdofs_;
119  Array<int> non_k_bdr_;
120 
121  std::map<std::string,socketstream*> socks_;
122 };
123 
125 {
126 public:
129  Array<int> & kbcs, Array<int> & vbcs, Vector & vbcv);
130  ~SurfaceCurrent();
131 
132  void InitSolver() const;
133 
135 
136  void Update();
137 
138  ParGridFunction * GetPsi() { return psi_; }
139 
140 private:
141  int myid_;
142 
143  ParFiniteElementSpace * H1FESpace_;
144  ParDiscreteGradOperator * grad_;
145  Array<int> * kbcs_;
146  Array<int> * vbcs_;
147  Vector * vbcv_;
148 
149  ParBilinearForm * s0_;
150  ParGridFunction * psi_;
151  ParGridFunction * rhs_;
152 
153  HypreParMatrix * S0_;
154  mutable Vector Psi_;
155  mutable Vector RHS_;
156 
157  mutable HypreBoomerAMG * amg_;
158  mutable HyprePCG * pcg_;
159 
160  Array<int> ess_bdr_, ess_bdr_tdofs_;
161  Array<int> non_k_bdr_;
162 };
163 
164 } // namespace electromagnetics
165 
166 } // namespace mfem
167 
168 #endif // MFEM_USE_MPI
169 
170 #endif // MFEM_TESLA_SOLVER
Abstract parallel finite element space.
Definition: pfespace.hpp:31
void GetErrorEstimates(Vector &errors)
The BoomerAMG solver in hypre.
Definition: hypre.hpp:783
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:643
Base class Coefficient that may optionally depend on time.
Definition: coefficient.hpp:31
Class for parallel bilinear form using different test and trial FE spaces.
const ParGridFunction & GetVectorPotential()
Class for parallel bilinear form.
void ComputeSurfaceCurrent(ParGridFunction &k)
Vector data type.
Definition: vector.hpp:41
Class for parallel grid function.
Definition: pgridfunc.hpp:32
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:175
Class for parallel meshes.
Definition: pmesh.hpp:29
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 &))