error_estimator.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 // Local include files
20 #include "libmesh/libmesh_common.h"
22 #include "libmesh/error_vector.h"
24 #include "libmesh/parallel.h"
26 #include "libmesh/int_range.h"
27 
28 namespace libMesh
29 {
30 
31 // ErrorEstimator functions
32 void ErrorEstimator::reduce_error (std::vector<ErrorVectorReal> & error_per_cell,
33  const Parallel::Communicator & comm) const
34 {
35  // This function must be run on all processors at once
36  // parallel_object_only();
37 
38  // Each processor has now computed the error contributions
39  // for its local elements. We may need to sum the vector to
40  // recover the error for each element.
41 
42  comm.sum(error_per_cell);
43 }
44 
45 
46 
47 void ErrorEstimator::estimate_errors(const EquationSystems & equation_systems,
48  ErrorVector & error_per_cell,
49  const std::map<const System *, SystemNorm> & error_norms,
50  const std::map<const System *, const NumericVector<Number> *> * solution_vectors,
51  bool estimate_parent_error)
52 {
53  SystemNorm old_error_norm = this->error_norm;
54 
55  // Sum the error values from each system
56  for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
57  {
58  ErrorVector system_error_per_cell;
59  const System & sys = equation_systems.get_system(s);
60  if (error_norms.find(&sys) == error_norms.end())
61  this->error_norm = old_error_norm;
62  else
63  this->error_norm = error_norms.find(&sys)->second;
64 
65  const NumericVector<Number> * solution_vector = nullptr;
66  if (solution_vectors &&
67  solution_vectors->find(&sys) != solution_vectors->end())
68  solution_vector = solution_vectors->find(&sys)->second;
69 
70  this->estimate_error(sys, system_error_per_cell,
71  solution_vector, estimate_parent_error);
72 
73  if (s)
74  {
75  libmesh_assert_equal_to (error_per_cell.size(), system_error_per_cell.size());
76  for (auto i : index_range(error_per_cell))
77  error_per_cell[i] += system_error_per_cell[i];
78  }
79  else
80  error_per_cell = system_error_per_cell;
81  }
82 
83  // Restore our old state before returning
84  this->error_norm = old_error_norm;
85 }
86 
87 
88 
93 void ErrorEstimator::estimate_errors(const EquationSystems & equation_systems,
94  ErrorMap & errors_per_cell,
95  const std::map<const System *, const NumericVector<Number> *> * solution_vectors,
96  bool estimate_parent_error)
97 {
98  SystemNorm old_error_norm = this->error_norm;
99 
100  // Find the requested error values from each system
101  for (unsigned int s = 0; s != equation_systems.n_systems(); ++s)
102  {
103  const System & sys = equation_systems.get_system(s);
104 
105  unsigned int n_vars = sys.n_vars();
106 
107  for (unsigned int v = 0; v != n_vars; ++v)
108  {
109  // Only fill in ErrorVectors the user asks for
110  if (errors_per_cell.find(std::make_pair(&sys, v)) ==
111  errors_per_cell.end())
112  continue;
113 
114  // Calculate error in only one variable
115  std::vector<Real> weights(n_vars, 0.0);
116  weights[v] = 1.0;
117  this->error_norm =
118  SystemNorm(std::vector<FEMNormType>(n_vars, old_error_norm.type(v)),
119  weights);
120 
121  const NumericVector<Number> * solution_vector = nullptr;
122  if (solution_vectors &&
123  solution_vectors->find(&sys) != solution_vectors->end())
124  solution_vector = solution_vectors->find(&sys)->second;
125 
126  this->estimate_error
127  (sys, *errors_per_cell[std::make_pair(&sys, v)],
128  solution_vector, estimate_parent_error);
129  }
130  }
131 
132  // Restore our old state before returning
133  this->error_norm = old_error_norm;
134 }
135 
136 } // namespace libMesh
virtual void estimate_errors(const EquationSystems &equation_systems, ErrorVector &error_per_cell, const std::map< const System *, SystemNorm > &error_norms, const std::map< const System *, const NumericVector< Number > *> *solution_vectors=nullptr, bool estimate_parent_error=false)
Manages multiples systems of equations.
unsigned int n_systems() const
const T_sys & get_system(const std::string &name) const
IntRange< std::size_t > index_range(const std::vector< T > &vec)
Definition: int_range.h:104
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0
const unsigned int n_vars
Definition: tecplot_io.C:69
Manages consistently variables, degrees of freedom, and coefficient vectors.
Definition: system.h:92
FEMNormType type(unsigned int var) const
Definition: system_norm.C:110
std::map< std::pair< const System *, unsigned int >, ErrorVector * > ErrorMap
void reduce_error(std::vector< ErrorVectorReal > &error_per_cell, const Parallel::Communicator &comm) const
unsigned int n_vars() const
Definition: system.h:2105