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

Classes

class  EllipsoidalChunkGeometry
 

Public Member Functions

virtual void initialize ()
 
virtual void create_coarse_mesh (parallel::distributed::Triangulation< dim > &coarse_grid) const
 
virtual double length_scale () const
 
virtual double depth (const Point< dim > &position) const
 
virtual double height_above_reference_surface (const Point< dim > &position) const
 
virtual Point< dim > representative_point (const double depth) const
 
virtual bool point_is_in_domain (const Point< dim > &point) const
 
virtual double maximal_depth () const
 
virtual std::set< types::boundary_idget_used_boundary_indicators () const
 
virtual std::map< std::string, types::boundary_idget_symbolic_boundary_names_map () const
 
virtual aspect::Utilities::Coordinates::CoordinateSystem natural_coordinate_system () const
 
virtual std::array< double, dim > cartesian_to_natural_coordinates (const Point< dim > &position) const
 
virtual Point< dim > natural_to_cartesian_coordinates (const std::array< double, dim > &position) const
 
virtual void parse_parameters (ParameterHandler &prm)
 
double get_radius (const Point< dim > &position) const
 
double get_semi_minor_axis_b () const
 
double get_semi_major_axis_a () const
 
double get_eccentricity () const
 
const std::vector< Point< 2 > > & get_corners () const
 
EllipsoidalChunkGeometry get_manifold () const
 
- Public Member Functions inherited from aspect::GeometryModel::Interface< dim >
virtual ~Interface ()
 
Utilities::NaturalCoordinate< dim > cartesian_to_other_coordinates (const Point< dim > &position, const Utilities::Coordinates::CoordinateSystem &coordinate_system) const
 
types::boundary_id translate_symbolic_boundary_name_to_id (const std::string &name) const
 
std::vector< types::boundary_idtranslate_symbolic_boundary_names_to_ids (const std::vector< std::string > &names) const
 
std::string translate_id_to_symbol_name (const types::boundary_id boundary_id) const
 
virtual std::set< std::pair< std::pair< types::boundary_id, types::boundary_id >, unsigned int > > get_periodic_boundary_pairs () const
 
virtual bool has_curved_elements () const
 
- Public Member Functions inherited from aspect::SimulatorAccess< dim >
 SimulatorAccess ()
 
 SimulatorAccess (const Simulator< dim > &simulator_object)
 
virtual ~SimulatorAccess ()
 
virtual void initialize_simulator (const Simulator< dim > &simulator_object)
 
template<typename PostprocessorType >
PostprocessorType * find_postprocessor () const
 
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
 
TimerOutputget_computing_timer () const
 
const ConditionalOStreamget_pcout () const
 
double get_time () const
 
double get_timestep () const
 
double get_old_timestep () const
 
unsigned int get_timestep_number () 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
 
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
 
void compute_material_model_input_values (const LinearAlgebra::BlockVector &input_solution, const FEValuesBase< dim, dim > &input_finite_element_values, const typename DoFHandler< dim >::active_cell_iterator &cell, const bool compute_strainrate, MaterialModel::MaterialModelInputs< dim > &material_model_inputs) 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
 
DEAL_II_DEPRECATED const BoundaryTemperature::Interface< dim > & get_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
 
DEAL_II_DEPRECATED const BoundaryComposition::Interface< dim > & get_boundary_composition () const
 
const BoundaryComposition::Manager< dim > & get_boundary_composition_manager () const
 
const std::map< types::boundary_id, std::unique_ptr< BoundaryTraction::Interface< dim > > > & get_boundary_traction () const
 
DEAL_II_DEPRECATED const InitialTemperature::Interface< dim > & get_initial_temperature () const
 
const InitialTemperature::Manager< dim > & get_initial_temperature_manager () const
 
DEAL_II_DEPRECATED const InitialComposition::Interface< dim > & get_initial_composition () 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 WorldBuilder::World & get_world_builder () const
 
const MeshDeformation::MeshDeformationHandler< dim > & get_mesh_deformation_handler () const
 
const LateralAveraging< dim > & get_lateral_averaging () const
 
const ConstraintMatrix & 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
 
TableHandlerget_statistics_object () const
 
template<typename PostprocessorType >
DEAL_II_DEPRECATED PostprocessorType * find_postprocessor () const
 
const Postprocess::Manager< dim > & get_postprocess_manager () const
 

Static Public Member Functions

static void declare_parameters (ParameterHandler &prm)
 
- Static Public Member Functions inherited from aspect::GeometryModel::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_boundary_ids (parallel::distributed::Triangulation< dim > &coarse_grid) const
 

Private Attributes

std::vector< Point< 2 > > corners
 
double semi_major_axis_a
 
double eccentricity
 
double semi_minor_axis_b
 
double rot_para_to_para_angle
 
double para_to_rect_angle
 
double rotation_longitude
 
double rotation_latitude
 
double bottom_depth
 
double westLongitude
 
double eastLongitude
 
double northLatitude
 
double southLatitude
 
unsigned int EW_subdiv
 
unsigned int NS_subdiv
 
unsigned int depth_subdiv
 
EllipsoidalChunkGeometry manifold
 

Detailed Description

template<int dim>
class aspect::GeometryModel::EllipsoidalChunk< dim >

A class that describes a geometry for an ellipsoid such as the WGS84 model of the earth.

This geometry model implements a (3d) ellipsoidal chunk geometry where two of the axis have the same length. The ellipsoidal chunk can be a non-coordinate parallel part of the ellipsoid.

Author
This plugin is a joined effort of Menno Fraters, D. Sarah Stamps and Wolfgang Bangerth

Definition at line 46 of file ellipsoidal_chunk.h.

Member Function Documentation

§ initialize()

template<int dim>
virtual void aspect::GeometryModel::EllipsoidalChunk< dim >::initialize ( )
virtual

Initialize function

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

§ create_coarse_mesh()

template<int dim>
virtual void aspect::GeometryModel::EllipsoidalChunk< dim >::create_coarse_mesh ( parallel::distributed::Triangulation< dim > &  coarse_grid) const
virtual

Generate a coarse mesh for the geometry described by this class.

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

§ length_scale()

template<int dim>
virtual double aspect::GeometryModel::EllipsoidalChunk< dim >::length_scale ( ) const
virtual

Return the typical length scale one would expect of features in this geometry, assuming realistic parameters.

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

§ depth()

template<int dim>
virtual double aspect::GeometryModel::EllipsoidalChunk< dim >::depth ( const Point< dim > &  position) const
virtual

Return the depth that corresponds to the given position. The documentation of the base class (see GeometryModel::Interface::depth()) describes in detail how "depth" is interpreted in general.

Computing a depth requires a geometry model to define a "vertical" direction. The current class considers the radial vector away from the origin as vertical and considers the "outer" boundary as the "surface". In almost all cases one will use a gravity model that also matches these definitions.

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

§ height_above_reference_surface()

template<int dim>
virtual double aspect::GeometryModel::EllipsoidalChunk< dim >::height_above_reference_surface ( const Point< dim > &  position) const
virtual

Placeholder for a function returning the height of the given position relative to the reference model surface.

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

§ representative_point()

template<int dim>
virtual Point<dim> aspect::GeometryModel::EllipsoidalChunk< dim >::representative_point ( const double  depth) const
virtual

Returns a point in the center of the domain.

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

§ point_is_in_domain()

template<int dim>
virtual bool aspect::GeometryModel::EllipsoidalChunk< dim >::point_is_in_domain ( const Point< dim > &  point) const
virtual

Return whether the given point lies within the domain specified by the geometry. This function does not take into account initial or dynamic surface topography.

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

§ maximal_depth()

template<int dim>
virtual double aspect::GeometryModel::EllipsoidalChunk< dim >::maximal_depth ( ) const
virtual

Returns the bottom depth which was used to create the geometry and which is defined by the depth parameter.

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

§ get_used_boundary_indicators()

template<int dim>
virtual std::set<types::boundary_id> aspect::GeometryModel::EllipsoidalChunk< dim >::get_used_boundary_indicators ( ) const
virtual

Return the set of boundary indicators that are used by this model. This information is used to determine what boundary indicators can be used in the input file.

The box model uses boundary indicators zero through 2*dim-1, with the first two being the faces perpendicular to the x-axis, the next two perpendicular to the y-axis, etc.

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

§ get_symbolic_boundary_names_map()

template<int dim>
virtual std::map<std::string,types::boundary_id> aspect::GeometryModel::EllipsoidalChunk< dim >::get_symbolic_boundary_names_map ( ) const
virtual

Return a mapping from symbolic names of each part of the boundary to the corresponding boundary indicator. This allows users to specify names, not just numbers in their input files when describing which parts of the boundary have to satisfy which boundary conditions.

An example would be that the "box" geometry returns a map of the form {{"left"->0}, {"right"->1}, {"bottom"->2}, {"top"->3}} in 2d.

The default implementation of this function returns an empty map. This still allows the use of a geometry model that does not implement this function but forces the user to identify parts of the boundary by their boundary indicator number, rather than using a symbolic name.

Note
Names may contain spaces and numbers, but they may not contain special characters and they should not equal the text representation of numbers (e.g., a name "10" is ill-advised).
Since in practice boundary indicators can be provided either via number or symbolic name, the mapping from something given in the input is not entirely trivial – in particular, because a function also has to do some error checking that a given string in fact matches any known boundary indicator. To this end, use GeometryModel::translate_boundary_indicator() and GeometryModel::translate_boundary_indicators().
Returns
A map from symbolic names to boundary indicators. The map should provide a symbolic name for each used boundary indicator as returned by get_user_boundary_indicators(). In the end, however, a geometry model may define multiple symbolic names for the same boundary or not define any.

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

§ natural_coordinate_system()

template<int dim>
virtual aspect::Utilities::Coordinates::CoordinateSystem aspect::GeometryModel::EllipsoidalChunk< dim >::natural_coordinate_system ( ) const
virtual

Returns what the natural coordinate system for this geometry model is.

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

§ cartesian_to_natural_coordinates()

template<int dim>
virtual std::array<double,dim> aspect::GeometryModel::EllipsoidalChunk< dim >::cartesian_to_natural_coordinates ( const Point< dim > &  position) const
virtual

Takes the Cartesian points (x,z or x,y,z) and returns standardized coordinates which are most 'natural' to the geometry model. For a ellispoidal chunk this is (radius, longitude) in 2d and (radius, longitude, latitude) in 3d. Note that internally the coordinates are stored in longitude, latitude, depth.

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

§ natural_to_cartesian_coordinates()

template<int dim>
virtual Point<dim> aspect::GeometryModel::EllipsoidalChunk< dim >::natural_to_cartesian_coordinates ( const std::array< double, dim > &  position) const
virtual

Undoes the action of cartesian_to_natural_coordinates, and turns the coordinate system which is most 'natural' to the geometry model into Cartesian coordinates.

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

§ declare_parameters()

template<int dim>
static void aspect::GeometryModel::EllipsoidalChunk< dim >::declare_parameters ( ParameterHandler prm)
static

Declare the parameters this class takes through input files.

§ parse_parameters()

template<int dim>
virtual void aspect::GeometryModel::EllipsoidalChunk< dim >::parse_parameters ( ParameterHandler prm)
virtual

Read the parameters this class declares from the parameter file.

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

§ get_radius()

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::get_radius ( const Point< dim > &  position) const

Calculate radius at current position.

§ get_semi_minor_axis_b()

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::get_semi_minor_axis_b ( ) const

Retrieve the semi minor axis b value.

§ get_semi_major_axis_a()

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::get_semi_major_axis_a ( ) const

Retrieve the semi major axis a value.

§ get_eccentricity()

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::get_eccentricity ( ) const

Retrieve the value of the eccentricity.

§ get_corners()

template<int dim>
const std::vector<Point<2> >& aspect::GeometryModel::EllipsoidalChunk< dim >::get_corners ( ) const

Retrieve the corners used to create the ellipsoid. This variable contains four vectors with each two elements. Each set of two elements represents a longitude and latitude value. The four vectors represent respectively the point in the North-East, North-West, South-West and South-East.

§ get_manifold()

template<int dim>
EllipsoidalChunkGeometry aspect::GeometryModel::EllipsoidalChunk< dim >::get_manifold ( ) const

Retrieve the manifold object.

§ set_boundary_ids()

template<int dim>
void aspect::GeometryModel::EllipsoidalChunk< dim >::set_boundary_ids ( parallel::distributed::Triangulation< dim > &  coarse_grid) const
private

Member Data Documentation

§ corners

template<int dim>
std::vector<Point<2> > aspect::GeometryModel::EllipsoidalChunk< dim >::corners
private

Definition at line 325 of file ellipsoidal_chunk.h.

§ semi_major_axis_a

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::semi_major_axis_a
private

Definition at line 326 of file ellipsoidal_chunk.h.

§ eccentricity

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::eccentricity
private

Definition at line 327 of file ellipsoidal_chunk.h.

§ semi_minor_axis_b

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::semi_minor_axis_b
private

Definition at line 328 of file ellipsoidal_chunk.h.

§ rot_para_to_para_angle

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::rot_para_to_para_angle
private

Definition at line 329 of file ellipsoidal_chunk.h.

§ para_to_rect_angle

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::para_to_rect_angle
private

Definition at line 330 of file ellipsoidal_chunk.h.

§ rotation_longitude

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::rotation_longitude
private

Definition at line 331 of file ellipsoidal_chunk.h.

§ rotation_latitude

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::rotation_latitude
private

Definition at line 332 of file ellipsoidal_chunk.h.

§ bottom_depth

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::bottom_depth
private

Definition at line 333 of file ellipsoidal_chunk.h.

§ westLongitude

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::westLongitude
private

Definition at line 334 of file ellipsoidal_chunk.h.

§ eastLongitude

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::eastLongitude
private

Definition at line 335 of file ellipsoidal_chunk.h.

§ northLatitude

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::northLatitude
private

Definition at line 336 of file ellipsoidal_chunk.h.

§ southLatitude

template<int dim>
double aspect::GeometryModel::EllipsoidalChunk< dim >::southLatitude
private

Definition at line 337 of file ellipsoidal_chunk.h.

§ EW_subdiv

template<int dim>
unsigned int aspect::GeometryModel::EllipsoidalChunk< dim >::EW_subdiv
private

Definition at line 339 of file ellipsoidal_chunk.h.

§ NS_subdiv

template<int dim>
unsigned int aspect::GeometryModel::EllipsoidalChunk< dim >::NS_subdiv
private

Definition at line 340 of file ellipsoidal_chunk.h.

§ depth_subdiv

template<int dim>
unsigned int aspect::GeometryModel::EllipsoidalChunk< dim >::depth_subdiv
private

Definition at line 341 of file ellipsoidal_chunk.h.

§ manifold

template<int dim>
EllipsoidalChunkGeometry aspect::GeometryModel::EllipsoidalChunk< dim >::manifold
private

Construct manifold object Pointer to an object that describes the geometry.

Definition at line 348 of file ellipsoidal_chunk.h.


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