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

Classes

struct  OutputHistory
 

Public Member Functions

 Visualization ()
 
std::pair< std::string, std::string > execute (TableHandler &statistics) override
 
void update () override
 
std::list< std::string > required_other_postprocessors () const override
 
void parse_parameters (ParameterHandler &prm) 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)
 
bool output_pointwise_stress_and_strain () const
 
 DeclException1 (ExcPostprocessorNameNotFound, std::string,<< "Could not find entry <"<< arg1<< "> among the names of registered postprocessors.")
 
- Public Member Functions inherited from aspect::Postprocess::Interface< dim >
virtual ~Interface ()=default
 
virtual void initialize ()
 
- 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 register_visualization_postprocessor (const std::string &name, const std::string &description, void(*declare_parameters_function)(ParameterHandler &), std::unique_ptr< VisualizationPostprocessors::Interface< dim >>(*factory_function)())
 
static void declare_parameters (ParameterHandler &prm)
 
static void write_plugin_graph (std::ostream &output_stream)
 
- 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 mesh_changed_signal ()
 
template<typename DataOutType >
void write_master_files (const DataOutType &data_out, const std::string &solution_file_prefix, const std::vector< std::string > &filenames, OutputHistory &output_history) const
 
template<typename DataOutType >
std::string write_data_out_data (DataOutType &data_out, OutputHistory &output_history, const std::map< std::string, std::string > &visualization_field_names_and_units) const
 

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 maximum_timesteps_between_outputs
 
unsigned int last_output_timestep
 
unsigned int output_file_number
 
std::string output_format
 
unsigned int group_files
 
std::string temporary_output_location
 
bool interpolate_output
 
bool filter_output
 
bool pointwise_stress_and_strain
 
bool write_higher_order_output
 
bool output_mesh_velocity
 
bool output_mesh_displacement
 
bool output_undeformed_mesh
 
bool write_in_background_thread
 
std::list< std::unique_ptr< VisualizationPostprocessors::Interface< dim > > > postprocessors
 
OutputHistory cell_output_history
 
OutputHistory face_output_history
 

Detailed Description

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

A postprocessor that generates graphical output in periodic intervals or every time step. The time interval between generating graphical output is obtained from the parameter file.

While this class acts as a plugin, i.e. as a postprocessor that can be registered with the postprocessing manager class aspect::Postprocess::Manager, this class at the same time also acts as a manager for plugins itself, namely for classes derived from the VisualizationPostprocessors::Interface class that are used to output different aspects of the solution, such as for example a computed seismic wave speed from temperature, velocity and pressure.

Definition at line 360 of file visualization.h.

Constructor & Destructor Documentation

§ Visualization()

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

Constructor.

Member Function Documentation

§ execute()

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

Generate graphical output from the current solution.

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

§ update()

template<int dim>
void aspect::Postprocess::Visualization< dim >::update ( )
overridevirtual

Update any temporary information needed by the visualization postprocessor.

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

§ register_visualization_postprocessor()

template<int dim>
static void aspect::Postprocess::Visualization< dim >::register_visualization_postprocessor ( const std::string &  name,
const std::string &  description,
void(*)(ParameterHandler &)  declare_parameters_function,
std::unique_ptr< VisualizationPostprocessors::Interface< dim >>(*)()  factory_function 
)
static

A function that is used to register visualization postprocessor objects in such a way that the Manager can deal with all of them without having to know them by name. This allows the files in which individual postprocessors are implement to register these postprocessors, rather than also having to modify the Manage class by adding the new postprocessor class.

Parameters
nameThe name under which this visualization postprocessor is to be called in parameter files.
descriptionA text description of what this visualization plugin does and that will be listed in the documentation of the parameter file.
declare_parameters_functionA pointer to a function that declares the parameters for this postprocessor.
factory_functionA pointer to a function that creates such a postprocessor object and returns a pointer to it.

§ required_other_postprocessors()

template<int dim>
std::list<std::string> aspect::Postprocess::Visualization< dim >::required_other_postprocessors ( ) const
overridevirtual

A function that is used to indicate to the postprocessor manager which other postprocessor(s) the current one depends upon.

For the current class, we simply loop over all of the visualization postprocessors and collect what they want.

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

§ declare_parameters()

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

Declare the parameters this class takes through input files.

§ parse_parameters()

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

Read the parameters this class declares from the parameter file.

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

§ save()

template<int dim>
void aspect::Postprocess::Visualization< 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::Visualization< 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::Visualization< 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.

§ write_plugin_graph()

template<int dim>
static void aspect::Postprocess::Visualization< dim >::write_plugin_graph ( std::ostream &  output_stream)
static

For the current plugin subsystem, write a connection graph of all of the plugins we know about, in the format that the programs dot and neato understand. This allows for a visualization of how all of the plugins that ASPECT knows about are interconnected, and connect to other parts of the ASPECT code.

Parameters
output_streamThe stream to write the output to.

§ output_pointwise_stress_and_strain()

template<int dim>
bool aspect::Postprocess::Visualization< dim >::output_pointwise_stress_and_strain ( ) const

Return the value of the parameter pointwise_stress_and_strain that is controlled by the parameter "Point-wise stress and strain".

§ DeclException1()

template<int dim>
aspect::Postprocess::Visualization< dim >::DeclException1 ( ExcPostprocessorNameNotFound  ,
std::string   
)

Exception.

§ set_last_output_time()

template<int dim>
void aspect::Postprocess::Visualization< 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.

§ mesh_changed_signal()

template<int dim>
void aspect::Postprocess::Visualization< dim >::mesh_changed_signal ( )
private

Record that the mesh changed. This helps some output writers avoid writing the same mesh multiple times.

§ writer()

template<int dim>
static void aspect::Postprocess::Visualization< 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.

§ write_master_files()

template<int dim>
template<typename DataOutType >
void aspect::Postprocess::Visualization< dim >::write_master_files ( const DataOutType &  data_out,
const std::string &  solution_file_prefix,
const std::vector< std::string > &  filenames,
OutputHistory output_history 
) const
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.
output_historyThe OutputHistory object to fill.

§ write_data_out_data()

template<int dim>
template<typename DataOutType >
std::string aspect::Postprocess::Visualization< dim >::write_data_out_data ( DataOutType &  data_out,
OutputHistory output_history,
const std::map< std::string, std::string > &  visualization_field_names_and_units 
) const
private

A function, called from the execute() function, that takes a DataOut object, does some preliminary work with it, and then writes the result out to files via the writer() function (in the case of VTU output) or through the XDMF facilities.

The function returns the base name of the output files produced, which can then be used for the statistics file and screen output.

Member Data Documentation

§ output_interval

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

Interval between the generation of graphical output. This parameter is read from the input file and consequently is not part of the state that needs to be saved and restored.

Definition at line 479 of file visualization.h.

§ last_output_time

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

A time (in seconds) at which the last graphical output was supposed to be produced. Used to check for the next necessary output time.

Definition at line 485 of file visualization.h.

§ maximum_timesteps_between_outputs

template<int dim>
unsigned int aspect::Postprocess::Visualization< dim >::maximum_timesteps_between_outputs
private

Maximum number of steps between the generation of graphical output. This parameter is read from the input file and consequently is not part of the state that needs to be saved and restored.

Definition at line 493 of file visualization.h.

§ last_output_timestep

template<int dim>
unsigned int aspect::Postprocess::Visualization< dim >::last_output_timestep
private

Timestep at which the last graphical output was produced Used to check for the next necessary output time.

Definition at line 499 of file visualization.h.

§ output_file_number

template<int dim>
unsigned int aspect::Postprocess::Visualization< 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 505 of file visualization.h.

§ output_format

template<int dim>
std::string aspect::Postprocess::Visualization< dim >::output_format
private

Graphical output format.

Definition at line 510 of file visualization.h.

§ group_files

template<int dim>
unsigned int aspect::Postprocess::Visualization< 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 518 of file visualization.h.

§ temporary_output_location

template<int dim>
std::string aspect::Postprocess::Visualization< 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 527 of file visualization.h.

§ interpolate_output

template<int dim>
bool aspect::Postprocess::Visualization< dim >::interpolate_output
private

deal.II offers the possibility to linearly interpolate output fields of higher order elements to a finer resolution. This somewhat compensates the fact that most visualization software only offers linear interpolation between grid points and therefore the output file is a very coarse representation of the actual solution field. Activating this option increases the spatial resolution in each dimension by a factor equal to the polynomial degree used for the velocity finite element (usually 2).

Definition at line 539 of file visualization.h.

§ filter_output

template<int dim>
bool aspect::Postprocess::Visualization< dim >::filter_output
private

deal.II offers the possibility to filter duplicate vertices in HDF5 output files. This merges the vertices of adjacent cells and therefore saves disk space, but misrepresents discontinuous output properties. Activating this function reduces the disk space by about a factor of \(2^{dim}\) for hdf5 output.

Definition at line 548 of file visualization.h.

§ pointwise_stress_and_strain

template<int dim>
bool aspect::Postprocess::Visualization< dim >::pointwise_stress_and_strain
private

If true, return quantities related to stresses and strain with point-wise values. Otherwise the values will be averaged on each cell.

Definition at line 555 of file visualization.h.

§ write_higher_order_output

template<int dim>
bool aspect::Postprocess::Visualization< dim >::write_higher_order_output
private

deal.II offers the possibility to write vtu files with higher order representations of the output data. This means each cell will correctly show the higher order representation of the output data instead of the linear interpolation between vertices that ParaView and VisIt usually show. Note that activating this option is safe and recommended, but requires that (i) Output format'' is set tovtu'', (ii) ``Interpolate output'' is set to true, (iii) you use a sufficiently new version of Paraview or VisIt to read the files (Paraview version 5.5 or newer, and VisIt version to be determined), and (iv) you use deal.II version 9.1.0 or newer.

Definition at line 568 of file visualization.h.

§ output_mesh_velocity

template<int dim>
bool aspect::Postprocess::Visualization< dim >::output_mesh_velocity
private

For mesh deformation computations ASPECT uses an Arbitrary-Lagrangian-Eulerian formulation to handle deforming the domain, so the mesh has its own velocity field. This may be written as an output field by setting output_mesh_velocity to true.

Definition at line 576 of file visualization.h.

§ output_mesh_displacement

template<int dim>
bool aspect::Postprocess::Visualization< dim >::output_mesh_displacement
private

For mesh deformation computations ASPECT uses an Arbitrary-Lagrangian-Eulerian formulation to handle deforming the domain, so the mesh has a field that determines the displacement from the reference configuration. This may be written as an output field by setting this flag to true.

Definition at line 585 of file visualization.h.

§ output_undeformed_mesh

template<int dim>
bool aspect::Postprocess::Visualization< dim >::output_undeformed_mesh
private

For mesh deformation computations ASPECT uses an Arbitrary-Lagrangian-Eulerian formulation to handle deforming the domain, and we output the mesh in its deformed state by default. If this flag is set to true, the mesh is written undeformed.

Definition at line 593 of file visualization.h.

§ write_in_background_thread

template<int dim>
bool aspect::Postprocess::Visualization< 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 601 of file visualization.h.

§ postprocessors

template<int dim>
std::list<std::unique_ptr<VisualizationPostprocessors::Interface<dim> > > aspect::Postprocess::Visualization< dim >::postprocessors
private

A list of postprocessor objects that have been requested in the parameter file.

Definition at line 634 of file visualization.h.

§ cell_output_history

template<int dim>
OutputHistory aspect::Postprocess::Visualization< dim >::cell_output_history
private

Information about the history of writing graphical output for cells (via DataOut).

Definition at line 712 of file visualization.h.

§ face_output_history

template<int dim>
OutputHistory aspect::Postprocess::Visualization< dim >::face_output_history
private

Information about the history of writing graphical output for faces (via DataOutFaces).

Definition at line 718 of file visualization.h.


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