ASPECT
Public Member Functions | Private Member Functions | Private Attributes | List of all members
aspect::Utilities::StructuredDataLookup< dim > Class Template Reference

Public Member Functions

 StructuredDataLookup (const unsigned int n_components, const double scale_factor)
 
 StructuredDataLookup (const double scale_factor)
 
void reinit (const std::vector< std::string > &column_names, std::vector< std::vector< double >> &&coordinate_values, std::vector< Table< dim, double >> &&data_table, const MPI_Comm mpi_communicator=MPI_COMM_SELF, const unsigned int root_process=numbers::invalid_unsigned_int)
 
void load_ascii (const std::string &filename, const MPI_Comm communicator)
 
void load_netcdf (const std::string &filename, const std::vector< std::string > &data_column_names={})
 
void load_file (const std::string &filename, const MPI_Comm communicator)
 
double get_data (const Point< dim > &position, const unsigned int component) const
 
Tensor< 1, dim > get_gradients (const Point< dim > &position, const unsigned int component)
 
std::vector< std::string > get_column_names () const
 
bool has_equidistant_coordinates () const
 
const std::vector< double > & get_interpolation_point_coordinates (const unsigned int dimension) const
 
unsigned int get_column_index_from_name (const std::string &column_name) const
 
std::string get_column_name_from_index (const unsigned int column_index) const
 
double get_maximum_component_value (const unsigned int component) const
 

Private Member Functions

TableIndices< dim > compute_table_indices (const TableIndices< dim > &sizes, const std::size_t idx) const
 

Private Attributes

unsigned int n_components
 
std::vector< std::string > data_component_names
 
std::vector< std::unique_ptr< Function< dim > > > data
 
std::array< std::vector< double >, dim > coordinate_values
 
std::vector< double > maximum_component_value
 
TableIndices< dim > table_points
 
const double scale_factor
 
bool coordinate_values_are_equidistant
 

Detailed Description

template<int dim>
class aspect::Utilities::StructuredDataLookup< dim >

StructuredDataLookup (formerly AsciiDataLookup) represents structured data that can be read from files including in ascii format.

For ascii files the files need to be formatted as follows: Note the required format of the input data: The first lines may contain any number of comments if they begin with '#', but one of these lines needs to contain the number of grid points in each dimension as for example '# POINTS: 3 3'. The comments can optionally be followed by a single line, which does not start with '#', containing the names of the data columns. The order of the following data columns has to be 'coordinates data' with dim coordinate columns and components data columns. Note that the data in the input files need to be sorted in a specific order: the first coordinate needs to ascend first, followed by the second and so on in order to assign the correct data to the prescribed coordinates. The coordinates do not need to be equidistant.

Definition at line 56 of file structured_data.h.

Constructor & Destructor Documentation

§ StructuredDataLookup() [1/2]

template<int dim>
aspect::Utilities::StructuredDataLookup< dim >::StructuredDataLookup ( const unsigned int  n_components,
const double  scale_factor 
)

Constructor that explicitly prescribes the number of data columns in the data file. If a list of data components is provided in the data file it is checked that the length of this list is consistent with this number of components. This constructor is mostly provided for backwards compatibility. Not prescribing the number of components and instead reading them from the input file allows for more flexible files.

§ StructuredDataLookup() [2/2]

template<int dim>
aspect::Utilities::StructuredDataLookup< dim >::StructuredDataLookup ( const double  scale_factor)
explicit

This constructor relies on the list of column names at the beginning of the model file to determine the number of data components, therefore when using this constructor it is necessary to provide this list in the first uncommented line of the data file.

Member Function Documentation

§ reinit()

template<int dim>
void aspect::Utilities::StructuredDataLookup< dim >::reinit ( const std::vector< std::string > &  column_names,
std::vector< std::vector< double >> &&  coordinate_values,
std::vector< Table< dim, double >> &&  data_table,
const MPI_Comm  mpi_communicator = MPI_COMM_SELF,
const unsigned int  root_process = numbers::invalid_unsigned_int 
)

Replace the data stored in this class by the data given to this function.

Note
This is a manual way to fill the data. Consider load_file() if your data is stored in a txt/csv file.

The data consists of n_components (implicitly given by the size of column_names) specified at coordinate_values[d] points in each of the dim coordinate directions d.

The data in data_table consists of a Table for each of the n_components components.

The coordinate_values and data_table arguments are rvalue references, and the function will move the data so provided into another storage location. In other words, after the call, the variables passed as the second and third arguments may be empty or otherwise altered.

This class is able to share data between processes located on the same machine if the last argument root_process is set to anything other than numbers::invalid_unsigned_int. Then only the indicated root process within the given MPI communicator needs to pass in valid arguments whereas on all other processes the values of column_names, coordinate_values, and data_table are ignored. Instead, the indicated root process will make sure that the other processes obtain valid data, for example by sharing data tables in shared memory spaces. This reduces both the computational effort in reading data from disk and parsing it on all processes, and can vastly reduce the amount of memory required on machines where many processor cores have access to shared memory resources.

If root_process equals numbers::invalid_unsigned_int, then every process needs to pass data for all arguments and the mpi_communicator argument is ignored.

§ load_ascii()

template<int dim>
void aspect::Utilities::StructuredDataLookup< dim >::load_ascii ( const std::string &  filename,
const MPI_Comm  communicator 
)

Loads a data text file replacing the current data.

Throws an exception if the file does not exist, if the data file format is incorrect, or if the file data like the number and names of columns change (if this class contains existing data for example when used with a collection of ascii files changing over the simulation time).

This function supports normal text / ASCII files that can optionally be compressed using .gz. If the given filename is a URL and libDAB support is enabled, the data is downloaded and parsed from the ASCII data received.

This function uses the given MPI communicator to only load the data on rank 0 and broadcasting the information to all other ranks.

§ load_netcdf()

template<int dim>
void aspect::Utilities::StructuredDataLookup< dim >::load_netcdf ( const std::string &  filename,
const std::vector< std::string > &  data_column_names = {} 
)

Fill the current object with data read from a NetCDF file with filename filename. This call will fail if ASPECT is not configured with NetCDF support.

data_column_names specifies the list of data columns to load (in the specified order). If an empty vector is passed, all columns will be loaded.

§ load_file()

template<int dim>
void aspect::Utilities::StructuredDataLookup< dim >::load_file ( const std::string &  filename,
const MPI_Comm  communicator 
)

Loads data from a file replacing the current data.

Throws an exception if the file does not exist, if the data file format is incorrect, or if the file data like the number and names of columns change (if this class contains existing data for example when used with a collection of ascii files changing over the simulation time).

The following formats are currently supported:

  • ASCII files (typically ending in .txt)
  • gzip compressed ASCII files (ending in .gz)
  • URLs starting with "http" (handled by libDAB)

§ get_data()

template<int dim>
double aspect::Utilities::StructuredDataLookup< dim >::get_data ( const Point< dim > &  position,
const unsigned int  component 
) const

Returns the computed data (velocity, temperature, etc. - according to the used plugin) in Cartesian coordinates.

Parameters
positionThe current position to compute the data (velocity, temperature, etc.)
componentThe index (starting at 0) of the data column to be returned. The index is therefore less than the number of data columns in the data file (or specified in the constructor).

§ get_gradients()

template<int dim>
Tensor<1,dim> aspect::Utilities::StructuredDataLookup< dim >::get_gradients ( const Point< dim > &  position,
const unsigned int  component 
)

Returns the gradient of the function based on the bilinear interpolation of the data (velocity, temperature, etc. - according to the used plugin) in Cartesian coordinates.

Parameters
positionThe current position to compute the data (velocity, temperature, etc.)
componentThe index of the data column to be returned.

§ get_column_names()

template<int dim>
std::vector<std::string> aspect::Utilities::StructuredDataLookup< dim >::get_column_names ( ) const

Returns a vector that contains the names of all data columns in the order of their appearance in the data file (and their order in the memory data table). Returns an empty vector if no names are provided or the file is not read in yet.

§ has_equidistant_coordinates()

template<int dim>
bool aspect::Utilities::StructuredDataLookup< dim >::has_equidistant_coordinates ( ) const

Returns whether the stored coordinates are equidistant. If coordinates are equidistant the lookup is more efficient. Returns false if no coordinates are loaded at the moment.

§ get_interpolation_point_coordinates()

template<int dim>
const std::vector<double>& aspect::Utilities::StructuredDataLookup< dim >::get_interpolation_point_coordinates ( const unsigned int  dimension) const

Returns the coordinates of the interpolation points at which data is stored. This function can be used to determine the number of data points in each of the coordinate directions, or to query data only at exactly the positions at which they are available (avoiding interpolation).

Parameters
dimensionThe spatial direction for which to return the data coordinates, e.g. 0 for \(x\)-direction, 1 for \(y\)-direction, or equivalent values if your data coordinates are other dimensions such as temperature, pressure.

§ get_column_index_from_name()

template<int dim>
unsigned int aspect::Utilities::StructuredDataLookup< dim >::get_column_index_from_name ( const std::string &  column_name) const

Returns the column index of a column with the given name column_name. Throws an exception if no such column exists or no names were provided in the file.

§ get_column_name_from_index()

template<int dim>
std::string aspect::Utilities::StructuredDataLookup< dim >::get_column_name_from_index ( const unsigned int  column_index) const

Returns a string that contains the name of the column with index column_index. Throws an exception if no such column exists or no name was provided in the file.

§ get_maximum_component_value()

template<int dim>
double aspect::Utilities::StructuredDataLookup< dim >::get_maximum_component_value ( const unsigned int  component) const

Return the maximum value of the component values.

§ compute_table_indices()

template<int dim>
TableIndices<dim> aspect::Utilities::StructuredDataLookup< dim >::compute_table_indices ( const TableIndices< dim > &  sizes,
const std::size_t  idx 
) const
private

Computes the table indices given the size sizes of the entry with index idx.

Member Data Documentation

§ n_components

template<int dim>
unsigned int aspect::Utilities::StructuredDataLookup< dim >::n_components
private

The number of data components read in (=columns in the data file).

Definition at line 255 of file structured_data.h.

§ data_component_names

template<int dim>
std::vector<std::string> aspect::Utilities::StructuredDataLookup< dim >::data_component_names
private

The names of the data components in the columns of the read file. Does not contain any strings if none are provided in the first uncommented line of the file.

Definition at line 262 of file structured_data.h.

§ data

template<int dim>
std::vector<std::unique_ptr<Function<dim> > > aspect::Utilities::StructuredDataLookup< dim >::data
private

Interpolation functions to access the data. Either InterpolatedUniformGridData or InterpolatedTensorProductGridData; the type is determined from the grid specified in the data file.

Definition at line 269 of file structured_data.h.

§ coordinate_values

template<int dim>
std::array<std::vector<double>,dim> aspect::Utilities::StructuredDataLookup< dim >::coordinate_values
private

The coordinate values in each direction as specified in the data file.

Definition at line 274 of file structured_data.h.

§ maximum_component_value

template<int dim>
std::vector<double> aspect::Utilities::StructuredDataLookup< dim >::maximum_component_value
private

The maximum value of each component

Definition at line 279 of file structured_data.h.

§ table_points

template<int dim>
TableIndices<dim> aspect::Utilities::StructuredDataLookup< dim >::table_points
private

Number of points in the data grid as specified in the data file.

Definition at line 284 of file structured_data.h.

§ scale_factor

template<int dim>
const double aspect::Utilities::StructuredDataLookup< dim >::scale_factor
private

Scales the data boundary condition by a scalar factor. Can be used to transform the unit of the data.

Definition at line 290 of file structured_data.h.

§ coordinate_values_are_equidistant

template<int dim>
bool aspect::Utilities::StructuredDataLookup< dim >::coordinate_values_are_equidistant
private

Stores whether the coordinate values are equidistant or not, this determines the type of data function stored.

Definition at line 296 of file structured_data.h.


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