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

#include <simulation.h>

Collaboration diagram for GRINS::Simulation:
Collaboration graph
[legend]

Public Member Functions

 Simulation (const GetPot &input, SimulationBuilder &sim_builder, const libMesh::Parallel::Communicator &comm LIBMESH_CAN_DEFAULT_TO_COMMWORLD)
 
 Simulation (const GetPot &input, GetPot &command_line, SimulationBuilder &sim_builder, const libMesh::Parallel::Communicator &comm LIBMESH_CAN_DEFAULT_TO_COMMWORLD)
 
virtual ~Simulation ()
 
void run ()
 
void print_sim_info ()
 
SharedPtr< libMesh::EquationSystems > get_equation_system ()
 
MultiphysicsSystemget_multiphysics_system ()
 
const MultiphysicsSystemget_multiphysics_system () const
 
libMesh::Number get_qoi_value (unsigned int qoi_index) const
 
const std::string & get_multiphysics_system_name () const
 

Protected Member Functions

void read_restart (const GetPot &input)
 
void init_multiphysics_system (const GetPot &input)
 Helper function. More...
 
void init_qois (const GetPot &input, SimulationBuilder &sim_builder)
 Helper function. More...
 
void init_params (const GetPot &input, SimulationBuilder &sim_builder)
 Helper function. More...
 
void init_restart (const GetPot &input, SimulationBuilder &sim_builder, const libMesh::Parallel::Communicator &comm)
 Helper function. More...
 
bool check_for_adjoint_solve (const GetPot &input) const
 Helper function. More...
 
void init_adjoint_solve (const GetPot &input, bool output_adjoint)
 Helper function. More...
 
void build_error_estimator (const GetPot &input)
 
void build_solver (const GetPot &input)
 

Protected Attributes

SharedPtr< libMesh::UnstructuredMesh > _mesh
 
SharedPtr< libMesh::EquationSystems > _equation_system
 
libMesh::UniquePtr< GRINS::Solver_solver
 
std::string _system_name
 GRINS::Multiphysics system name. More...
 
GRINS::MultiphysicsSystem_multiphysics_system
 
SharedPtr< GRINS::Visualization_vis
 
SharedPtr< PostProcessedQuantities< libMesh::Real > > _postprocessing
 
bool _print_mesh_info
 
bool _print_log_info
 
bool _print_equation_system_info
 
bool _print_constraint_info
 
bool _print_perflog
 
bool _print_scalars
 
SharedPtr< QoIOutput_qoi_output
 
bool _output_vis
 
bool _output_adjoint
 
bool _output_residual
 
bool _output_residual_sensitivities
 
bool _output_solution_sensitivities
 
unsigned int _timesteps_per_vis
 
unsigned int _timesteps_per_perflog
 
ErrorEstimatorOptions _error_estimator_options
 
SharedPtr< libMesh::ErrorEstimator > _error_estimator
 
ParameterManager _adjoint_parameters
 
ParameterManager _forward_parameters
 
bool _do_adjoint_solve
 
bool _have_restart
 

Private Member Functions

 Simulation ()
 

Detailed Description

Definition at line 58 of file simulation.h.

Constructor & Destructor Documentation

GRINS::Simulation::Simulation ( const GetPot &  input,
SimulationBuilder sim_builder,
const libMesh::Parallel::Communicator &comm  LIBMESH_CAN_DEFAULT_TO_COMMWORLD 
)
GRINS::Simulation::Simulation ( const GetPot &  input,
GetPot &  command_line,
SimulationBuilder sim_builder,
const libMesh::Parallel::Communicator &comm  LIBMESH_CAN_DEFAULT_TO_COMMWORLD 
)
virtual GRINS::Simulation::~Simulation ( )
inlinevirtual

Definition at line 73 of file simulation.h.

73 {};
GRINS::Simulation::Simulation ( )
private

Member Function Documentation

void GRINS::Simulation::build_error_estimator ( const GetPot &  input)
protected

Definition at line 491 of file simulation.C.

References _error_estimator, _error_estimator_options, _multiphysics_system, GRINS::FactoryAbstract< libMesh::ErrorEstimator >::build(), GRINS::ErrorEstimatorOptions::estimator_type(), GRINS::ErrorEstimatorFactoryBase::set_estimator_options(), GRINS::FactoryWithGetPot< libMesh::ErrorEstimator >::set_getpot(), and GRINS::ErrorEstimatorFactoryBase::set_system().

492  {
493  std::string estimator_type = _error_estimator_options.estimator_type();
494  if( estimator_type != std::string("none") )
495  {
499  _error_estimator.reset( (ErrorEstimatorFactoryBase::build(estimator_type)).release() );
500  }
501  }
static void set_estimator_options(const ErrorEstimatorOptions &estimator_options)
const std::string & estimator_type() const
static libMesh::UniquePtr< libMesh::ErrorEstimator > build(const std::string &name)
Use this method to build objects of type Base.
SharedPtr< libMesh::ErrorEstimator > _error_estimator
Definition: simulation.h:153
static void set_system(MultiphysicsSystem &system)
ErrorEstimatorOptions _error_estimator_options
Definition: simulation.h:152
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
void GRINS::Simulation::build_solver ( const GetPot &  input)
protected

Definition at line 503 of file simulation.C.

References _solver, GRINS::FactoryAbstract< Solver >::build(), GRINS::FactoryWithGetPot< Solver >::set_getpot(), and GRINS::SolverParsing::solver_type().

504  {
506 
507  std::string solver_name = SolverParsing::solver_type(input);
508 
509  _solver = SolverFactoryAbstract::build(solver_name);
510  }
static std::string solver_type(const GetPot &input)
static void set_getpot(const GetPot &input)
static libMesh::UniquePtr< Solver > build(const std::string &name)
Use this method to build objects of type Base.
libMesh::UniquePtr< GRINS::Solver > _solver
Definition: simulation.h:119
bool GRINS::Simulation::check_for_adjoint_solve ( const GetPot &  input) const
protected

Helper function.

Todo:
This option name WILL change at some point.

Definition at line 469 of file simulation.C.

References _adjoint_parameters, _error_estimator_options, GRINS::StrategiesParsing::do_adjoint_solve(), GRINS::ErrorEstimatorOptions::estimator_requires_adjoint(), GRINS::ErrorEstimatorOptions::estimator_type(), and GRINS::ParameterManager::parameter_vector.

Referenced by init_adjoint_solve().

470  {
471  std::string error_estimator = _error_estimator_options.estimator_type();
472 
473  bool do_adjoint_solve = false;
474 
475  // First check if the error estimator requires an adjoint solve
477  do_adjoint_solve = true;
478 
479  // Next check if parameter sensitivies require an adjoint solve
481  do_adjoint_solve = true;
482 
483  // Now check if the user requested to do an adjoint solve regardless
486  do_adjoint_solve = true;
487 
488  return do_adjoint_solve;
489  }
static bool do_adjoint_solve(const GetPot &input)
Option to let user manually trigger adjoint solve.
const std::string & estimator_type() const
libMesh::ParameterVector parameter_vector
ParameterManager _adjoint_parameters
Definition: simulation.h:155
ErrorEstimatorOptions _error_estimator_options
Definition: simulation.h:152
SharedPtr< libMesh::EquationSystems > GRINS::Simulation::get_equation_system ( )

Definition at line 393 of file simulation.C.

References _equation_system.

Referenced by main(), and run().

394  {
395  return _equation_system;
396  }
SharedPtr< libMesh::EquationSystems > _equation_system
Definition: simulation.h:117
MultiphysicsSystem * GRINS::Simulation::get_multiphysics_system ( )
inline

Definition at line 179 of file simulation.h.

References _multiphysics_system.

Referenced by main().

180  {
181  return this->_multiphysics_system;
182  }
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
const MultiphysicsSystem * GRINS::Simulation::get_multiphysics_system ( ) const
inline

Definition at line 172 of file simulation.h.

References _multiphysics_system.

173  {
174  return this->_multiphysics_system;
175  }
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
const std::string & GRINS::Simulation::get_multiphysics_system_name ( ) const
inline

Definition at line 185 of file simulation.h.

References _system_name.

Referenced by main().

186  {
187  return this->_system_name;
188  }
std::string _system_name
GRINS::Multiphysics system name.
Definition: simulation.h:122
libMesh::Number GRINS::Simulation::get_qoi_value ( unsigned int  qoi_index) const

Definition at line 398 of file simulation.C.

References _multiphysics_system, and GRINS::CompositeQoI::get_qoi_value().

Referenced by main().

399  {
400  const CompositeQoI* qoi = libMesh::cast_ptr<const CompositeQoI*>(this->_multiphysics_system->get_qoi());
401  return qoi->get_qoi_value(qoi_index);
402  }
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
void GRINS::Simulation::init_adjoint_solve ( const GetPot &  input,
bool  output_adjoint 
)
protected

Helper function.

Definition at line 446 of file simulation.C.

References _do_adjoint_solve, _multiphysics_system, and check_for_adjoint_solve().

447  {
448  // Check if we're doing an adjoint solve
450 
451  const libMesh::DifferentiableQoI* raw_qoi = _multiphysics_system->get_qoi();
452  const CompositeQoI* qoi = dynamic_cast<const CompositeQoI*>( raw_qoi );
453 
454  // If we are trying to do an adjoint solve without a QoI, that's an error
455  // If there are no QoIs, the CompositeQoI never gets built and qoi will be NULL
456  if( _do_adjoint_solve && !qoi )
457  {
458  libmesh_error_msg("Error: Adjoint solve requested, but no QoIs detected.");
459  }
460 
461  /* If we are not doing an adjoint solve, but adjoint output was requested:
462  that's an error. */
463  if( !_do_adjoint_solve && output_adjoint )
464  {
465  libmesh_error_msg("Error: Adjoint output requested, but no adjoint solve requested.");
466  }
467  }
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
bool check_for_adjoint_solve(const GetPot &input) const
Helper function.
Definition: simulation.C:469
void GRINS::Simulation::init_multiphysics_system ( const GetPot &  input)
protected

Helper function.

Definition at line 142 of file simulation.C.

References _equation_system, _multiphysics_system, _postprocessing, _print_log_info, _solver, GRINS::MultiphysicsSystem::attach_physics_list(), GRINS::PhysicsBuilder::build_physics_map(), GRINS::VariableBuilder::build_variables(), GRINS::StrategiesParsing::extra_quadrature_order(), GRINS::MultiphysicsSystem::read_input_options(), and GRINS::MultiphysicsSystem::register_postprocessing_vars().

143  {
144  // Only print libMesh logging info if the user requests it
145  libMesh::perflog.disable_logging();
146  if( this->_print_log_info ) libMesh::perflog.enable_logging();
147 
149 
150  PhysicsList physics_list = PhysicsBuilder::build_physics_map(input);
151 
152  _multiphysics_system->attach_physics_list( physics_list );
153 
155 
157 
158  /* Postprocessing needs to be initialized before the solver since that's
159  where equation_system gets init'ed */
161 
162  _solver->initialize( input, _equation_system, _multiphysics_system );
163 
164  // Useful for debugging
165  if( input("screen-options/print_dof_constraints", false ) )
166  {
167  _multiphysics_system->get_dof_map().print_dof_constraints();
168  }
169 
170  // This *must* be done after equation_system->init in order to get variable indices
171  // Set any extra quadrature order the user requested. By default, is 0.
172  _multiphysics_system->extra_quadrature_order = StrategiesParsing::extra_quadrature_order(input);
173 
174  if (input("Strategies/MeshAdaptivity/disable_two_step_refinement",false))
175  _equation_system->disable_refine_in_reinit();
176  }
static PhysicsList build_physics_map(const GetPot &input)
Returns container of constructed Physics objects.
void attach_physics_list(PhysicsList physics_list)
PhysicsList gets built by GRINS::PhysicsFactory and attached here.
void register_postprocessing_vars(const GetPot &input, PostProcessedQuantities< libMesh::Real > &postprocessing)
Each Physics will register their postprocessed quantities with this call.
SharedPtr< PostProcessedQuantities< libMesh::Real > > _postprocessing
Definition: simulation.h:129
static int extra_quadrature_order(const GetPot &input)
static void build_variables(const GetPot &input, MultiphysicsSystem &system)
libMesh::UniquePtr< GRINS::Solver > _solver
Definition: simulation.h:119
std::map< std::string, SharedPtr< GRINS::Physics > > PhysicsList
Container for GRINS::Physics object pointers.
Definition: var_typedefs.h:59
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
virtual void read_input_options(const GetPot &input)
Reads input options for this class and all physics that are enabled.
SharedPtr< libMesh::EquationSystems > _equation_system
Definition: simulation.h:117
void GRINS::Simulation::init_params ( const GetPot &  input,
SimulationBuilder sim_builder 
)
protected

Helper function.

Definition at line 199 of file simulation.C.

References _adjoint_parameters, _forward_parameters, _multiphysics_system, and GRINS::ParameterManager::initialize().

201  {
202  unsigned int n_adjoint_parameters =
203  input.vector_variable_size("QoI/adjoint_sensitivity_parameters");
204 
205  unsigned int n_forward_parameters =
206  input.vector_variable_size("QoI/forward_sensitivity_parameters");
207 
208  // If the user actually asks for parameter sensitivities, then we
209  // set up the parameter vectors to use.
210  if ( n_adjoint_parameters )
211  {
212  // If we're doing adjoint sensitivities, dq/dp only makes
213  // sense if we have q
214  CompositeQoI* qoi =
215  libMesh::cast_ptr<CompositeQoI*>
216  (this->_multiphysics_system->get_qoi());
217 
218  if (!qoi)
219  {
220  std::cout <<
221  "Error: adjoint_sensitivity_parameters are specified but\n"
222  << "no QoIs have been specified.\n" << std::endl;
223  libmesh_error();
224  }
225 
227  (input, "QoI/adjoint_sensitivity_parameters",
228  *this->_multiphysics_system, qoi);
229  }
230 
231  if ( n_forward_parameters )
232  {
233  // If we're doing forward sensitivities, du/dp can make
234  // sense even with no q defined
235  CompositeQoI* qoi =
236  dynamic_cast<CompositeQoI*>
237  (this->_multiphysics_system->get_qoi());
238 
239  // dynamic_cast returns NULL if our QoI isn't a CompositeQoI;
240  // i.e. if there were no QoIs that made us bother setting up
241  // the CompositeQoI object. Passing NULL tells
242  // ParameterManager not to bother asking for qoi registration
243  // of parameters.
244 
246  (input, "QoI/forward_sensitivity_parameters",
247  *this->_multiphysics_system, qoi);
248  }
249  }
ParameterManager _forward_parameters
Definition: simulation.h:157
ParameterManager _adjoint_parameters
Definition: simulation.h:155
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
virtual void initialize(const GetPot &input, const std::string &parameters_varname, GRINS::MultiphysicsSystem &system, GRINS::CompositeQoI *qoi)
void GRINS::Simulation::init_qois ( const GetPot &  input,
SimulationBuilder sim_builder 
)
protected

Helper function.

Definition at line 178 of file simulation.C.

References _multiphysics_system, _qoi_output, and GRINS::SimulationBuilder::build_qoi().

179  {
180  // If the user actually asks for a QoI, then we add it.
181  SharedPtr<CompositeQoI> qois = sim_builder.build_qoi( input );
182  if( qois->n_qois() > 0 )
183  {
184  // This *must* be done after equation_system->init in order to get variable indices
185  qois->init(input, *_multiphysics_system );
186 
187  /* Note that we are effectively transfering ownership of the qoi pointer because
188  it will be cloned in _multiphysics_system and all the calculations are done there. */
189  _multiphysics_system->attach_qoi( qois.get() );
190  }
191  else if (_qoi_output->output_qoi_set())
192  {
193  std::cout << "Error: print_qoi is specified but\n" <<
194  "no QoIs have been specified.\n" << std::endl;
195  libmesh_error();
196  }
197  }
SharedPtr< QoIOutput > _qoi_output
Definition: simulation.h:140
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
void GRINS::Simulation::init_restart ( const GetPot &  input,
SimulationBuilder sim_builder,
const libMesh::Parallel::Communicator &  comm 
)
protected

Helper function.

Definition at line 252 of file simulation.C.

References _equation_system, _mesh, GRINS::MeshBuilder::do_mesh_refinement_from_input(), GRINS::SimulationBuilder::mesh_builder(), and read_restart().

254  {
255  this->read_restart( input );
256 
257  /* We do this here only if there's a restart file. Otherwise, this was done
258  at mesh construction time */
259  sim_builder.mesh_builder().do_mesh_refinement_from_input( input, comm, *_mesh );
260 
261  /* \todo Any way to tell if the mesh got refined so we don't unnecessarily
262  call reinit()? */
263  _equation_system->reinit();
264  }
void read_restart(const GetPot &input)
Definition: simulation.C:404
SharedPtr< libMesh::UnstructuredMesh > _mesh
Definition: simulation.h:115
SharedPtr< libMesh::EquationSystems > _equation_system
Definition: simulation.h:117
void GRINS::Simulation::print_sim_info ( )

Definition at line 376 of file simulation.C.

References _equation_system, _mesh, _multiphysics_system, _print_constraint_info, _print_equation_system_info, and _print_mesh_info.

Referenced by run().

377  {
378  // Print mesh info if the user wants it
379  if( this->_print_mesh_info )
380  this->_mesh->print_info();
381 
382  // Print EquationSystems info if requested
383  if( this->_print_equation_system_info )
384  this->_equation_system->print_info();
385 
386  // Print DofMap constraint info if requested
387  if( this->_print_constraint_info )
388  this->_multiphysics_system->get_dof_map().print_dof_constraints();
389 
390  return;
391  }
bool _print_constraint_info
Definition: simulation.h:135
bool _print_equation_system_info
Definition: simulation.h:134
SharedPtr< libMesh::UnstructuredMesh > _mesh
Definition: simulation.h:115
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
SharedPtr< libMesh::EquationSystems > _equation_system
Definition: simulation.h:117
void GRINS::Simulation::read_restart ( const GetPot &  input)
protected

Definition at line 404 of file simulation.C.

References _equation_system, _have_restart, and GRINS::SimulationParsing::restart_file().

Referenced by init_restart().

405  {
406  const std::string restart_file = SimulationParsing::restart_file(input);
407 
408  // Most of this was pulled from FIN-S
409  if (restart_file != "none")
410  {
411  std::cout << " ====== Restarting from " << restart_file << std::endl;
412 
413  // Must have correct file type to restart
414  if (restart_file.rfind(".xdr") < restart_file.size())
415  _equation_system->read(restart_file,GRINSEnums::DECODE,
416  //libMesh::EquationSystems::READ_HEADER | // Allow for thermochemistry upgrades
417  libMesh::EquationSystems::READ_DATA |
418  libMesh::EquationSystems::READ_ADDITIONAL_DATA);
419 
420  else if (restart_file.rfind(".xda") < restart_file.size())
421  _equation_system->read(restart_file,GRINSEnums::READ,
422  //libMesh::EquationSystems::READ_HEADER | // Allow for thermochemistry upgrades
423  libMesh::EquationSystems::READ_DATA |
424  libMesh::EquationSystems::READ_ADDITIONAL_DATA);
425 
426  else
427  {
428  std::cerr << "Error: Restart filename must have .xdr or .xda extension!" << std::endl;
429  libmesh_error();
430  }
431 
432  this->_have_restart = true;
433 
434  const std::string system_name = input("screen-options/system_name", "GRINS" );
435 
436  MultiphysicsSystem& system =
437  _equation_system->get_system<MultiphysicsSystem>(system_name);
438 
439  // Update the old data
440  system.update();
441  }
442 
443  return;
444  }
static std::string restart_file(const GetPot &input)
SharedPtr< libMesh::EquationSystems > _equation_system
Definition: simulation.h:117
void GRINS::Simulation::run ( )

Definition at line 267 of file simulation.C.

References _adjoint_parameters, _do_adjoint_solve, _equation_system, _error_estimator, _forward_parameters, _have_restart, _multiphysics_system, _output_adjoint, _output_residual, _output_residual_sensitivities, _output_solution_sensitivities, _output_vis, _postprocessing, _print_log_info, _print_scalars, _qoi_output, _solver, _timesteps_per_perflog, _timesteps_per_vis, _vis, GRINS::SolverContext::do_adjoint_solve, GRINS::SolverContext::equation_system, GRINS::SolverContext::error_estimator, GRINS::SolverContext::have_restart, GRINS::SolverContext::output_adjoint, GRINS::SolverContext::output_residual, GRINS::SolverContext::output_residual_sensitivities, GRINS::SolverContext::output_solution_sensitivities, GRINS::SolverContext::output_vis, GRINS::ParameterManager::parameter_vector, GRINS::SolverContext::postprocessing, GRINS::SolverContext::print_perflog, GRINS::SolverContext::print_scalars, print_sim_info(), GRINS::SolverContext::qoi_output, GRINS::SolverContext::system, GRINS::SolverContext::timesteps_per_perflog, GRINS::SolverContext::timesteps_per_vis, and GRINS::SolverContext::vis.

268  {
269  this->print_sim_info();
270 
271  SolverContext context;
272  context.system = _multiphysics_system;
273  context.equation_system = _equation_system;
274  context.vis = _vis;
275  context.output_adjoint = _output_adjoint;
276  context.timesteps_per_vis = _timesteps_per_vis;
277  context.timesteps_per_perflog = _timesteps_per_perflog;
278  context.output_vis = _output_vis;
279  context.output_residual = _output_residual;
280  context.output_residual_sensitivities = _output_residual_sensitivities;
281  context.output_solution_sensitivities = _output_solution_sensitivities;
282  context.print_scalars = _print_scalars;
283  context.print_perflog = _print_log_info;
284  context.postprocessing = _postprocessing;
285  context.error_estimator = _error_estimator;
286  context.qoi_output = _qoi_output;
287  context.do_adjoint_solve = _do_adjoint_solve;
288  context.have_restart = _have_restart;
289 
292  {
293  std::cout <<
294  "Error: output_residual_sensitivities is specified but\n" <<
295  "no forward sensitivity parameters have been specified.\n" <<
296  std::endl;
297  libmesh_error();
298  }
299 
302  {
303  std::cout <<
304  "Error: output_solution_sensitivities is specified but\n" <<
305  "no forward sensitivity parameters have been specified.\n" <<
306  std::endl;
307  libmesh_error();
308  }
309 
310  _solver->solve( context );
311 
312  if (_qoi_output->output_qoi_set())
313  {
314  _multiphysics_system->assemble_qoi();
315  const CompositeQoI * my_qoi = libMesh::cast_ptr<const CompositeQoI*>(this->_multiphysics_system->get_qoi());
316  _qoi_output->output_qois(*my_qoi, this->_multiphysics_system->comm() );
317  }
318 
320  {
321  // Default: "calculate sensitivities for all QoIs"
322  libMesh::QoISet qois;
323 
324  const libMesh::ParameterVector & params =
326 
327  libMesh::SensitivityData sensitivities
328  (qois, *this->_multiphysics_system, params);
329 
330  _solver->adjoint_qoi_parameter_sensitivity
331  (context, qois, params, sensitivities);
332 
333  std::cout << "Adjoint sensitivities:" << std::endl;
334 
335  for (unsigned int q=0;
336  q != this->_multiphysics_system->qoi.size(); ++q)
337  {
338  for (unsigned int p=0; p != params.size(); ++p)
339  {
340  std::cout << "dq" << q << "/dp" << p << " = " <<
341  sensitivities[q][p] << std::endl;
342  }
343  }
344  }
345 
347  {
348  // Default: "calculate sensitivities for all QoIs"
349  libMesh::QoISet qois;
350 
351  const libMesh::ParameterVector & params =
353 
354  libMesh::SensitivityData sensitivities
355  (qois, *this->_multiphysics_system, params);
356 
357  _solver->forward_qoi_parameter_sensitivity
358  (context, qois, params, sensitivities);
359 
360  std::cout << "Forward sensitivities:" << std::endl;
361 
362  for (unsigned int q=0;
363  q != this->_multiphysics_system->qoi.size(); ++q)
364  {
365  for (unsigned int p=0; p != params.size(); ++p)
366  {
367  std::cout << "dq" << q << "/dp" << p << " = " <<
368  sensitivities[q][p] << std::endl;
369  }
370  }
371  }
372 
373  return;
374  }
ParameterManager _forward_parameters
Definition: simulation.h:157
unsigned int _timesteps_per_vis
Definition: simulation.h:149
SharedPtr< PostProcessedQuantities< libMesh::Real > > _postprocessing
Definition: simulation.h:129
SharedPtr< QoIOutput > _qoi_output
Definition: simulation.h:140
SharedPtr< libMesh::ErrorEstimator > _error_estimator
Definition: simulation.h:153
libMesh::ParameterVector parameter_vector
libMesh::UniquePtr< GRINS::Solver > _solver
Definition: simulation.h:119
ParameterManager _adjoint_parameters
Definition: simulation.h:155
bool _output_residual_sensitivities
Definition: simulation.h:146
bool _output_solution_sensitivities
Definition: simulation.h:147
GRINS::MultiphysicsSystem * _multiphysics_system
Definition: simulation.h:125
unsigned int _timesteps_per_perflog
Definition: simulation.h:150
SharedPtr< GRINS::Visualization > _vis
Definition: simulation.h:127
void print_sim_info()
Definition: simulation.C:376
SharedPtr< libMesh::EquationSystems > _equation_system
Definition: simulation.h:117

Member Data Documentation

ParameterManager GRINS::Simulation::_adjoint_parameters
protected

Definition at line 155 of file simulation.h.

Referenced by check_for_adjoint_solve(), init_params(), and run().

bool GRINS::Simulation::_do_adjoint_solve
protected

Definition at line 160 of file simulation.h.

Referenced by init_adjoint_solve(), and run().

SharedPtr<libMesh::EquationSystems> GRINS::Simulation::_equation_system
protected
SharedPtr<libMesh::ErrorEstimator> GRINS::Simulation::_error_estimator
protected

Definition at line 153 of file simulation.h.

Referenced by build_error_estimator(), and run().

ErrorEstimatorOptions GRINS::Simulation::_error_estimator_options
protected

Definition at line 152 of file simulation.h.

Referenced by build_error_estimator(), and check_for_adjoint_solve().

ParameterManager GRINS::Simulation::_forward_parameters
protected

Definition at line 157 of file simulation.h.

Referenced by init_params(), and run().

bool GRINS::Simulation::_have_restart
protected

Definition at line 162 of file simulation.h.

Referenced by read_restart(), and run().

SharedPtr<libMesh::UnstructuredMesh> GRINS::Simulation::_mesh
protected

Definition at line 115 of file simulation.h.

Referenced by init_restart(), and print_sim_info().

GRINS::MultiphysicsSystem* GRINS::Simulation::_multiphysics_system
protected
bool GRINS::Simulation::_output_adjoint
protected

Definition at line 144 of file simulation.h.

Referenced by run().

bool GRINS::Simulation::_output_residual
protected

Definition at line 145 of file simulation.h.

Referenced by run().

bool GRINS::Simulation::_output_residual_sensitivities
protected

Definition at line 146 of file simulation.h.

Referenced by run().

bool GRINS::Simulation::_output_solution_sensitivities
protected

Definition at line 147 of file simulation.h.

Referenced by run().

bool GRINS::Simulation::_output_vis
protected

Definition at line 143 of file simulation.h.

Referenced by run().

SharedPtr<PostProcessedQuantities<libMesh::Real> > GRINS::Simulation::_postprocessing
protected

Definition at line 129 of file simulation.h.

Referenced by init_multiphysics_system(), and run().

bool GRINS::Simulation::_print_constraint_info
protected

Definition at line 135 of file simulation.h.

Referenced by print_sim_info().

bool GRINS::Simulation::_print_equation_system_info
protected

Definition at line 134 of file simulation.h.

Referenced by print_sim_info().

bool GRINS::Simulation::_print_log_info
protected

Definition at line 133 of file simulation.h.

Referenced by init_multiphysics_system(), and run().

bool GRINS::Simulation::_print_mesh_info
protected

Definition at line 132 of file simulation.h.

Referenced by print_sim_info().

bool GRINS::Simulation::_print_perflog
protected

Definition at line 136 of file simulation.h.

bool GRINS::Simulation::_print_scalars
protected

Definition at line 137 of file simulation.h.

Referenced by run().

SharedPtr<QoIOutput> GRINS::Simulation::_qoi_output
protected

Definition at line 140 of file simulation.h.

Referenced by init_qois(), and run().

libMesh::UniquePtr<GRINS::Solver> GRINS::Simulation::_solver
protected

Definition at line 119 of file simulation.h.

Referenced by build_solver(), init_multiphysics_system(), and run().

std::string GRINS::Simulation::_system_name
protected

GRINS::Multiphysics system name.

Definition at line 122 of file simulation.h.

Referenced by get_multiphysics_system_name().

unsigned int GRINS::Simulation::_timesteps_per_perflog
protected

Definition at line 150 of file simulation.h.

Referenced by run().

unsigned int GRINS::Simulation::_timesteps_per_vis
protected

Definition at line 149 of file simulation.h.

Referenced by run().

SharedPtr<GRINS::Visualization> GRINS::Simulation::_vis
protected

Definition at line 127 of file simulation.h.

Referenced by run().


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

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