GRINS-0.8.0
unsteady_mesh_adaptive_solver.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // GRINS - General Reacting Incompressible Navier-Stokes
5 //
6 // Copyright (C) 2014-2017 Paul T. Bauman, Roy H. Stogner
7 // Copyright (C) 2010-2013 The PECOS Development Team
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the Version 2.1 GNU Lesser General
11 // Public License as published by the Free Software Foundation.
12 //
13 // This library is distributed in the hope that it will be useful,
14 // but WITHOUT ANY WARRANTY; without even the implied warranty of
15 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 // Lesser General Public License for more details.
17 //
18 // You should have received a copy of the GNU Lesser General Public
19 // License along with this library; if not, write to the Free Software
20 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
21 // Boston, MA 02110-1301 USA
22 //
23 //-----------------------------------------------------------------------el-
24 
25 
26 // This class
28 
29 // GRINS
30 #include "grins/common.h"
31 #include "grins/solver_context.h"
32 #include "grins/multiphysics_sys.h"
33 
34 // libMesh
35 #include "libmesh/getpot.h"
36 #include "libmesh/error_vector.h"
37 
38 namespace GRINS
39 {
41  : UnsteadySolver(input),
42  MeshAdaptiveSolverBase( input )
43  {}
44 
46  {
47  context.system->deltat = this->_deltat;
48 
49  libMesh::Real sim_time;
50 
51  if( context.output_vis )
52  {
53  context.postprocessing->update_quantities( *(context.equation_system) );
54  context.vis->output( context.equation_system );
55  }
56 
57  // Setup MeshRefinement
58  libMesh::MeshBase& mesh = context.equation_system->get_mesh();
59  this->build_mesh_refinement( mesh );
60 
61  std::time_t first_wall_time = std::time(NULL);
62 
63  // Now we begin the timestep loop to compute the time-accurate
64  // solution of the equations.
65  for (unsigned int t_step=0; t_step < this->_n_timesteps; t_step++)
66  {
67  std::time_t latest_wall_time = std::time(NULL);
68 
69  std::cout << "==========================================================" << std::endl
70  << " Beginning time step " << t_step <<
71  ", t = " << context.system->time <<
72  ", dt = " << context.system->deltat <<
73  ", runtime = " << (latest_wall_time - first_wall_time) <<
74  std::endl
75  << "==========================================================" << std::endl;
76 
77  // If we have any solution-dependent Dirichlet boundaries, we
78  // need to update them with the current solution.
79  this->update_dirichlet_bcs(context);
80 
81  for ( unsigned int r_step = 0; r_step < _mesh_adaptivity_options.max_refinement_steps(); r_step++ )
82  {
83  std::cout << "==========================================================" << std::endl
84  << "Adaptive Refinement Step " << r_step << std::endl
85  << "==========================================================" << std::endl;
86 
87  context.system->solve();
88 
89  sim_time = context.system->time;
90 
91  if( context.output_vis && !((t_step+1)%context.timesteps_per_vis) )
92  {
93  context.postprocessing->update_quantities( *(context.equation_system) );
94  context.vis->output( context.equation_system, t_step, sim_time );
95  }
96 
97  if( context.output_residual && !((t_step+1)%context.timesteps_per_vis) )
98  context.vis->output_residual( context.equation_system, context.system,
99  t_step, sim_time );
100 
101  if ( context.print_perflog && context.timesteps_per_perflog
102  && !((t_step+1)%context.timesteps_per_perflog) )
103  libMesh::perflog.print_log();
104 
105  if ( context.print_scalars )
106  this->print_scalar_vars(context);
107 
108 
109  // Now we construct the data structures for the mesh refinement process
110  libMesh::ErrorVector error;
111  this->estimate_error_for_amr(context,error);
112 
113  // Check for convergence of error
114  bool converged = this->check_for_convergence( context, error );
115 
116  if( converged )
117  {
118  // Break out of adaptive loop
119  std::cout << "==========================================================" << std::endl
120  << "Convergence detected!" << std::endl
121  << "==========================================================" << std::endl;
122  break;
123  }
124  else
125  {
126  // Only bother refining if we're not on the last step.
128  this->perform_amr(context,error);
129  }
130 
131  } // End mesh adaptive loop
132 
133  // Advance to the next timestep
134  context.system->time_solver->advance_timestep();
135 
136  } // End time step loop
137 
138  std::time_t final_wall_time = std::time(NULL);
139  std::cout << "==========================================================" << std::endl
140  << " Ending time stepping, t = " << context.system->time <<
141  ", runtime = " << (final_wall_time - first_wall_time) <<
142  std::endl
143  << "==========================================================" << std::endl;
144 
145  // Print out the QoI, but only do it if the user asks for it
146  if(context.qoi_output->output_qoi_set())
147  this->print_qoi(context);
148  }
149 
150 } // namespace GRINS
unsigned int _n_timesteps
SharedPtr< libMesh::EquationSystems > equation_system
SharedPtr< PostProcessedQuantities< libMesh::Real > > postprocessing
void build_mesh_refinement(libMesh::MeshBase &mesh)
virtual void solve(SolverContext &context)
void perform_amr(SolverContext &context, const libMesh::ErrorVector &error)
void estimate_error_for_amr(SolverContext &context, libMesh::ErrorVector &error)
bool check_for_convergence(SolverContext &context, const libMesh::ErrorVector &error) const
GRINS namespace.
unsigned int timesteps_per_perflog
SharedPtr< QoIOutput > qoi_output
unsigned int timesteps_per_vis
GRINS::MultiphysicsSystem * system
SharedPtr< GRINS::Visualization > vis
void print_scalar_vars(SolverContext &context)
Definition: solver.C:133
Simple class to hold objects passed to Solver::solve.
void print_qoi(SolverContext &context)
Definition: solver.C:155
unsigned int max_refinement_steps() const
void update_dirichlet_bcs(SolverContext &context)
Updates Dirichlet boundary conditions.
MeshAdaptivityOptions _mesh_adaptivity_options

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