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

Public Member Functions

DEAL_II_DEPRECATED std::vector< std::vector< double > > get_averages (const unsigned int n_slices, const std::vector< std::string > &property_names) const
 
std::vector< std::vector< double > > compute_lateral_averages (const unsigned int n_slices, const std::vector< std::string > &property_names) const
 
std::vector< std::vector< double > > compute_lateral_averages (const std::vector< double > &depth_bounds, const std::vector< std::string > &property_names) const
 
std::vector< std::vector< double > > compute_lateral_averages (const std::vector< double > &depth_bounds, std::vector< std::unique_ptr< internal::FunctorBase< dim >>> &functors) const
 
void get_temperature_averages (std::vector< double > &values) const
 
void get_composition_averages (const unsigned int composition_index, std::vector< double > &values) const
 
void get_viscosity_averages (std::vector< double > &values) const
 
void get_log_viscosity_averages (std::vector< double > &values) const
 
void get_velocity_magnitude_averages (std::vector< double > &values) const
 
void get_sinking_velocity_averages (std::vector< double > &values) const
 
void get_rising_velocity_averages (std::vector< double > &values) const
 
void get_Vs_averages (std::vector< double > &values) const
 
void get_Vp_averages (std::vector< double > &values) const
 
void get_vertical_heat_flux_averages (std::vector< double > &values) const
 
void get_vertical_mass_flux_averages (std::vector< double > &values) 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
 

Additional Inherited Members

- 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)
 

Detailed Description

template<int dim>
class aspect::LateralAveraging< dim >

LateralAveraging is a class that performs various averaging operations on the solution. The functions of this class take a vector as an argument. The model is divided into depth slices where the number of slices is the length of the output vector. Each function averages a specific quantity (as specified by their name), and that quantity is averaged laterally for each depth slice.

Plugins may access the LateralAveraging plugin through the SimulatorAccess function get_lateral_averaging(), and then query that for the desired averaged quantity.

Definition at line 110 of file lateral_averaging.h.

Member Function Documentation

§ get_averages()

template<int dim>
DEAL_II_DEPRECATED std::vector<std::vector<double> > aspect::LateralAveraging< dim >::get_averages ( const unsigned int  n_slices,
const std::vector< std::string > &  property_names 
) const
Deprecated:
: This function is deprecated and only maintained for backward compatibility. Use the function compute_lateral_averages() with the same arguments instead.

§ compute_lateral_averages() [1/3]

template<int dim>
std::vector<std::vector<double> > aspect::LateralAveraging< dim >::compute_lateral_averages ( const unsigned int  n_slices,
const std::vector< std::string > &  property_names 
) const

Return a depth profile of lateral averages of the selected property_names. This function is a convenience interface for the other functions of the same name and is more efficient for multiple properties than calling multiple of the single property functions sequentially. This version of the function creates n_slices equidistant depth slices to compute averages in.

Parameters
n_slicesThe number of equidistant depth slices to perform the averaging in.
property_namesNames of the available quantities to average. Check the implementation of this function for available names.
Returns
The output vector of laterally averaged values. Each vector has the same size of n_slices, and there are as many vectors returned as names in property_names.

§ compute_lateral_averages() [2/3]

template<int dim>
std::vector<std::vector<double> > aspect::LateralAveraging< dim >::compute_lateral_averages ( const std::vector< double > &  depth_bounds,
const std::vector< std::string > &  property_names 
) const

Return a depth profile of lateral averages of the selected property_names. This function is a convenience interface for the other functions of the same name and is more efficient for multiple properties than calling multiple of the single property functions sequentially. This version of the function uses depth_bounds to determine the upper and lower boundary of each depth slice, and returns one averaged value per slice (= one less than the number of depth bounds).

Parameters
depth_boundsThe boundaries of the depth slices to compute averages in. It is expected to be a vector of monotonically increasing depth values with consecutive entries representing the minimum and maximum depth of each depth slice. All depths smaller than entry 0 or larger than the last entry are ignored, and depths between entries 0 and 1 fall into depth slice 0, and so on.
property_namesNames of the available quantities to average. Check the implementation of this function for available names.
Returns
The output vector of laterally averaged values. The length of each vector is one less than the number of depth_bounds, and there are as many vectors returned as names in property_names.

§ compute_lateral_averages() [3/3]

template<int dim>
std::vector<std::vector<double> > aspect::LateralAveraging< dim >::compute_lateral_averages ( const std::vector< double > &  depth_bounds,
std::vector< std::unique_ptr< internal::FunctorBase< dim >>> &  functors 
) const

Return a depth profile of lateral averages. This function is the low-level implementation for the other functions of the same name and is more efficient for multiple properties than calling multiple of the single property functions sequentially. This version of the function uses depth_bounds to determine the upper and lower boundary of each depth slice, and returns one averaged value per slice (= one less than the number of depth bounds). The vector of functors functors must contain one or more objects of classes that are derived from FunctorBase and are used to fill the values vectors. All of the other functions in this class use this low-level implementation for the actual computation. Using this function allows to provide user-defined functors to average properties not already implemented (e.g. to average additional material model outputs).

Parameters
depth_boundsThe boundaries of the depth slices to compute averages in. It is expected to be a vector of monotonically increasing depth values with consecutive entries representing the minimum and maximum depth of each depth slice. All depths smaller than entry 0 or larger than the last entry are ignored, and depths between entries 0 and 1 fall into depth slice 0, and so on.
functorsInstances of a class derived from FunctorBase that are used to compute the averaged properties.
Returns
The output vectors of depth averaged values. The function returns one vector of doubles per property and uses depth_bounds to determine the depth extents of each slice. Each returned vector has the same size (one entry less than the number of depth_bounds).

§ get_temperature_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_temperature_averages ( std::vector< double > &  values) const

Fill the argument with a set of lateral averages of the current temperature field. The function fills a vector that contains average field values over slices of the domain of same depth.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_composition_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_composition_averages ( const unsigned int  composition_index,
std::vector< double > &  values 
) const

Fill the argument with a set of lateral averages of the current compositional fields.

Parameters
composition_indexThe index of the compositional field whose matrix we want to assemble (0 <= composition_index < number of compositional fields in this problem).
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_viscosity_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_viscosity_averages ( std::vector< double > &  values) const

Compute a lateral average of the current viscosity.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_log_viscosity_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_log_viscosity_averages ( std::vector< double > &  values) const

Compute a lateral average of the log10 of the current viscosity.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_velocity_magnitude_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_velocity_magnitude_averages ( std::vector< double > &  values) const

Compute a lateral average of the current velocity magnitude.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_sinking_velocity_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_sinking_velocity_averages ( std::vector< double > &  values) const

Compute a lateral average of the current sinking velocity.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_rising_velocity_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_rising_velocity_averages ( std::vector< double > &  values) const

Compute a lateral average of the current rising velocity.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_Vs_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_Vs_averages ( std::vector< double > &  values) const

Compute a lateral average of the seismic shear wave speed: Vs.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_Vp_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_Vp_averages ( std::vector< double > &  values) const

Compute a lateral average of the seismic pressure wave speed: Vp.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_vertical_heat_flux_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_vertical_heat_flux_averages ( std::vector< double > &  values) const

Compute a lateral average of the vertical heat flux, with the sign convention of positive heat flux when it flows upwards.

Parameters
valuesThe output vector of laterally averaged values. The function takes the pre-existing size of this vector as the number of depth slices.

§ get_vertical_mass_flux_averages()

template<int dim>
void aspect::LateralAveraging< dim >::get_vertical_mass_flux_averages ( std::vector< double > &  values) const

Compute a lateral average of the vertical mass flux. Note that while get_vertical_heat_flux_averages() computes the average vertical heat flux (positive or negative), this function computes the average of the total mass flux through a certain layer (both down- and upward motion counted positively).

Parameters
valuesThe output vector of laterally averaged values.

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