ASPECT
structured_data.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2014 - 2023 by the authors of the ASPECT code.
3 
4  This file is part of ASPECT.
5 
6  ASPECT is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2, or (at your option)
9  any later version.
10 
11  ASPECT is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with ASPECT; see the file LICENSE. If not see
18  <http://www.gnu.org/licenses/>.
19 */
20 
21 
22 #ifndef _aspect_structured_data_h
23 #define _aspect_structured_data_h
24 
25 #include <aspect/global.h>
27 
28 #include <array>
29 
30 namespace aspect
31 {
32  namespace Utilities
33  {
34  using namespace dealii;
35 
55  template <int dim>
57  {
58  public:
68  StructuredDataLookup(const unsigned int n_components,
69  const double scale_factor);
70 
77  explicit StructuredDataLookup(const double scale_factor);
78 
114  void reinit(const std::vector<std::string> &column_names,
115  std::vector<std::vector<double>> &&coordinate_values,
116  std::vector<Table<dim,double>> &&data_table,
117  const MPI_Comm &mpi_communicator = MPI_COMM_SELF,
118  const unsigned int root_process = numbers::invalid_unsigned_int);
119 
125  void
126  load_file(const std::string &filename,
127  const MPI_Comm &communicator);
128 
137  void load_netcdf(const std::string &filename, const std::vector<std::string> &data_column_names = {});
138 
149  double
150  get_data(const Point<dim> &position,
151  const unsigned int component) const;
152 
162  Tensor<1,dim>
163  get_gradients(const Point<dim> &position,
164  const unsigned int component);
165 
172  std::vector<std::string>
173  get_column_names() const;
174 
180  bool
181  has_equidistant_coordinates() const;
182 
195  const std::vector<double> &
196  get_interpolation_point_coordinates(const unsigned int dimension) const;
197 
203  unsigned int
204  get_column_index_from_name(const std::string &column_name) const;
205 
211  std::string
212  get_column_name_from_index(const unsigned int column_index) const;
213 
217  double get_maximum_component_value(const unsigned int component) const;
218 
219  private:
223  unsigned int n_components;
224 
230  std::vector<std::string> data_component_names;
231 
237  std::vector<std::unique_ptr<Function<dim>>> data;
238 
242  std::array<std::vector<double>,dim> coordinate_values;
243 
247  std::vector<double> maximum_component_value;
248 
252  TableIndices<dim> table_points;
253 
258  const double scale_factor;
259 
265 
270  TableIndices<dim>
271  compute_table_indices(const TableIndices<dim> &sizes, const std::size_t idx) const;
272 
273  };
274 
279  template <int dim>
281  {
282  public:
286  AsciiDataBase();
287 
291  static
292  void
293  declare_parameters (ParameterHandler &prm,
294  const std::string &default_directory,
295  const std::string &default_filename,
296  const std::string &subsection_name = "Ascii data model");
297 
301  void
302  parse_parameters (ParameterHandler &prm,
303  const std::string &subsection_name = "Ascii data model");
304 
308  std::string data_directory;
309 
315  std::string data_file_name;
316 
323  double scale_factor;
324  };
325 
330  template <int dim>
332  {
333  public:
338 
343  virtual
344  void
345  initialize (const std::set<types::boundary_id> &boundary_ids,
346  const unsigned int components);
347 
354  void
355  update();
356 
360  double
361  get_data_component (const types::boundary_id boundary_indicator,
362  const Point<dim> &position,
363  const unsigned int component) const;
364 
368  double
369  get_maximum_component_value (const types::boundary_id boundary_indicator,
370  const unsigned int component) const;
371 
375  Tensor<1,dim-1>
376  vector_gradient(const types::boundary_id boundary_indicator,
377  const Point<dim> &p,
378  const unsigned int component) const;
379 
396  static
397  void
398  declare_parameters (ParameterHandler &prm,
399  const std::string &default_directory,
400  const std::string &default_filename,
401  const std::string &subsection_name = "Ascii data model",
402  const bool declare_time_dependent_parameters = true);
403 
417  void
418  parse_parameters (ParameterHandler &prm,
419  const std::string &subsection_name = "Ascii data model",
420  const bool parse_time_dependent_parameters = true);
421 
422  protected:
428 
434 
442 
448 
453  double time_weight;
454 
461 
466  std::map<types::boundary_id,
468 
472  std::map<types::boundary_id,
474 
478  void
479  update_data (const types::boundary_id boundary_id,
480  const bool reload_both_files);
481 
486  void
487  end_time_dependence ();
488 
492  std::string
493  create_filename (const int filenumber,
494  const types::boundary_id boundary_id) const;
495  };
496 
497 
498 
503  template <int dim>
505  {
506  public:
511 
516  virtual
517  void
518  initialize (const unsigned int components);
519 
520 
524  double
525  get_data_component (const Point<dim> &position,
526  const unsigned int component) const;
527 
531  static
532  void
533  declare_parameters (ParameterHandler &prm,
534  const std::string &default_directory,
535  const std::string &default_filename,
536  const std::string &subsection_name = "Ascii data model");
537 
541  void
542  parse_parameters (ParameterHandler &prm,
543  const std::string &subsection_name = "Ascii data model");
544 
545  protected:
550  std::unique_ptr<aspect::Utilities::StructuredDataLookup<dim>> lookup;
551 
557  std::unique_ptr<aspect::Utilities::StructuredDataLookup<3>> slice_lookup;
558 
564 
570  Tensor<2,3> rotation_matrix;
571  };
572 
573 
578  template <int dim>
580  {
581  public:
586 
591  virtual
592  void
593  initialize (const unsigned int components);
594 
595 
599  double
600  get_data_component (const Point<dim> &position,
601  const unsigned int component) const;
602 
603 
607  static
608  void
609  declare_parameters (ParameterHandler &prm,
610  const std::string &default_directory,
611  const std::string &default_filename,
612  const std::string &subsection_name = "Ascii data model");
613 
617  void
618  parse_parameters (ParameterHandler &prm,
619  const std::string &subsection_name = "Ascii data model");
620 
621  protected:
626  std::vector<std::unique_ptr<aspect::Utilities::StructuredDataLookup<dim-1>>> lookups;
627 
628  private:
629 
633  std::string data_directory;
634 
638  std::vector<std::string> data_file_names;
639 
644 
648  std::string interpolation_scheme;
649 
650 
651  };
652 
653 
657  template <int dim>
659  {
660  public:
665 
670  virtual
671  void
672  initialize (const MPI_Comm &communicator);
673 
674 
678  double
679  get_data_component (const Point<1> &position,
680  const unsigned int component) const;
681 
688  std::vector<std::string>
689  get_column_names() const;
690 
705  const std::vector<double> &
706  get_interpolation_point_coordinates() const;
707 
713  unsigned int
714  get_column_index_from_name(const std::string &column_name) const;
715 
721  unsigned int
722  maybe_get_column_index_from_name(const std::string &column_name) const;
723 
729  std::string
730  get_column_name_from_index(const unsigned int column_index) const;
731  protected:
736  std::unique_ptr<aspect::Utilities::StructuredDataLookup<1>> lookup;
737  };
738 
739 
740 
741  template<int dim>
743  }
744 }
745 
746 #endif
std::unique_ptr< aspect::Utilities::StructuredDataLookup< 3 > > slice_lookup
std::map< types::boundary_id, std::unique_ptr< aspect::Utilities::StructuredDataLookup< dim-1 > > > lookups
std::vector< std::unique_ptr< Function< dim > > > data
std::vector< std::string > data_component_names
void declare_parameters(ParameterHandler &prm)
std::vector< std::string > data_file_names
std::vector< double > maximum_component_value
std::vector< std::unique_ptr< aspect::Utilities::StructuredDataLookup< dim-1 > > > lookups
std::unique_ptr< aspect::Utilities::StructuredDataLookup< 1 > > lookup
std::map< types::boundary_id, std::unique_ptr< aspect::Utilities::StructuredDataLookup< dim-1 > > > old_lookups
std::array< std::vector< double >, dim > coordinate_values
Definition: compat.h:42
std::unique_ptr< aspect::Utilities::StructuredDataLookup< dim > > lookup