libMesh
Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
libMesh::AdjointResidualErrorEstimator Class Reference

This class implements a goal oriented error indicator, by weighting residual-based estimates from the primal problem against estimates from the adjoint problem. More...

#include <adjoint_residual_error_estimator.h>

Inheritance diagram for libMesh::AdjointResidualErrorEstimator:
[legend]

Public Types

typedef std::map< std::pair< const System *, unsigned int >, std::unique_ptr< ErrorVector > > ErrorMap
 When calculating many error vectors at once, we need a data structure to hold them all. More...
 

Public Member Functions

 AdjointResidualErrorEstimator ()
 Constructor. More...
 
 AdjointResidualErrorEstimator (const AdjointResidualErrorEstimator &)=delete
 This class cannot be (default) copy constructed/assigned because it has unique_ptr members. More...
 
AdjointResidualErrorEstimatoroperator= (const AdjointResidualErrorEstimator &)=delete
 
 AdjointResidualErrorEstimator (AdjointResidualErrorEstimator &&)=default
 Defaulted move ctor, move assignment operator, and destructor. More...
 
AdjointResidualErrorEstimatoroperator= (AdjointResidualErrorEstimator &&)=default
 
virtual ~AdjointResidualErrorEstimator ()=default
 
std::unique_ptr< ErrorEstimator > & primal_error_estimator ()
 Access to the "subestimator" (default: PatchRecovery) to use on the primal/forward solution. More...
 
std::unique_ptr< ErrorEstimator > & dual_error_estimator ()
 Access to the "subestimator" (default: PatchRecovery) to use on the dual/adjoint solution. More...
 
QoISetqoi_set ()
 Access to the QoISet (default: weight all QoIs equally) to use when computing errors. More...
 
const QoISetqoi_set () const
 Access to the QoISet (default: weight all QoIs equally) to use when computing errors. More...
 
virtual void estimate_error (const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false) override
 Compute the adjoint-weighted error on each element and place it in the error_per_cell vector. More...
 
virtual ErrorEstimatorType type () const override
 
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)
 This virtual function can be redefined in derived classes, but by default computes the sum of the error_per_cell for each system in the equation_systems. More...
 
virtual void estimate_errors (const EquationSystems &equation_systems, ErrorMap &errors_per_cell, const std::map< const System *, const NumericVector< Number > *> *solution_vectors=nullptr, bool estimate_parent_error=false)
 This virtual function can be redefined in derived classes, but by default it calls estimate_error repeatedly to calculate the requested error vectors. More...
 

Public Attributes

std::string error_plot_suffix
 To aid in investigating error estimator behavior, set this string to a suffix with which to plot (prefixed by "primal_" or "dual_") the subestimator results. More...
 
SystemNorm error_norm
 When estimating the error in a single system, the error_norm is used to control the scaling and norm choice for each variable. More...
 

Protected Member Functions

void reduce_error (std::vector< ErrorVectorReal > &error_per_cell, const Parallel::Communicator &comm) const
 This method takes the local error contributions in error_per_cell from each processor and combines them to get the global error vector. More...
 

Protected Attributes

std::unique_ptr< ErrorEstimator_primal_error_estimator
 An error estimator for the forward problem. More...
 
std::unique_ptr< ErrorEstimator_dual_error_estimator
 An error estimator for the adjoint problem. More...
 
QoISet _qoi_set
 A QoISet to handle cases with multiple QoIs available. More...
 

Detailed Description

This class implements a goal oriented error indicator, by weighting residual-based estimates from the primal problem against estimates from the adjoint problem.

This is based on a trick suggested by Brian Carnes, (first proposed by Babuska and Miller in 1984) but if it doesn't actually work then the misunderstanding or misimplementation will be the fault of Roy Stogner. It's also Roy's fault there's no literature reference here yet.

Author
Roy H. Stogner
Date
2009

Definition at line 49 of file adjoint_residual_error_estimator.h.

Member Typedef Documentation

◆ ErrorMap

typedef std::map<std::pair<const System *, unsigned int>, std::unique_ptr<ErrorVector> > libMesh::ErrorEstimator::ErrorMap
inherited

When calculating many error vectors at once, we need a data structure to hold them all.

Definition at line 121 of file error_estimator.h.

Constructor & Destructor Documentation

◆ AdjointResidualErrorEstimator() [1/3]

libMesh::AdjointResidualErrorEstimator::AdjointResidualErrorEstimator ( )

Constructor.

Responsible for picking default subestimators.

Definition at line 42 of file adjoint_residual_error_estimator.C.

42  :
45  _primal_error_estimator(std::make_unique<PatchRecoveryErrorEstimator>()),
46  _dual_error_estimator(std::make_unique<PatchRecoveryErrorEstimator>()),
47  _qoi_set(QoISet())
48 {
49 }
QoISet _qoi_set
A QoISet to handle cases with multiple QoIs available.
std::string error_plot_suffix
To aid in investigating error estimator behavior, set this string to a suffix with which to plot (pre...
ErrorEstimator()=default
Constructor.
std::unique_ptr< ErrorEstimator > _dual_error_estimator
An error estimator for the adjoint problem.
std::unique_ptr< ErrorEstimator > _primal_error_estimator
An error estimator for the forward problem.

◆ AdjointResidualErrorEstimator() [2/3]

libMesh::AdjointResidualErrorEstimator::AdjointResidualErrorEstimator ( const AdjointResidualErrorEstimator )
delete

This class cannot be (default) copy constructed/assigned because it has unique_ptr members.

Explicitly deleting these functions is the best way to document this fact.

◆ AdjointResidualErrorEstimator() [3/3]

libMesh::AdjointResidualErrorEstimator::AdjointResidualErrorEstimator ( AdjointResidualErrorEstimator &&  )
default

Defaulted move ctor, move assignment operator, and destructor.

◆ ~AdjointResidualErrorEstimator()

virtual libMesh::AdjointResidualErrorEstimator::~AdjointResidualErrorEstimator ( )
virtualdefault

Member Function Documentation

◆ dual_error_estimator()

std::unique_ptr<ErrorEstimator>& libMesh::AdjointResidualErrorEstimator::dual_error_estimator ( )
inline

Access to the "subestimator" (default: PatchRecovery) to use on the dual/adjoint solution.

Definition at line 83 of file adjoint_residual_error_estimator.h.

References _dual_error_estimator.

83 { return _dual_error_estimator; }
std::unique_ptr< ErrorEstimator > _dual_error_estimator
An error estimator for the adjoint problem.

◆ estimate_error()

void libMesh::AdjointResidualErrorEstimator::estimate_error ( const System system,
ErrorVector error_per_cell,
const NumericVector< Number > *  solution_vector = nullptr,
bool  estimate_parent_error = false 
)
overridevirtual

Compute the adjoint-weighted error on each element and place it in the error_per_cell vector.

Note
this->error_norm is ignored; the error estimate is in the seminorm given by the absolute value of the error in the quantity of interest functional. The primal and dual subestimator error_norm values are used, and should be chosen appropriately for your model.

Implements libMesh::ErrorEstimator.

Definition at line 61 of file adjoint_residual_error_estimator.C.

References _dual_error_estimator, _primal_error_estimator, _qoi_set, libMesh::System::adjoint_solve(), libMesh::SystemNorm::calculate_norm(), libMesh::ErrorEstimator::error_norm, error_plot_suffix, libMesh::ErrorVectorReal, libMesh::System::get_adjoint_solution(), libMesh::System::get_equation_systems(), libMesh::System::get_mesh(), libMesh::QoISet::has_index(), libMesh::index_range(), libMesh::System::is_adjoint_already_solved(), libMesh::SystemNorm::is_identity(), libMesh::libmesh_assert(), mesh, libMesh::System::n_qois(), n_vars, libMesh::System::n_vars(), libMesh::ErrorVector::plot_error(), libMesh::Real, libMesh::System::solution, and libMesh::QoISet::weight().

65 {
66  LOG_SCOPE("estimate_error()", "AdjointResidualErrorEstimator");
67 
68  // The current mesh
69  const MeshBase & mesh = _system.get_mesh();
70 
71  // Resize the error_per_cell vector to be
72  // the number of elements, initialize it to 0.
73  error_per_cell.resize (mesh.max_elem_id());
74  std::fill (error_per_cell.begin(), error_per_cell.end(), 0.);
75 
76  // Get the number of variables in the system
77  unsigned int n_vars = _system.n_vars();
78 
79  // We need to make a map of the pointer to the solution vector
80  std::map<const System *, const NumericVector<Number> *>solutionvecs;
81  solutionvecs[&_system] = _system.solution.get();
82 
83  // Solve the dual problem if we have to
84  if (!_system.is_adjoint_already_solved())
85  {
86  // FIXME - we'll need to change a lot of APIs to make this trick
87  // work with a const System...
88  System & system = const_cast<System &>(_system);
89  system.adjoint_solve(_qoi_set);
90  }
91 
92  // Flag to check whether we have not been asked to weight the variable error contributions in any specific manner
93  bool error_norm_is_identity = error_norm.is_identity();
94 
95  // Create an ErrorMap/ErrorVector to store the primal, dual and total_dual variable errors
96  ErrorMap primal_errors_per_cell;
97  ErrorMap dual_errors_per_cell;
98  ErrorMap total_dual_errors_per_cell;
99  // Allocate ErrorVectors to this map if we're going to use it
100  if (!error_norm_is_identity)
101  for (unsigned int v = 0; v < n_vars; v++)
102  {
103  primal_errors_per_cell[std::make_pair(&_system, v)] = std::make_unique<ErrorVector>();
104  dual_errors_per_cell[std::make_pair(&_system, v)] = std::make_unique<ErrorVector>();
105  total_dual_errors_per_cell[std::make_pair(&_system, v)] = std::make_unique<ErrorVector>();
106  }
107  ErrorVector primal_error_per_cell;
108  ErrorVector dual_error_per_cell;
109  ErrorVector total_dual_error_per_cell;
110 
111  // Have we been asked to weight the variable error contributions in any specific manner
112  if (!error_norm_is_identity) // If we do
113  {
114  // Estimate the primal problem error for each variable
115  _primal_error_estimator->estimate_errors
116  (_system.get_equation_systems(), primal_errors_per_cell, &solutionvecs, estimate_parent_error);
117  }
118  else // If not
119  {
120  // Just get the combined error estimate
121  _primal_error_estimator->estimate_error
122  (_system, primal_error_per_cell, solution_vector, estimate_parent_error);
123  }
124 
125  // Sum and weight the dual error estimate based on our QoISet
126  for (unsigned int i = 0, n_qois = _system.n_qois(); i != n_qois; ++i)
127  {
128  if (_qoi_set.has_index(i))
129  {
130  // Get the weight for the current QoI
131  Real error_weight = _qoi_set.weight(i);
132 
133  // We need to make a map of the pointer to the adjoint solution vector
134  std::map<const System *, const NumericVector<Number> *>adjointsolutionvecs;
135  adjointsolutionvecs[&_system] = &_system.get_adjoint_solution(i);
136 
137  // Have we been asked to weight the variable error contributions in any specific manner
138  if (!error_norm_is_identity) // If we have
139  {
140  _dual_error_estimator->estimate_errors
141  (_system.get_equation_systems(), dual_errors_per_cell, &adjointsolutionvecs,
142  estimate_parent_error);
143  }
144  else // If not
145  {
146  // Just get the combined error estimate
147  _dual_error_estimator->estimate_error
148  (_system, dual_error_per_cell, &(_system.get_adjoint_solution(i)), estimate_parent_error);
149  }
150 
151  std::size_t error_size;
152 
153  // Get the size of the first ErrorMap vector; this will give us the number of elements
154  if (!error_norm_is_identity) // If in non default weights case
155  {
156  error_size = dual_errors_per_cell[std::make_pair(&_system, 0)]->size();
157  }
158  else // If in the standard default weights case
159  {
160  error_size = dual_error_per_cell.size();
161  }
162 
163  // Resize the ErrorVector(s)
164  if (!error_norm_is_identity)
165  {
166  // Loop over variables
167  for (unsigned int v = 0; v < n_vars; v++)
168  {
169  libmesh_assert(!total_dual_errors_per_cell[std::make_pair(&_system, v)]->size() ||
170  total_dual_errors_per_cell[std::make_pair(&_system, v)]->size() == error_size) ;
171  total_dual_errors_per_cell[std::make_pair(&_system, v)]->resize(error_size);
172  }
173  }
174  else
175  {
176  libmesh_assert(!total_dual_error_per_cell.size() ||
177  total_dual_error_per_cell.size() == error_size);
178  total_dual_error_per_cell.resize(error_size);
179  }
180 
181  for (std::size_t e = 0; e != error_size; ++e)
182  {
183  // Have we been asked to weight the variable error contributions in any specific manner
184  if (!error_norm_is_identity) // If we have
185  {
186  // Loop over variables
187  for (unsigned int v = 0; v < n_vars; v++)
188  {
189  // Now fill in total_dual_error ErrorMap with the weight
190  (*total_dual_errors_per_cell[std::make_pair(&_system, v)])[e] +=
191  static_cast<ErrorVectorReal>
192  (error_weight *
193  (*dual_errors_per_cell[std::make_pair(&_system, v)])[e]);
194  }
195  }
196  else // If not
197  {
198  total_dual_error_per_cell[e] +=
199  static_cast<ErrorVectorReal>(error_weight * dual_error_per_cell[e]);
200  }
201  }
202  }
203  }
204 
205  // Do some debugging plots if requested
206  if (!error_plot_suffix.empty())
207  {
208  if (!error_norm_is_identity) // If we have
209  {
210  // Loop over variables
211  for (unsigned int v = 0; v < n_vars; v++)
212  {
213  std::ostringstream primal_out;
214  std::ostringstream dual_out;
215  primal_out << "primal_" << error_plot_suffix << ".";
216  dual_out << "dual_" << error_plot_suffix << ".";
217 
218  primal_out << std::setw(1)
219  << std::setprecision(0)
220  << std::setfill('0')
221  << std::right
222  << v;
223 
224  dual_out << std::setw(1)
225  << std::setprecision(0)
226  << std::setfill('0')
227  << std::right
228  << v;
229 
230  primal_errors_per_cell[std::make_pair(&_system, v)]->plot_error(primal_out.str(), _system.get_mesh());
231  total_dual_errors_per_cell[std::make_pair(&_system, v)]->plot_error(dual_out.str(), _system.get_mesh());
232 
233  primal_out.clear();
234  dual_out.clear();
235  }
236  }
237  else // If not
238  {
239  std::ostringstream primal_out;
240  std::ostringstream dual_out;
241  primal_out << "primal_" << error_plot_suffix ;
242  dual_out << "dual_" << error_plot_suffix ;
243 
244  primal_error_per_cell.plot_error(primal_out.str(), _system.get_mesh());
245  total_dual_error_per_cell.plot_error(dual_out.str(), _system.get_mesh());
246 
247  primal_out.clear();
248  dual_out.clear();
249  }
250  }
251 
252  // Weight the primal error by the dual error using the system norm object
253  // FIXME: we ought to thread this
254  for (auto i : index_range(error_per_cell))
255  {
256  // Have we been asked to weight the variable error contributions in any specific manner
257  if (!error_norm_is_identity) // If we do
258  {
259  // Create Error Vectors to pass to calculate_norm
260  std::vector<Real> cell_primal_error;
261  std::vector<Real> cell_dual_error;
262 
263  for (unsigned int v = 0; v < n_vars; v++)
264  {
265  cell_primal_error.push_back((*primal_errors_per_cell[std::make_pair(&_system, v)])[i]);
266  cell_dual_error.push_back((*total_dual_errors_per_cell[std::make_pair(&_system, v)])[i]);
267  }
268 
269  error_per_cell[i] =
270  static_cast<ErrorVectorReal>
271  (error_norm.calculate_norm(cell_primal_error, cell_dual_error));
272  }
273  else // If not
274  {
275  error_per_cell[i] = primal_error_per_cell[i]*total_dual_error_per_cell[i];
276  }
277  }
278 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
MeshBase & mesh
QoISet _qoi_set
A QoISet to handle cases with multiple QoIs available.
std::map< std::pair< const System *, unsigned int >, std::unique_ptr< ErrorVector > > ErrorMap
When calculating many error vectors at once, we need a data structure to hold them all...
DIE A HORRIBLE DEATH HERE typedef float ErrorVectorReal
bool has_index(std::size_t) const
Return whether or not this index is in the set to be calculated.
Definition: qoi_set.h:224
std::string error_plot_suffix
To aid in investigating error estimator behavior, set this string to a suffix with which to plot (pre...
unsigned int n_vars
libmesh_assert(ctx)
std::unique_ptr< ErrorEstimator > _dual_error_estimator
An error estimator for the adjoint problem.
Real weight(std::size_t) const
Get the weight for this index (default 1.0)
Definition: qoi_set.h:243
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Real calculate_norm(const std::vector< Real > &v)
Definition: system_norm.C:232
std::unique_ptr< ErrorEstimator > _primal_error_estimator
An error estimator for the forward problem.
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ estimate_errors() [1/2]

void libMesh::ErrorEstimator::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 
)
virtualinherited

This virtual function can be redefined in derived classes, but by default computes the sum of the error_per_cell for each system in the equation_systems.

Currently this function ignores the error_norm member variable, and uses the function argument error_norms instead.

This function is named estimate_errors instead of estimate_error because otherwise C++ can get confused.

Reimplemented in libMesh::UniformRefinementEstimator.

Definition at line 47 of file error_estimator.C.

References libMesh::ErrorEstimator::error_norm, libMesh::ErrorEstimator::estimate_error(), libMesh::EquationSystems::get_system(), libMesh::index_range(), libMesh::make_range(), and libMesh::EquationSystems::n_systems().

52 {
53  SystemNorm old_error_norm = this->error_norm;
54 
55  // Sum the error values from each system
56  for (auto s : make_range(equation_systems.n_systems()))
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 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
template class LIBMESH_EXPORT NumericVector< Number >
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ estimate_errors() [2/2]

void libMesh::ErrorEstimator::estimate_errors ( const EquationSystems equation_systems,
ErrorMap errors_per_cell,
const std::map< const System *, const NumericVector< Number > *> *  solution_vectors = nullptr,
bool  estimate_parent_error = false 
)
virtualinherited

This virtual function can be redefined in derived classes, but by default it calls estimate_error repeatedly to calculate the requested error vectors.

FIXME: This is a default implementation - derived classes should reimplement it for efficiency.

Currently this function ignores the error_norm.weight() values because it calculates each variable's error individually, unscaled.

The user selects which errors get computed by filling a map with error vectors: If errors_per_cell[&system][v] exists, it will be filled with the error values in variable v of system

Reimplemented in libMesh::UniformRefinementEstimator.

Definition at line 93 of file error_estimator.C.

References libMesh::ErrorEstimator::error_norm, libMesh::ErrorEstimator::estimate_error(), libMesh::EquationSystems::get_system(), libMesh::make_range(), libMesh::EquationSystems::n_systems(), n_vars, libMesh::System::n_vars(), and libMesh::SystemNorm::type().

97 {
98  SystemNorm old_error_norm = this->error_norm;
99 
100  // Find the requested error values from each system
101  for (auto s : make_range(equation_systems.n_systems()))
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.count(std::make_pair(&sys, v)))
111  continue;
112 
113  // Calculate error in only one variable
114  std::vector<Real> weights(n_vars, 0.0);
115  weights[v] = 1.0;
116  this->error_norm =
117  SystemNorm(std::vector<FEMNormType>(n_vars, old_error_norm.type(v)),
118  weights);
119 
120  const NumericVector<Number> * solution_vector = nullptr;
121  if (solution_vectors &&
122  solution_vectors->find(&sys) != solution_vectors->end())
123  solution_vector = solution_vectors->find(&sys)->second;
124 
125  this->estimate_error
126  (sys, *errors_per_cell[std::make_pair(&sys, v)],
127  solution_vector, estimate_parent_error);
128  }
129  }
130 
131  // Restore our old state before returning
132  this->error_norm = old_error_norm;
133 }
SystemNorm error_norm
When estimating the error in a single system, the error_norm is used to control the scaling and norm ...
virtual void estimate_error(const System &system, ErrorVector &error_per_cell, const NumericVector< Number > *solution_vector=nullptr, bool estimate_parent_error=false)=0
This pure virtual function must be redefined in derived classes to compute the error for each active ...
unsigned int n_vars
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
template class LIBMESH_EXPORT NumericVector< Number >

◆ operator=() [1/2]

AdjointResidualErrorEstimator& libMesh::AdjointResidualErrorEstimator::operator= ( const AdjointResidualErrorEstimator )
delete

◆ operator=() [2/2]

AdjointResidualErrorEstimator& libMesh::AdjointResidualErrorEstimator::operator= ( AdjointResidualErrorEstimator &&  )
default

◆ primal_error_estimator()

std::unique_ptr<ErrorEstimator>& libMesh::AdjointResidualErrorEstimator::primal_error_estimator ( )
inline

Access to the "subestimator" (default: PatchRecovery) to use on the primal/forward solution.

Definition at line 77 of file adjoint_residual_error_estimator.h.

References _primal_error_estimator.

77 { return _primal_error_estimator; }
std::unique_ptr< ErrorEstimator > _primal_error_estimator
An error estimator for the forward problem.

◆ qoi_set() [1/2]

QoISet& libMesh::AdjointResidualErrorEstimator::qoi_set ( )
inline

Access to the QoISet (default: weight all QoIs equally) to use when computing errors.

Definition at line 89 of file adjoint_residual_error_estimator.h.

References _qoi_set.

89 { return _qoi_set; }
QoISet _qoi_set
A QoISet to handle cases with multiple QoIs available.

◆ qoi_set() [2/2]

const QoISet& libMesh::AdjointResidualErrorEstimator::qoi_set ( ) const
inline

Access to the QoISet (default: weight all QoIs equally) to use when computing errors.

Definition at line 95 of file adjoint_residual_error_estimator.h.

References _qoi_set.

95 { return _qoi_set; }
QoISet _qoi_set
A QoISet to handle cases with multiple QoIs available.

◆ reduce_error()

void libMesh::ErrorEstimator::reduce_error ( std::vector< ErrorVectorReal > &  error_per_cell,
const Parallel::Communicator comm 
) const
protectedinherited

This method takes the local error contributions in error_per_cell from each processor and combines them to get the global error vector.

Definition at line 32 of file error_estimator.C.

References TIMPI::Communicator::sum().

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), and libMesh::ExactErrorEstimator::estimate_error().

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 }

◆ type()

ErrorEstimatorType libMesh::AdjointResidualErrorEstimator::type ( ) const
overridevirtual
Returns
The type for the ErrorEstimator subclass.

Implements libMesh::ErrorEstimator.

Definition at line 54 of file adjoint_residual_error_estimator.C.

References libMesh::ADJOINT_RESIDUAL.

Member Data Documentation

◆ _dual_error_estimator

std::unique_ptr<ErrorEstimator> libMesh::AdjointResidualErrorEstimator::_dual_error_estimator
protected

An error estimator for the adjoint problem.

Definition at line 133 of file adjoint_residual_error_estimator.h.

Referenced by dual_error_estimator(), and estimate_error().

◆ _primal_error_estimator

std::unique_ptr<ErrorEstimator> libMesh::AdjointResidualErrorEstimator::_primal_error_estimator
protected

An error estimator for the forward problem.

Definition at line 128 of file adjoint_residual_error_estimator.h.

Referenced by estimate_error(), and primal_error_estimator().

◆ _qoi_set

QoISet libMesh::AdjointResidualErrorEstimator::_qoi_set
protected

A QoISet to handle cases with multiple QoIs available.

Definition at line 138 of file adjoint_residual_error_estimator.h.

Referenced by estimate_error(), and qoi_set().

◆ error_norm

SystemNorm libMesh::ErrorEstimator::error_norm
inherited

When estimating the error in a single system, the error_norm is used to control the scaling and norm choice for each variable.

Not all estimators will support all norm choices. The default scaling is for all variables to be weighted equally. The default norm choice depends on the error estimator.

Part of this functionality was supported via component_scale and sobolev_order in older libMesh versions, and a small part was supported via component_mask in even older versions. Hopefully the encapsulation here will allow us to avoid changing this API again.

Definition at line 158 of file error_estimator.h.

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::AdjointRefinementEstimator::AdjointRefinementEstimator(), libMesh::DiscontinuityMeasure::boundary_side_integration(), libMesh::KellyErrorEstimator::boundary_side_integration(), libMesh::DiscontinuityMeasure::DiscontinuityMeasure(), libMesh::JumpErrorEstimator::estimate_error(), estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::ErrorEstimator::estimate_errors(), libMesh::ExactErrorEstimator::ExactErrorEstimator(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::LaplacianErrorEstimator::init_context(), libMesh::DiscontinuityMeasure::init_context(), libMesh::KellyErrorEstimator::init_context(), libMesh::LaplacianErrorEstimator::internal_side_integration(), libMesh::DiscontinuityMeasure::internal_side_integration(), libMesh::KellyErrorEstimator::internal_side_integration(), libMesh::KellyErrorEstimator::KellyErrorEstimator(), libMesh::LaplacianErrorEstimator::LaplacianErrorEstimator(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::PatchRecoveryErrorEstimator(), libMesh::JumpErrorEstimator::reinit_sides(), and libMesh::UniformRefinementEstimator::UniformRefinementEstimator().

◆ error_plot_suffix

std::string libMesh::AdjointResidualErrorEstimator::error_plot_suffix

To aid in investigating error estimator behavior, set this string to a suffix with which to plot (prefixed by "primal_" or "dual_") the subestimator results.

The suffix should end with a file extension (e.g. ".gmv") that the ErrorVector::plot_error recognizes.

Definition at line 104 of file adjoint_residual_error_estimator.h.

Referenced by estimate_error().


The documentation for this class was generated from the following files: