MFEM  v3.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
volta_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_VOLTA_SOLVER
13 #define MFEM_VOLTA_SOLVER
14 
15 #include "../../config/config.hpp"
16 
17 #ifdef MFEM_USE_MPI
18 
19 #include "../common/pfem_extras.hpp"
20 #include <string>
21 #include <map>
22 
23 namespace mfem
24 {
25 
26 using miniapps::H1_ParFESpace;
27 using miniapps::ND_ParFESpace;
28 using miniapps::RT_ParFESpace;
29 using miniapps::ParDiscreteGradOperator;
30 
31 namespace electromagnetics
32 {
33 
34 // Physical Constants
35 // Permittivity of Free Space (units F/m)
36 static double epsilon0_ = 8.8541878176e-12;
37 
39 {
40 public:
41  VoltaSolver(ParMesh & pmesh, int order,
42  Array<int> & dbcs, Vector & dbcv,
43  Array<int> & nbcs, Vector & nbcv,
44  double (*eps )(const Vector&),
45  double (*phi_bc )(const Vector&),
46  double (*rho_src)(const Vector&),
47  void (*p_src )(const Vector&, Vector&));
48  ~VoltaSolver();
49 
50  HYPRE_Int GetProblemSize();
51 
52  void PrintSizes();
53 
54  void Update();
55 
56  void Solve();
57 
58  void GetErrorEstimates(Vector & errors);
59 
61 
62  void WriteVisItFields(int it = 0);
63 
64  void InitializeGLVis();
65 
66  void DisplayToGLVis();
67 
68  const ParGridFunction & GetVectorPotential() { return *phi_; }
69 
70 private:
71 
72  int myid_;
73  int num_procs_;
74  int order_;
75 
76  ParMesh * pmesh_;
77 
78  Array<int> * dbcs_;
79  Vector * dbcv_;
80  Array<int> * nbcs_;
81  Vector * nbcv_;
82 
83  VisItDataCollection * visit_dc_;
84 
85  H1_ParFESpace * H1FESpace_;
86  ND_ParFESpace * HCurlFESpace_;
87  RT_ParFESpace * HDivFESpace_;
88 
89  ParBilinearForm * divEpsGrad_;
90  ParBilinearForm * h1Mass_;
91  ParBilinearForm * h1SurfMass_;
92  ParBilinearForm * hCurlMass_;
93  ParBilinearForm * hDivMass_;
94  ParMixedBilinearForm * hCurlHDivEps_;
95  ParMixedBilinearForm * hCurlHDiv_;
96 
98 
99  ParGridFunction * phi_;
100  ParGridFunction * rho_;
101  ParGridFunction * sigma_;
102  ParGridFunction * e_;
103  ParGridFunction * d_;
104  ParGridFunction * p_;
105 
106  Coefficient * epsCoef_;
107  Coefficient * phiBCCoef_;
108  Coefficient * rhoCoef_;
109  VectorCoefficient * pCoef_;
110 
111  double (*eps_ )(const Vector&);
112  double (*phi_bc_ )(const Vector&);
113  double (*rho_src_)(const Vector&);
114  void (*p_src_ )(const Vector&, Vector&);
115 
116  std::map<std::string,socketstream*> socks_;
117 
118  Array<int> ess_bdr_;
119 };
120 
121 } // namespace electromagnetics
122 
123 } // namespace mfem
124 
125 #endif // MFEM_USE_MPI
126 
127 #endif // MFEM_VOLTA_SOLVER
void RegisterVisItFields(VisItDataCollection &visit_dc)
void GetErrorEstimates(Vector &errors)
Data collection with VisIt I/O routines.
Base class Coefficient that may optionally depend on time.
Definition: coefficient.hpp:31
Class for parallel bilinear form.
const ParGridFunction & GetVectorPotential()
Class for parallel bilinear form.
Vector data type.
Definition: vector.hpp:33
Class for parallel grid function.
Definition: pgridfunc.hpp:31
Class for parallel meshes.
Definition: pmesh.hpp:28
VoltaSolver(ParMesh &pmesh, int order, Array< int > &dbcs, Vector &dbcv, Array< int > &nbcs, Vector &nbcv, double(*eps)(const Vector &), double(*phi_bc)(const Vector &), double(*rho_src)(const Vector &), void(*p_src)(const Vector &, Vector &))