22 #ifndef _aspect_simulator_h 23 #define _aspect_simulator_h 25 #include <deal.II/base/timer.h> 26 #include <deal.II/base/parameter_handler.h> 27 #include <deal.II/base/conditional_ostream.h> 28 #include <deal.II/base/symmetric_tensor.h> 30 DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
32 #include <deal.II/lac/affine_constraints.h> 34 #include <deal.II/distributed/tria.h> 36 #include <deal.II/dofs/dof_handler.h> 37 #include <deal.II/dofs/dof_tools.h> 39 #include <deal.II/fe/fe_system.h> 40 #include <deal.II/fe/mapping.h> 41 #include <deal.II/base/tensor_function.h> 43 DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
71 #include <boost/iostreams/tee.hpp> 72 #include <boost/iostreams/stream.hpp> 89 namespace StokesSolver
95 template <
int dim,
int velocity_degree>
98 namespace MeshDeformation
153 scalar_moment_of_inertia(numbers::signaling_nan<double>()),
154 scalar_angular_momentum(numbers::signaling_nan<double>()),
155 scalar_rotation(numbers::signaling_nan<double>()),
156 tensor_moment_of_inertia(numbers::signaling_nan<SymmetricTensor<2,dim>>()),
157 tensor_angular_momentum(numbers::signaling_nan<Tensor<1,dim>>()),
158 tensor_rotation(numbers::signaling_nan<Tensor<1,dim>>())
181 "Nonlinear solver failed to converge in the prescribed number of steps. " 182 "Consider changing `Max nonlinear iterations` or `Nonlinear solver failure " 209 Simulator (
const MPI_Comm mpi_communicator,
210 ParameterHandler &prm);
335 void setup_introspection ();
348 void set_initial_temperature_and_compositional_fields ();
365 void compute_initial_pressure_field ();
372 void compute_initial_velocity_boundary_constraints (AffineConstraints<double> &constraints);
379 void compute_current_velocity_boundary_constraints (AffineConstraints<double> &constraints);
393 void compute_current_constraints ();
408 double compute_pressure_scaling_factor ()
const;
421 void start_timestep ();
430 void solve_timestep ();
443 void solve_single_advection_single_stokes ();
457 void solve_no_advection_iterated_stokes ();
469 void solve_no_advection_single_stokes ();
482 void solve_first_timestep_only_single_stokes ();
496 void solve_iterated_advection_and_stokes ();
510 void solve_single_advection_iterated_stokes ();
524 void solve_no_advection_iterated_defect_correction_stokes ();
538 void solve_single_advection_iterated_defect_correction_stokes ();
553 void solve_iterated_advection_and_defect_correction_stokes ();
573 void solve_iterated_advection_and_newton_stokes (
bool use_newton_iterations);
594 void solve_single_advection_and_iterated_newton_stokes (
bool use_newton_iterations);
607 void solve_single_advection_no_stokes ();
620 void solve_no_advection_no_stokes ();
630 void build_stokes_preconditioner ();
639 void build_advection_preconditioner (
const AdvectionField &advection_field,
641 const double diagonal_strengthening);
649 void assemble_stokes_system ();
665 double assemble_and_solve_temperature (
const double &initial_residual = 0,
666 double *residual =
nullptr);
683 std::vector<double> assemble_and_solve_composition (
const std::vector<double> &initial_residual = {},
684 std::vector<double> *residual =
nullptr);
700 double assemble_and_solve_stokes (
const double &initial_nonlinear_residual = 0,
701 double *nonlinear_residual =
nullptr);
716 const bool use_picard);
725 void assemble_advection_system (
const AdvectionField &advection_field);
743 void interpolate_particle_properties (
const std::vector<AdvectionField> &advection_fields);
821 std::pair<double,double>
853 void refine_mesh (
const unsigned int max_grid_level);
870 unsigned int determine_last_good_snapshot()
const;
883 void create_snapshot();
898 void resume_from_snapshot();
906 template <
class Archive>
907 void serialize (Archive &ar,
const unsigned int version);
924 Table<2,DoFTools::Coupling>
925 setup_system_matrix_coupling ()
const;
934 void setup_system_matrix (
const std::vector<IndexSet> &system_partitioning);
947 void setup_system_preconditioner (
const std::vector<IndexSet> &system_partitioning);
980 void set_assemblers ();
992 void set_advection_assemblers ();
1003 void set_stokes_assemblers ();
1011 void assemble_stokes_preconditioner ();
1021 local_assemble_stokes_preconditioner (
const typename DoFHandler<dim>::active_cell_iterator &cell,
1043 local_assemble_stokes_system (
const typename DoFHandler<dim>::active_cell_iterator &cell,
1065 local_assemble_advection_face_terms(
const AdvectionField &advection_field,
1066 const typename DoFHandler<dim>::active_cell_iterator &cell,
1077 local_assemble_advection_system (
const AdvectionField &advection_field,
1078 const Vector<double> &viscosity_per_cell,
1079 const typename DoFHandler<dim>::active_cell_iterator &cell,
1091 copy_local_to_global_advection_system (
const AdvectionField &advection_field,
1129 template <
typename T>
1130 void get_artificial_viscosity (Vector<T> &viscosity_per_cell,
1132 const bool skip_interior_cells =
false)
const;
1207 void denormalize_pressure(
const double pressure_adjustment,
1218 void apply_limiter_to_dg_solutions (
const AdvectionField &advection_field);
1236 void compute_unique_advection_support_points (
const std::vector<AdvectionField> &advection_fields,
1237 std::vector<Point<dim>> &unique_support_points,
1238 std::vector<std::vector<unsigned int>> &support_point_index_by_field)
const;
1262 void compute_reactions ();
1273 void update_solution_vectors_with_reaction_results (
const unsigned int block_index,
1287 void initialize_current_linearization_point ();
1308 void interpolate_material_output_into_advection_field (
const std::vector<AdvectionField> &adv_field);
1318 void interpolate_onto_velocity_system(
const TensorFunction<1,dim> &func,
1336 void setup_nullspace_constraints(AffineConstraints<double> &constraints);
1369 compute_net_angular_momentum(
const bool use_constant_density,
1371 const bool limit_to_top_faces =
false)
const;
1389 void remove_net_angular_momentum(
const bool use_constant_density,
1392 const bool limit_to_top_faces =
false)
const;
1401 void replace_outflow_boundary_ids(
const unsigned int boundary_id_offset,
1402 const bool is_composition,
1403 const unsigned int composition_index);
1412 void restore_outflow_boundary_ids(
const unsigned int boundary_id_offset);
1427 void remove_net_linear_momentum(
const bool use_constant_density,
1452 double get_entropy_variation (
const double average_field,
1463 std::pair<double,double>
1464 get_extrapolated_advection_field_range (
const AdvectionField &advection_field)
const;
1474 void exchange_refinement_flags();
1485 void maybe_write_timing_output ()
const;
1494 bool maybe_write_checkpoint (
const time_t last_checkpoint_time,
1495 const bool force_writing_checkpoint);
1513 bool maybe_do_initial_refinement (
const unsigned int max_refinement_level);
1523 void maybe_refine_mesh (
const double new_time_step,
1524 unsigned int &max_refinement_level);
1533 void advance_time (
const double step_size);
1544 const double global_u_infty,
1545 const double global_field_variation,
1546 const double average_field,
1547 const double global_entropy_variation,
1548 const double cell_diameter,
1561 const double average_field,
1563 double &max_residual,
1564 double &max_velocity,
1565 double &max_density,
1566 double &max_specific_heat,
1567 double &conductivity)
const;
1584 stokes_matrix_depends_on_solution ()
const;
1596 stokes_A_block_is_symmetric ()
const;
1611 check_consistency_of_formulation ();
1619 select_default_solver_and_averaging ();
1632 check_consistency_of_boundary_conditions ()
const;
1638 compute_initial_newton_residual ();
1649 compute_Eisenstat_Walker_linear_tolerance(
const bool EisenstatWalkerChoiceOne,
1650 const double maximum_linear_stokes_solver_tolerance,
1651 const double linear_stokes_solver_tolerance,
1652 const double stokes_residual,
1653 const double newton_residual,
1654 const double newton_residual_old);
1665 void output_statistics();
1679 compute_initial_stokes_residual();
1730 using TeeDevice = boost::iostreams::tee_device<std::ostream, std::ofstream>;
1833 #ifdef ASPECT_WITH_WORLD_BUILDER 1845 std::shared_ptr<WorldBuilder::World> world_builder;
2042 friend class boost::serialization::access;
The NullspaceRemoval struct.
unsigned int nonlinear_iteration
BoundaryVelocity::Manager< dim > boundary_velocity_manager
const std::unique_ptr< AdiabaticConditions::Interface< dim > > adiabatic_conditions
std::shared_ptr< InitialComposition::Manager< dim > > initial_composition_manager
void write_plugin_graph(std::ostream &output_stream)
parallel::distributed::Triangulation< dim > triangulation
TimerOutput computing_timer
BoundaryTemperature::Manager< dim > boundary_temperature_manager
boost::iostreams::tee_device< std::ostream, std::ofstream > TeeDevice
std::unique_ptr< StokesSolver::Direct< dim > > stokes_direct
::TrilinosWrappers::MPI::BlockVector BlockVector
double scalar_moment_of_inertia
Tensor< 1, dim > tensor_rotation
Tensor< 1, dim > tensor_angular_momentum
AffineConstraints< double > constraints
bool assemble_newton_stokes_system
double scalar_angular_momentum
LinearAlgebra::BlockVector old_solution
const IntermediaryConstructorAction post_geometry_model_creation_action
std::unique_ptr< VolumeOfFluidHandler< dim > > volume_of_fluid_handler
bool rebuild_stokes_matrix
Parameters< dim > parameters
SymmetricTensor< 2, dim > tensor_moment_of_inertia
AffineConstraints< double > current_constraints
std::unique_ptr< MeshDeformation::MeshDeformationHandler< dim > > mesh_deformation
LinearAlgebra::BlockVector system_rhs
std::pair< double, double > stokes_residuals
std::vector< Particle::Manager< dim > > particle_managers
std::thread output_statistics_thread
LateralAveraging< dim > lateral_averaging
DeclExceptionMsg(ExcNonlinearSolverNoConvergence, "Nonlinear solver failed to converge in the prescribed number of steps. " "Consider changing `Max nonlinear iterations` or `Nonlinear solver failure " "strategy`.")
TeeStream iostream_tee_stream
::TrilinosWrappers::MPI::Vector Vector
MeshRefinement::Manager< dim > mesh_refinement_manager
std::size_t statistics_last_write_size
const std::unique_ptr< MaterialModel::Interface< dim > > material_model
unsigned int nonlinear_solver_failures
double last_pressure_normalization_adjustment
::TrilinosWrappers::BlockSparseMatrix BlockSparseMatrix
double total_walltime_until_last_snapshot
void declare_parameters(ParameterHandler &prm)
const std::unique_ptr< GeometryModel::Interface< dim > > geometry_model
LinearAlgebra::BlockVector operator_split_reaction_vector
DoFHandler< dim > dof_handler
std::unique_ptr< StokesMatrixFreeHandler< dim > > stokes_matrix_free
unsigned int last_checkpoint_id
std::unique_ptr< Assemblers::Manager< dim > > assemblers
bool rebuild_sparsity_and_matrices
const IntermediaryConstructorAction post_signal_creation
double global_Omega_diameter
unsigned int timestep_number
bool assemble_newton_stokes_matrix
unsigned int pre_refinement_step
TimeStepping::Manager< dim > time_stepping_manager
std::size_t statistics_last_hash
Introspection< dim > introspection
MPI_Comm mpi_communicator
typename Parameters< dim >::NullspaceRemoval NullspaceRemoval
LinearAlgebra::BlockSparseMatrix system_preconditioner_matrix
bool do_pressure_rhs_compatibility_modification
LinearAlgebra::BlockVector current_linearization_point
BoundaryConvectiveHeating::Manager< dim > boundary_convective_heating_manager
LinearAlgebra::BlockVector old_old_solution
BoundaryTraction::Manager< dim > boundary_traction_manager
LinearAlgebra::BlockVector pressure_shape_function_integrals
std::ofstream log_file_stream
std::unique_ptr< LinearAlgebra::PreconditionAMG > Amg_preconditioner
std::shared_ptr< InitialTemperature::Manager< dim > > initial_temperature_manager
TeeDevice iostream_tee_device
boost::iostreams::stream< TeeDevice > TeeStream
LinearAlgebra::BlockSparseMatrix system_matrix
std::unique_ptr< MeltHandler< dim > > melt_handler
LinearAlgebra::BlockVector inverse_lumped_mass_matrix
bool rebuild_stokes_preconditioner
std::unique_ptr< Mapping< dim > > mapping
const std::unique_ptr< PrescribedStokesSolution::Interface< dim > > prescribed_stokes_solution
double switch_initial_residual
HeatingModel::Manager< dim > heating_model_manager
SimulatorSignals< dim > signals
const std::unique_ptr< BoundaryHeatFlux::Interface< dim > > boundary_heat_flux
const FESystem< dim > finite_element
bool simulator_is_past_initialization
const std::unique_ptr< InitialTopographyModel::Interface< dim > > initial_topography_model
LinearAlgebra::BlockVector solution
std::unique_ptr< NewtonHandler< dim > > newton_handler
typename Parameters< dim >::NonlinearSolver NonlinearSolver
Postprocess::Manager< dim > postprocess_manager
const std::unique_ptr< GravityModel::Interface< dim > > gravity_model
double newton_residual_for_derivative_scaling_factor
BoundaryComposition::Manager< dim > boundary_composition_manager
std::unique_ptr< LinearAlgebra::PreconditionBase > Mp_preconditioner
::TrilinosWrappers::PreconditionILU PreconditionILU