MFEM  v4.4.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
convergence.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, 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_CONVERGENCE
13 #define MFEM_CONVERGENCE
14 
15 #include "../linalg/linalg.hpp"
16 #include "gridfunc.hpp"
17 #ifdef MFEM_USE_MPI
18 #include "pgridfunc.hpp"
19 #endif
20 
21 namespace mfem
22 {
23 
24 /** @brief Class to compute error and convergence rates.
25  It supports H1, H(curl) (ND elements), H(div) (RT elements) and L2 (DG).
26 
27  For "smooth enough" solutions the Galerkin error measured in the appropriate
28  norm satisfies || u - u_h || ~ h^k
29 
30  Here, k is called the asymptotic rate of convergence
31 
32  For successive h-refinements the rate can be estimated by
33  k = log(||u - u_h|| / ||u - u_{h/2}||)/(1/dim * log(N_{h/2}/N_{h})
34 */
36 {
37 private:
38  // counters for solutions/derivatives
39  int counter=0;
40  int dcounter=0;
41  int fcounter=0;
42 
43  // space continuity type
44  int cont_type=-1;
45 
46  // printing flag for helpful for MPI calls
47  int print_flag=1;
48 
49  // exact solution and derivatives
50  double CoeffNorm;
51  double CoeffDNorm;
52 
53  // Arrays to store error/rates
54  Array<double> L2Errors, DGFaceErrors, DErrors, EnErrors;
55  Array<double> L2Rates, DGFaceRates, DRates, EnRates;
56  Array<int> ndofs;
57 
58  void AddL2Error(GridFunction *gf, Coefficient *scalar_u,
59  VectorCoefficient *vector_u);
60  void AddGf(GridFunction *gf, Coefficient *scalar_u,
61  VectorCoefficient *grad=nullptr,
62  Coefficient *ell_coeff=nullptr,
63  JumpScaling jump_scaling = {1.0, JumpScaling::ONE_OVER_H});
64  void AddGf(GridFunction *gf, VectorCoefficient *vector_u,
65  VectorCoefficient *curl, Coefficient *div);
66  // returns the L2-norm of scalar_u or vector_u
67  double GetNorm(GridFunction *gf, Coefficient *scalar_u,
68  VectorCoefficient *vector_u);
69 
70 public:
71 
72  /// Clear any internal data
73  void Reset();
74 
75  /// Add L2 GridFunction, the exact solution and possibly its gradient and/or
76  /// DG face jumps parameters
78  VectorCoefficient *grad=nullptr,
79  Coefficient *ell_coeff=nullptr,
80  JumpScaling jump_scaling = {1.0, JumpScaling::ONE_OVER_H})
81  {
82  AddGf(gf, scalar_u, grad, ell_coeff, jump_scaling);
83  }
84 
85  /// Add H1 GridFunction, the exact solution and possibly its gradient
87  VectorCoefficient *grad=nullptr)
88  {
89  AddGf(gf, scalar_u, grad);
90  }
91 
92  /// Add H(curl) GridFunction, the exact solution and possibly its curl
94  VectorCoefficient *curl=nullptr)
95  {
96  AddGf(gf, vector_u, curl, nullptr);
97  }
98 
99  /// Add H(div) GridFunction, the exact solution and possibly its div
101  Coefficient *div=nullptr)
102  {
103  AddGf(gf,vector_u, nullptr, div);
104  }
105 
106  /// Get the L2 error at step n
107  double GetL2Error(int n)
108  {
109  MFEM_VERIFY( n <= counter,"Step out of bounds")
110  return L2Errors[n];
111  }
112 
113  /// Get all L2 errors
114  void GetL2Errors(Array<double> & L2Errors_)
115  {
116  L2Errors_ = L2Errors;
117  }
118 
119  /// Get the Grad/Curl/Div error at step n
120  double GetDError(int n)
121  {
122  MFEM_VERIFY(n <= dcounter,"Step out of bounds")
123  return DErrors[n];
124  }
125 
126  /// Get all Grad/Curl/Div errors
127  void GetDErrors(Array<double> & DErrors_)
128  {
129  DErrors_ = DErrors;
130  }
131 
132  /// Get the DGFaceJumps error at step n
133  double GetDGFaceJumpsError(int n)
134  {
135  MFEM_VERIFY(n<= fcounter,"Step out of bounds")
136  return DGFaceErrors[n];
137  }
138 
139  /// Get all DGFaceJumps errors
140  void GetDGFaceJumpsErrors(Array<double> & DGFaceErrors_)
141  {
142  DGFaceErrors_ = DGFaceErrors;
143  }
144 
145  /// Print rates and errors
146  void Print(bool relative = false, std::ostream &out = mfem::out);
147 };
148 
149 } // namespace mfem
150 
151 #endif // MFEM_CONVERGENCE
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
void Reset()
Clear any internal data.
Definition: convergence.cpp:9
Base class for vector Coefficients that optionally depend on time and space.
void AddHcurlGridFunction(GridFunction *gf, VectorCoefficient *vector_u, VectorCoefficient *curl=nullptr)
Add H(curl) GridFunction, the exact solution and possibly its curl.
Definition: convergence.hpp:93
void AddL2GridFunction(GridFunction *gf, Coefficient *scalar_u, VectorCoefficient *grad=nullptr, Coefficient *ell_coeff=nullptr, JumpScaling jump_scaling={1.0, JumpScaling::ONE_OVER_H})
Definition: convergence.hpp:77
void AddHdivGridFunction(GridFunction *gf, VectorCoefficient *vector_u, Coefficient *div=nullptr)
Add H(div) GridFunction, the exact solution and possibly its div.
void GetDGFaceJumpsErrors(Array< double > &DGFaceErrors_)
Get all DGFaceJumps errors.
double GetL2Error(int n)
Get the L2 error at step n.
void GetDErrors(Array< double > &DErrors_)
Get all Grad/Curl/Div errors.
double GetDGFaceJumpsError(int n)
Get the DGFaceJumps error at step n.
void GetL2Errors(Array< double > &L2Errors_)
Get all L2 errors.
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:39
void Print(bool relative=false, std::ostream &out=mfem::out)
Print rates and errors.
Class to compute error and convergence rates. It supports H1, H(curl) (ND elements), H(div) (RT elements) and L2 (DG).
Definition: convergence.hpp:35
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:66
void AddH1GridFunction(GridFunction *gf, Coefficient *scalar_u, VectorCoefficient *grad=nullptr)
Add H1 GridFunction, the exact solution and possibly its gradient.
Definition: convergence.hpp:86
double GetDError(int n)
Get the Grad/Curl/Div error at step n.