GRINS-0.8.0
List of all members | Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | Private Member Functions
GRINS::MeshAdaptiveSolverBase Class Reference

#include <mesh_adaptive_solver_base.h>

Inheritance diagram for GRINS::MeshAdaptiveSolverBase:
Inheritance graph
[legend]
Collaboration diagram for GRINS::MeshAdaptiveSolverBase:
Collaboration graph
[legend]

Public Member Functions

 MeshAdaptiveSolverBase (const GetPot &input)
 
virtual ~MeshAdaptiveSolverBase ()
 

Protected Types

enum  RefinementFlaggingType {
  INVALID = 0, ERROR_TOLERANCE, N_ELEM_TARGET, ERROR_FRACTION,
  ELEM_FRACTION, MEAN_STD_DEV
}
 

Protected Member Functions

void build_mesh_refinement (libMesh::MeshBase &mesh)
 
void set_refinement_type (const GetPot &input, const MeshAdaptivityOptions &mesh_adaptivity_options, RefinementFlaggingType &refinement_type)
 
bool check_for_convergence (SolverContext &context, const libMesh::ErrorVector &error) const
 
void flag_elements_for_refinement (const libMesh::ErrorVector &error)
 
void estimate_error_for_amr (SolverContext &context, libMesh::ErrorVector &error)
 
void perform_amr (SolverContext &context, const libMesh::ErrorVector &error)
 

Protected Attributes

ErrorEstimatorOptions _error_estimator_options
 
MeshAdaptivityOptions _mesh_adaptivity_options
 
RefinementFlaggingType _refinement_type
 
libMesh::UniquePtr< libMesh::MeshRefinement > _mesh_refinement
 

Private Member Functions

 MeshAdaptiveSolverBase ()
 

Detailed Description

Definition at line 52 of file mesh_adaptive_solver_base.h.

Member Enumeration Documentation

Constructor & Destructor Documentation

GRINS::MeshAdaptiveSolverBase::MeshAdaptiveSolverBase ( const GetPot &  input)

Definition at line 40 of file mesh_adaptive_solver_base.C.

References _mesh_adaptivity_options, _refinement_type, and set_refinement_type().

41  : _error_estimator_options(input),
44  {
46  }
void set_refinement_type(const GetPot &input, const MeshAdaptivityOptions &mesh_adaptivity_options, RefinementFlaggingType &refinement_type)
ErrorEstimatorOptions _error_estimator_options
MeshAdaptivityOptions _mesh_adaptivity_options
virtual GRINS::MeshAdaptiveSolverBase::~MeshAdaptiveSolverBase ( )
inlinevirtual

Definition at line 58 of file mesh_adaptive_solver_base.h.

58 {}
GRINS::MeshAdaptiveSolverBase::MeshAdaptiveSolverBase ( )
private

Member Function Documentation

void GRINS::MeshAdaptiveSolverBase::build_mesh_refinement ( libMesh::MeshBase &  mesh)
protected

Definition at line 48 of file mesh_adaptive_solver_base.C.

References _mesh_adaptivity_options, _mesh_refinement, GRINS::MeshAdaptivityOptions::absolute_global_tolerance(), GRINS::MeshAdaptivityOptions::coarsen_by_parents(), GRINS::MeshAdaptivityOptions::coarsen_fraction(), GRINS::MeshAdaptivityOptions::coarsen_threshold(), GRINS::MeshAdaptivityOptions::edge_level_mismatch_limit(), GRINS::MeshAdaptivityOptions::enforce_mismatch_limit_prior_to_refinement(), GRINS::MeshAdaptivityOptions::face_level_mismatch_limit(), GRINS::MeshAdaptivityOptions::max_h_level(), GRINS::MeshAdaptivityOptions::nelem_target(), GRINS::MeshAdaptivityOptions::node_level_mismatch_limit(), and GRINS::MeshAdaptivityOptions::refine_fraction().

Referenced by GRINS::SteadyMeshAdaptiveSolver::solve(), and GRINS::UnsteadyMeshAdaptiveSolver::solve().

49  {
50  _mesh_refinement.reset( new libMesh::MeshRefinement( mesh ) );
60  _mesh_refinement->enforce_mismatch_limit_prior_to_refinement() = _mesh_adaptivity_options.enforce_mismatch_limit_prior_to_refinement();
62  }
libMesh::Real coarsen_fraction() const
unsigned int edge_level_mismatch_limit() const
libMesh::Real absolute_global_tolerance() const
unsigned int face_level_mismatch_limit() const
libMesh::Real coarsen_threshold() const
libMesh::Real refine_fraction() const
bool enforce_mismatch_limit_prior_to_refinement() const
unsigned int node_level_mismatch_limit() const
libMesh::UniquePtr< libMesh::MeshRefinement > _mesh_refinement
MeshAdaptivityOptions _mesh_adaptivity_options
bool GRINS::MeshAdaptiveSolverBase::check_for_convergence ( SolverContext context,
const libMesh::ErrorVector &  error 
) const
protected

Definition at line 120 of file mesh_adaptive_solver_base.C.

References _mesh_adaptivity_options, GRINS::MeshAdaptivityOptions::absolute_global_tolerance(), and GRINS::SolverContext::do_adjoint_solve.

Referenced by GRINS::SteadyMeshAdaptiveSolver::solve(), and GRINS::UnsteadyMeshAdaptiveSolver::solve().

122  {
123  std::cout << "==========================================================" << std::endl
124  << "Checking convergence" << std::endl
125  << "==========================================================" << std::endl;
126 
127  bool converged = false;
128 
129  libMesh::Real error_estimate = 0.0;
130 
131  if( context.do_adjoint_solve )
132  {
133  error_estimate = std::accumulate( error.begin(), error.end(), 0.0 );
134  }
135  else
136  {
137  error_estimate = error.l2_norm();
138  }
139 
140  std::cout << "==========================================================" << std::endl
141  << "Error estimate = " << error_estimate << std::endl
142  << "==========================================================" << std::endl;
143 
144  // For now, we just check the norm
145  if( std::fabs(error_estimate) <= _mesh_adaptivity_options.absolute_global_tolerance() )
146  {
147  converged = true;
148  }
149 
150  return converged;
151  }
libMesh::Real absolute_global_tolerance() const
MeshAdaptivityOptions _mesh_adaptivity_options
void GRINS::MeshAdaptiveSolverBase::estimate_error_for_amr ( SolverContext context,
libMesh::ErrorVector &  error 
)
protected

Definition at line 205 of file mesh_adaptive_solver_base.C.

References _mesh_adaptivity_options, GRINS::SolverContext::equation_system, GRINS::SolverContext::error_estimator, GRINS::MeshAdaptivityOptions::error_plot_prefix(), GRINS::MeshAdaptivityOptions::plot_cell_errors(), and GRINS::SolverContext::system.

Referenced by GRINS::SteadyMeshAdaptiveSolver::solve(), and GRINS::UnsteadyMeshAdaptiveSolver::solve().

206  {
207  std::cout << "==========================================================" << std::endl
208  << "Estimating error" << std::endl
209  << "==========================================================" << std::endl;
210  context.error_estimator->estimate_error( *context.system, error );
211 
212  libMesh::MeshBase& mesh = context.equation_system->get_mesh();
213 
214  // Plot error vector
216  {
217  error.plot_error( _mesh_adaptivity_options.error_plot_prefix()+".exo", mesh );
218  }
219  }
const std::string & error_plot_prefix() const
MeshAdaptivityOptions _mesh_adaptivity_options
void GRINS::MeshAdaptiveSolverBase::flag_elements_for_refinement ( const libMesh::ErrorVector &  error)
protected

Definition at line 153 of file mesh_adaptive_solver_base.C.

References _mesh_refinement, _refinement_type, ELEM_FRACTION, ERROR_FRACTION, ERROR_TOLERANCE, INVALID, MEAN_STD_DEV, and N_ELEM_TARGET.

Referenced by perform_amr().

154  {
155  switch(_refinement_type)
156  {
157  case(ERROR_TOLERANCE):
158  {
159  _mesh_refinement->flag_elements_by_error_tolerance( error );
160  }
161  break;
162 
163  case(N_ELEM_TARGET):
164  {
165  _mesh_refinement->flag_elements_by_nelem_target( error );
166  }
167  break;
168 
169  case( ERROR_FRACTION ):
170  {
171  _mesh_refinement->flag_elements_by_error_fraction( error );
172  }
173  break;
174 
175  case( ELEM_FRACTION ):
176  {
177  _mesh_refinement->flag_elements_by_elem_fraction( error );
178  }
179  break;
180 
181  case( MEAN_STD_DEV ):
182  {
183  _mesh_refinement->flag_elements_by_mean_stddev( error );
184  }
185  break;
186 
187  case(INVALID):
188  {
189  std::cerr << "Error: Invalid refinement option!" << std::endl;
190  libmesh_error();
191  }
192  break;
193 
194  // Wat?!
195  default:
196  {
197  libmesh_error();
198  }
199 
200  } // switch(_refinement_type)
201 
202  return;
203  }
libMesh::UniquePtr< libMesh::MeshRefinement > _mesh_refinement
void GRINS::MeshAdaptiveSolverBase::perform_amr ( SolverContext context,
const libMesh::ErrorVector &  error 
)
protected

Definition at line 221 of file mesh_adaptive_solver_base.C.

References _mesh_refinement, GRINS::SolverContext::equation_system, flag_elements_for_refinement(), and GRINS::SolverContext::system.

Referenced by GRINS::SteadyMeshAdaptiveSolver::solve(), and GRINS::UnsteadyMeshAdaptiveSolver::solve().

222  {
223  libMesh::MeshBase& mesh = context.equation_system->get_mesh();
224 
225  std::cout << "==========================================================" << std::endl
226  << "Performing Mesh Refinement" << std::endl
227  << "==========================================================" << std::endl;
228 
229  this->flag_elements_for_refinement( error );
230  _mesh_refinement->refine_and_coarsen_elements();
231 
232  // Dont forget to reinit the system after each adaptive refinement!
233  context.equation_system->reinit();
234 
235  // This output cannot be toggled in the input file.
236  std::cout << "==========================================================" << std::endl
237  << "Refined mesh to " << std::setw(12) << mesh.n_active_elem()
238  << " active elements" << std::endl
239  << " " << std::setw(16) << context.system->n_active_dofs()
240  << " active dofs" << std::endl
241  << "==========================================================" << std::endl;
242  }
void flag_elements_for_refinement(const libMesh::ErrorVector &error)
libMesh::UniquePtr< libMesh::MeshRefinement > _mesh_refinement
void GRINS::MeshAdaptiveSolverBase::set_refinement_type ( const GetPot &  input,
const MeshAdaptivityOptions mesh_adaptivity_options,
MeshAdaptiveSolverBase::RefinementFlaggingType refinement_type 
)
protected

Definition at line 64 of file mesh_adaptive_solver_base.C.

References ELEM_FRACTION, ERROR_FRACTION, ERROR_TOLERANCE, MEAN_STD_DEV, N_ELEM_TARGET, and GRINS::MeshAdaptivityOptions::refinement_strategy().

Referenced by MeshAdaptiveSolverBase().

67  {
68  const std::string& refinement_stategy = mesh_adaptivity_options.refinement_strategy();
69 
70  if( refinement_stategy == std::string("error_tolerance") )
71  {
72  refinement_type = ERROR_TOLERANCE;
73  }
74 
75  else if( refinement_stategy == std::string("nelem_target") )
76  {
77  refinement_type = N_ELEM_TARGET;
78 
79  if( !input.have_variable("MeshAdaptivity/nelem_target") )
80  {
81  std::cerr << "Error: Must specify nelem_target for" << std::endl
82  << " for error_tolerance refinement strategy."
83  << std::endl;
84 
85  libmesh_error();
86  }
87  }
88 
89  else if( refinement_stategy == std::string("error_fraction") )
90  {
91  refinement_type = ERROR_FRACTION;
92  }
93 
94  else if( refinement_stategy == std::string("elem_fraction") )
95  {
96  refinement_type = ELEM_FRACTION;
97  }
98 
99  else if( refinement_stategy == std::string("mean_std_dev") )
100  {
101  refinement_type = MEAN_STD_DEV;
102  }
103 
104  else
105  {
106  std::cerr << "Error: Invalid refinement strategy " << refinement_stategy << std::endl
107  << " Valid refinement strategy options are: absolute_global_tolerance" << std::endl
108  << " error_tolerance" << std::endl
109  << " nelem_target" << std::endl
110  << " error_fraction" << std::endl
111  << " elem_fraction" << std::endl
112  << " mean_std_dev" << std::endl;
113 
114  libmesh_error();
115  }
116 
117  return;
118  }

Member Data Documentation

ErrorEstimatorOptions GRINS::MeshAdaptiveSolverBase::_error_estimator_options
protected
MeshAdaptivityOptions GRINS::MeshAdaptiveSolverBase::_mesh_adaptivity_options
protected
libMesh::UniquePtr<libMesh::MeshRefinement> GRINS::MeshAdaptiveSolverBase::_mesh_refinement
protected
RefinementFlaggingType GRINS::MeshAdaptiveSolverBase::_refinement_type
protected

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

Generated on Tue Dec 19 2017 12:47:31 for GRINS-0.8.0 by  doxygen 1.8.9.1