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

Public Member Functions

 CrystalPreferredOrientation ()
 
 ~CrystalPreferredOrientation ()
 
void initialize () override
 
std::pair< std::string, std::string > execute (TableHandler &statistics) override
 
std::list< std::string > required_other_postprocessors () const override
 
void parse_parameters (ParameterHandler &prm) override
 
- Public Member Functions inherited from aspect::Postprocess::Interface< dim >
virtual ~Interface ()=default
 
virtual void update ()
 
virtual void save (std::map< std::string, std::string > &status_strings) const
 
virtual void load (const std::map< std::string, std::string > &status_strings)
 
- 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 Types

enum  Output { Output::VolumeFraction, Output::RotationMatrix, Output::EulerAngles, Output::not_found }
 

Private Member Functions

Output string_to_output_enum (const std::string &string)
 
void set_last_output_time (const double current_time)
 

Static Private Member Functions

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

Private Attributes

double end_time
 
std::mt19937 random_number_generator
 
unsigned int random_number_seed
 
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_master
 
std::vector< std::pair< unsigned int, Output > > write_raw_cpo
 
bool compute_raw_euler_angles
 
std::thread background_thread_content_raw
 
std::vector< std::pair< unsigned int, Output > > write_draw_volume_weighted_cpo
 
bool compute_weighted_rotation_matrix
 
std::thread background_thread_content_draw_volume_weighting
 
bool compress_cpo_data_files
 

Detailed Description

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

A Postprocessor that writes out CPO specific particle data. It can write out the CPO data as it is stored (raw) and/or as a random draw volume weighted representation. The latter one is recommended for plotting against real data. For both representations the specific output fields and their order can be set.

Definition at line 45 of file crystal_preferred_orientation.h.

Member Enumeration Documentation

§ Output

template<int dim>
enum aspect::Postprocess::CrystalPreferredOrientation::Output
strongprivate

Enums specifying what information to write:

VolumeFraction: Write the volume fraction RotationMatrix: Write the whole rotation matrix of a grain (in deal.ii order) EulerAngles: Write out the Z-X-Z Euler angle of a grain not_found: Error enum, the requested output was not found

Enumerator
VolumeFraction 
RotationMatrix 
EulerAngles 
not_found 

Definition at line 109 of file crystal_preferred_orientation.h.

Constructor & Destructor Documentation

§ CrystalPreferredOrientation()

Constructor.

§ ~CrystalPreferredOrientation()

Destructor.

Member Function Documentation

§ initialize()

template<int dim>
void aspect::Postprocess::CrystalPreferredOrientation< dim >::initialize ( )
overridevirtual

Initialize function.

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

§ execute()

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

A Postprocessor that writes out CPO specific particle data. It can write out the CPO data as it is stored (raw) and/or as a random draw volume weighted representation. The latter one is recommended for plotting against real data. For both representations the specific output fields and their order can be set.

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

§ required_other_postprocessors()

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

This function ensures that the particle postprocessor is run before this postprocessor.

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

§ declare_parameters()

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

Declare the parameters this class takes through input files.

§ parse_parameters()

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

Read the parameters this class declares from the parameter file.

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

§ string_to_output_enum()

template<int dim>
Output aspect::Postprocess::CrystalPreferredOrientation< dim >::string_to_output_enum ( const std::string &  string)
private

Converts a string to an Output enum.

§ set_last_output_time()

template<int dim>
void aspect::Postprocess::CrystalPreferredOrientation< 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::CrystalPreferredOrientation< dim >::writer ( const std::string &  filename,
const std::string &  temporary_filename,
const std::string &  file_contents,
const bool  compress_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.

Member Data Documentation

§ end_time

template<int dim>
double aspect::Postprocess::CrystalPreferredOrientation< dim >::end_time
private

Stores the simulation end time, so that it always produces output at the last timestep.

Definition at line 99 of file crystal_preferred_orientation.h.

§ random_number_generator

template<int dim>
std::mt19937 aspect::Postprocess::CrystalPreferredOrientation< dim >::random_number_generator
mutableprivate

Random number generator used for random draw volume weighting.

Definition at line 122 of file crystal_preferred_orientation.h.

§ random_number_seed

template<int dim>
unsigned int aspect::Postprocess::CrystalPreferredOrientation< dim >::random_number_seed
private

Random number generator seed used to initialize the random number generator.

Definition at line 127 of file crystal_preferred_orientation.h.

§ output_interval

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

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

Definition at line 133 of file crystal_preferred_orientation.h.

§ last_output_time

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

Records time for next output to occur

Definition at line 138 of file crystal_preferred_orientation.h.

§ output_file_number

template<int dim>
unsigned int aspect::Postprocess::CrystalPreferredOrientation< 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 154 of file crystal_preferred_orientation.h.

§ output_formats

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

Graphical output format.

Definition at line 159 of file crystal_preferred_orientation.h.

§ times_and_pvtu_file_names

template<int dim>
std::vector<std::pair<double,std::string> > aspect::Postprocess::CrystalPreferredOrientation< 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 169 of file crystal_preferred_orientation.h.

§ times_and_vtu_file_names

template<int dim>
std::vector<std::pair<double,std::vector<std::string> > > aspect::Postprocess::CrystalPreferredOrientation< 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 176 of file crystal_preferred_orientation.h.

§ output_file_names_by_timestep

template<int dim>
std::vector<std::vector<std::string> > aspect::Postprocess::CrystalPreferredOrientation< 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 183 of file crystal_preferred_orientation.h.

§ xdmf_entries

template<int dim>
std::vector<XDMFEntry> aspect::Postprocess::CrystalPreferredOrientation< 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 191 of file crystal_preferred_orientation.h.

§ group_files

template<int dim>
unsigned int aspect::Postprocess::CrystalPreferredOrientation< 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 199 of file crystal_preferred_orientation.h.

§ temporary_output_location

template<int dim>
std::string aspect::Postprocess::CrystalPreferredOrientation< 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 208 of file crystal_preferred_orientation.h.

§ write_in_background_thread

template<int dim>
bool aspect::Postprocess::CrystalPreferredOrientation< 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 216 of file crystal_preferred_orientation.h.

§ background_thread_master

template<int dim>
std::thread aspect::Postprocess::CrystalPreferredOrientation< dim >::background_thread_master
private

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

Definition at line 222 of file crystal_preferred_orientation.h.

§ write_raw_cpo

template<int dim>
std::vector<std::pair<unsigned int,Output> > aspect::Postprocess::CrystalPreferredOrientation< dim >::write_raw_cpo
private

What "raw" CPO data to write out.

Definition at line 227 of file crystal_preferred_orientation.h.

§ compute_raw_euler_angles

template<int dim>
bool aspect::Postprocess::CrystalPreferredOrientation< dim >::compute_raw_euler_angles
private

Whether computing raw Euler angles is needed.

Definition at line 232 of file crystal_preferred_orientation.h.

§ background_thread_content_raw

template<int dim>
std::thread aspect::Postprocess::CrystalPreferredOrientation< dim >::background_thread_content_raw
private

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

Definition at line 238 of file crystal_preferred_orientation.h.

§ write_draw_volume_weighted_cpo

template<int dim>
std::vector<std::pair<unsigned int,Output> > aspect::Postprocess::CrystalPreferredOrientation< dim >::write_draw_volume_weighted_cpo
private

What "draw volume weighted" CPO data to write out. Draw volume weighted means the grain properties (size and rotation) are put in a list sorted based by volume and picked randomly, with large volumes having a higher chance of being picked.

Definition at line 245 of file crystal_preferred_orientation.h.

§ compute_weighted_rotation_matrix

template<int dim>
bool aspect::Postprocess::CrystalPreferredOrientation< dim >::compute_weighted_rotation_matrix
private

Whether computing weighted A matrix is needed.

Definition at line 250 of file crystal_preferred_orientation.h.

§ background_thread_content_draw_volume_weighting

template<int dim>
std::thread aspect::Postprocess::CrystalPreferredOrientation< dim >::background_thread_content_draw_volume_weighting
private

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

Definition at line 256 of file crystal_preferred_orientation.h.

§ compress_cpo_data_files

template<int dim>
bool aspect::Postprocess::CrystalPreferredOrientation< dim >::compress_cpo_data_files
private

Whether to compress the raw and weighed cpo data output files with zlib.

Definition at line 261 of file crystal_preferred_orientation.h.


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