ASPECT
Public Member Functions | Static Public Member Functions | Private Member Functions | Static Private Member Functions | Private Attributes | List of all members
aspect::Postprocess::Particles< dim > Class Template Reference
Inheritance diagram for aspect::Postprocess::Particles< dim >:
Inheritance graph
[legend]

Public Member Functions

 Particles ()
 
 ~Particles () override
 
std::pair< std::string, std::string > execute (TableHandler &statistics) override
 
void save (std::map< std::string, std::string > &status_strings) const override
 
void load (const std::map< std::string, std::string > &status_strings) override
 
template<class Archive >
void serialize (Archive &ar, const unsigned int version)
 
void parse_parameters (ParameterHandler &prm) override
 
- Public Member Functions inherited from aspect::Postprocess::Interface< dim >
virtual ~Interface ()=default
 
virtual void initialize ()
 
virtual void update ()
 
virtual std::list< std::string > required_other_postprocessors () const
 
- Public Member Functions inherited from aspect::SimulatorAccess< dim >
 SimulatorAccess ()
 
 SimulatorAccess (const Simulator< dim > &simulator_object)
 
virtual ~SimulatorAccess ()=default
 
virtual void initialize_simulator (const Simulator< dim > &simulator_object)
 
const Introspection< dim > & introspection () const
 
const Simulator< dim > & get_simulator () const
 
const Parameters< dim > & get_parameters () const
 
SimulatorSignals< dim > & get_signals () const
 
MPI_Comm get_mpi_communicator () const
 
TimerOutput & get_computing_timer () const
 
const ConditionalOStream & get_pcout () const
 
double get_time () const
 
double get_timestep () const
 
double get_old_timestep () const
 
unsigned int get_timestep_number () const
 
const TimeStepping::Manager< dim > & get_timestepping_manager () const
 
unsigned int get_nonlinear_iteration () const
 
const parallel::distributed::Triangulation< dim > & get_triangulation () const
 
double get_volume () const
 
const Mapping< dim > & get_mapping () const
 
std::string get_output_directory () const
 
bool include_adiabatic_heating () const
 
bool include_latent_heat () const
 
bool include_melt_transport () const
 
int get_stokes_velocity_degree () const
 
double get_adiabatic_surface_temperature () const
 
double get_surface_pressure () const
 
bool convert_output_to_years () const
 
unsigned int get_pre_refinement_step () const
 
unsigned int n_compositional_fields () const
 
double get_end_time () const
 
void get_refinement_criteria (Vector< float > &estimated_error_per_cell) const
 
void get_artificial_viscosity (Vector< float > &viscosity_per_cell, const bool skip_interior_cells=false) const
 
void get_artificial_viscosity_composition (Vector< float > &viscosity_per_cell, const unsigned int compositional_variable) const
 
const LinearAlgebra::BlockVectorget_current_linearization_point () const
 
const LinearAlgebra::BlockVectorget_solution () const
 
const LinearAlgebra::BlockVectorget_old_solution () const
 
const LinearAlgebra::BlockVectorget_old_old_solution () const
 
const LinearAlgebra::BlockVectorget_reaction_vector () const
 
const LinearAlgebra::BlockVectorget_mesh_velocity () const
 
const DoFHandler< dim > & get_dof_handler () const
 
const FiniteElement< dim > & get_fe () const
 
const LinearAlgebra::BlockSparseMatrixget_system_matrix () const
 
const LinearAlgebra::BlockSparseMatrixget_system_preconditioner_matrix () const
 
const MaterialModel::Interface< dim > & get_material_model () const
 
const GravityModel::Interface< dim > & get_gravity_model () const
 
const InitialTopographyModel::Interface< dim > & get_initial_topography_model () const
 
const GeometryModel::Interface< dim > & get_geometry_model () const
 
const AdiabaticConditions::Interface< dim > & get_adiabatic_conditions () const
 
bool has_boundary_temperature () const
 
const BoundaryTemperature::Manager< dim > & get_boundary_temperature_manager () const
 
const BoundaryHeatFlux::Interface< dim > & get_boundary_heat_flux () const
 
bool has_boundary_composition () const
 
const BoundaryComposition::Manager< dim > & get_boundary_composition_manager () const
 
const BoundaryTraction::Manager< dim > & get_boundary_traction_manager () const
 
std::shared_ptr< const InitialTemperature::Manager< dim > > get_initial_temperature_manager_pointer () const
 
const InitialTemperature::Manager< dim > & get_initial_temperature_manager () const
 
std::shared_ptr< const InitialComposition::Manager< dim > > get_initial_composition_manager_pointer () const
 
const InitialComposition::Manager< dim > & get_initial_composition_manager () const
 
const std::set< types::boundary_id > & get_fixed_temperature_boundary_indicators () const
 
const std::set< types::boundary_id > & get_fixed_heat_flux_boundary_indicators () const
 
const std::set< types::boundary_id > & get_fixed_composition_boundary_indicators () const
 
const std::set< types::boundary_id > & get_mesh_deformation_boundary_indicators () const
 
const BoundaryVelocity::Manager< dim > & get_boundary_velocity_manager () const
 
const HeatingModel::Manager< dim > & get_heating_model_manager () const
 
const MeshRefinement::Manager< dim > & get_mesh_refinement_manager () const
 
const MeltHandler< dim > & get_melt_handler () const
 
const VolumeOfFluidHandler< dim > & get_volume_of_fluid_handler () const
 
const NewtonHandler< dim > & get_newton_handler () const
 
const MeshDeformation::MeshDeformationHandler< dim > & get_mesh_deformation_handler () const
 
const LateralAveraging< dim > & get_lateral_averaging () const
 
const AffineConstraints< double > & get_current_constraints () const
 
bool simulator_is_past_initialization () const
 
double get_pressure_scaling () const
 
bool pressure_rhs_needs_compatibility_modification () const
 
bool model_has_prescribed_stokes_solution () const
 
TableHandler & get_statistics_object () const
 
const Postprocess::Manager< dim > & get_postprocess_manager () const
 
const Particle::World< dim > & get_particle_world () const
 
Particle::World< dim > & get_particle_world ()
 
bool is_stokes_matrix_free ()
 
const StokesMatrixFreeHandler< dim > & get_stokes_matrix_free () const
 
RotationProperties< dim > compute_net_angular_momentum (const bool use_constant_density, const LinearAlgebra::BlockVector &solution, const bool limit_to_top_faces=false) const
 

Static Public Member Functions

static void declare_parameters (ParameterHandler &prm)
 
- Static Public Member Functions inherited from aspect::Postprocess::Interface< dim >
static void declare_parameters (ParameterHandler &prm)
 
- Static Public Member Functions inherited from aspect::SimulatorAccess< dim >
static void get_composition_values_at_q_point (const std::vector< std::vector< double >> &composition_values, const unsigned int q, std::vector< double > &composition_values_at_q_point)
 

Private Member Functions

void set_last_output_time (const double current_time)
 
void write_master_files (const internal::ParticleOutput< dim > &data_out, const std::string &solution_file_prefix, const std::vector< std::string > &filenames)
 

Static Private Member Functions

static void writer (const std::string &filename, const std::string &temporary_filename, const std::string &file_contents)
 

Private Attributes

double output_interval
 
double last_output_time
 
unsigned int output_file_number
 
std::vector< std::string > output_formats
 
std::vector< std::pair< double, std::string > > times_and_pvtu_file_names
 
std::vector< std::pair< double, std::vector< std::string > > > times_and_vtu_file_names
 
std::vector< std::vector< std::string > > output_file_names_by_timestep
 
std::vector< XDMFEntry > xdmf_entries
 
unsigned int group_files
 
std::string temporary_output_location
 
bool write_in_background_thread
 
std::thread background_thread
 
std::vector< std::string > exclude_output_properties
 

Detailed Description

template<int dim>
class aspect::Postprocess::Particles< dim >

A Postprocessor that creates particles, which follow the velocity field of the simulation. The particles can be generated and propagated in various ways and they can carry a number of constant or time-varying properties. The postprocessor can write output positions and properties of all particles at chosen intervals, although this is not mandatory. It also allows other parts of the code to query the particles for information.

Definition at line 120 of file particles.h.

Constructor & Destructor Documentation

§ Particles()

template<int dim>
aspect::Postprocess::Particles< dim >::Particles ( )

Constructor.

§ ~Particles()

template<int dim>
aspect::Postprocess::Particles< dim >::~Particles ( )
override

Destructor.

Member Function Documentation

§ execute()

template<int dim>
std::pair<std::string,std::string> aspect::Postprocess::Particles< dim >::execute ( TableHandler &  statistics)
overridevirtual

Execute this postprocessor. Derived classes will implement this function to do whatever they want to do to evaluate the solution at the current time step.

Parameters
[in,out]statisticsAn object that contains statistics that are collected throughout the simulation and that will be written to an output file at the end of each time step. Postprocessors may deposit data in these tables for later visualization or further processing.
Returns
A pair of strings that will be printed to the screen after running the postprocessor in two columns; typically the first column contains a description of what the data is and the second contains a numerical value of this data. If there is nothing to print, simply return two empty strings.

Implements aspect::Postprocess::Interface< dim >.

§ save()

template<int dim>
void aspect::Postprocess::Particles< dim >::save ( std::map< std::string, std::string > &  status_strings) const
overridevirtual

Save the state of this object.

Reimplemented from aspect::Postprocess::Interface< dim >.

§ load()

template<int dim>
void aspect::Postprocess::Particles< dim >::load ( const std::map< std::string, std::string > &  status_strings)
overridevirtual

Restore the state of the object.

Reimplemented from aspect::Postprocess::Interface< dim >.

§ serialize()

template<int dim>
template<class Archive >
void aspect::Postprocess::Particles< dim >::serialize ( Archive &  ar,
const unsigned int  version 
)

Serialize the contents of this class as far as they are not read from input parameter files.

§ declare_parameters()

template<int dim>
static void aspect::Postprocess::Particles< dim >::declare_parameters ( ParameterHandler &  prm)
static

Declare the parameters this class takes through input files.

§ parse_parameters()

template<int dim>
void aspect::Postprocess::Particles< dim >::parse_parameters ( ParameterHandler &  prm)
overridevirtual

Read the parameters this class declares from the parameter file.

Reimplemented from aspect::Postprocess::Interface< dim >.

§ set_last_output_time()

template<int dim>
void aspect::Postprocess::Particles< dim >::set_last_output_time ( const double  current_time)
private

Set the time output was supposed to be written. In the simplest case, this is the previous last output time plus the interval, but in general we'd like to ensure that it is the largest supposed output time, which is smaller than the current time, to avoid falling behind with last_output_time and having to catch up once the time step becomes larger. This is done after every output.

§ writer()

template<int dim>
static void aspect::Postprocess::Particles< dim >::writer ( const std::string &  filename,
const std::string &  temporary_filename,
const std::string &  file_contents 
)
staticprivate

A function that writes the text in the second argument to a file with the name given in the first argument. The function is run on a separate thread to allow computations to continue even though writing data is still continuing. The function takes over ownership of these arguments and deletes them at the end of its work.

§ write_master_files()

template<int dim>
void aspect::Postprocess::Particles< dim >::write_master_files ( const internal::ParticleOutput< dim > &  data_out,
const std::string &  solution_file_prefix,
const std::vector< std::string > &  filenames 
)
private

Write the various master record files. The master files are used by visualization programs to identify which of the output files in a directory, possibly one file written by each processor, belong to a single time step and/or form the different time steps of a simulation. For Paraview, this is a .pvtu file per time step and a .pvd for all time steps. For VisIt it is a .visit file per time step and one for all time steps.

Parameters
data_outThe DataOut object that was used to write the solutions.
solution_file_prefixThe stem of the filename to be written.
filenamesList of filenames for the current output from all processors.

Member Data Documentation

§ output_interval

template<int dim>
double aspect::Postprocess::Particles< dim >::output_interval
private

Interval between output (in years if appropriate simulation parameter is set, otherwise seconds)

Definition at line 187 of file particles.h.

§ last_output_time

template<int dim>
double aspect::Postprocess::Particles< dim >::last_output_time
private

Records time for next output to occur

Definition at line 192 of file particles.h.

§ output_file_number

template<int dim>
unsigned int aspect::Postprocess::Particles< dim >::output_file_number
private

Consecutively counted number indicating the how-manyth time we will create output the next time we get to it.

Definition at line 208 of file particles.h.

§ output_formats

template<int dim>
std::vector<std::string> aspect::Postprocess::Particles< dim >::output_formats
private

Graphical output format.

Definition at line 213 of file particles.h.

§ times_and_pvtu_file_names

template<int dim>
std::vector<std::pair<double,std::string> > aspect::Postprocess::Particles< dim >::times_and_pvtu_file_names
private

A list of pairs (time, pvtu_filename) that have so far been written and that we will pass to DataOutInterface::write_pvd_record to create a master file that can make the association between simulation time and corresponding file name (this is done because there is no way to store the simulation time inside the .pvtu or .vtu files).

Definition at line 223 of file particles.h.

§ times_and_vtu_file_names

template<int dim>
std::vector<std::pair<double,std::vector<std::string> > > aspect::Postprocess::Particles< dim >::times_and_vtu_file_names
private

A corresponding variable that we use for the .visit files created by DataOutInterface::write_visit_record. The second part of a pair contains all files that together form a time step.

Definition at line 230 of file particles.h.

§ output_file_names_by_timestep

template<int dim>
std::vector<std::vector<std::string> > aspect::Postprocess::Particles< dim >::output_file_names_by_timestep
private

A list of list of filenames, sorted by timestep, that correspond to what has been created as output. This is used to create a master .visit file for the entire simulation.

Definition at line 237 of file particles.h.

§ xdmf_entries

template<int dim>
std::vector<XDMFEntry> aspect::Postprocess::Particles< dim >::xdmf_entries
private

A set of data related to XDMF file sections describing the HDF5 heavy data files created. These contain things such as the dimensions and names of data written at all steps during the simulation.

Definition at line 245 of file particles.h.

§ group_files

template<int dim>
unsigned int aspect::Postprocess::Particles< dim >::group_files
private

VTU file output supports grouping files from several CPUs into one file using MPI I/O when writing on a parallel filesystem. 0 means no grouping (and no parallel I/O). 1 will generate one big file containing the whole solution.

Definition at line 253 of file particles.h.

§ temporary_output_location

template<int dim>
std::string aspect::Postprocess::Particles< dim >::temporary_output_location
private

On large clusters it can be advantageous to first write the output to a temporary file on a local file system and later move this file to a network file system. If this variable is set to a non-empty string it will be interpreted as a temporary storage location.

Definition at line 262 of file particles.h.

§ write_in_background_thread

template<int dim>
bool aspect::Postprocess::Particles< dim >::write_in_background_thread
private

File operations can potentially take a long time, blocking the progress of the rest of the model run. Setting this variable to 'true' moves this process into a background thread, while the rest of the model continues.

Definition at line 270 of file particles.h.

§ background_thread

template<int dim>
std::thread aspect::Postprocess::Particles< dim >::background_thread
private

Handle to a thread that is used to write data in the background. The writer() function runs on this background thread.

Definition at line 276 of file particles.h.

§ exclude_output_properties

template<int dim>
std::vector<std::string> aspect::Postprocess::Particles< dim >::exclude_output_properties
private

Stores the particle property fields which are excluded from output to the visualization file.

Definition at line 282 of file particles.h.


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