World Builder  1.1.0-pre
A geodynamic initial conditions generator
uniform.cc
Go to the documentation of this file.
1 /*
2  Copyright (C) 2018-2024 by the authors of the World Builder code.
3 
4  This file is part of the World Builder.
5 
6  This program is free software: you can redistribute it and/or modify
7  it under the terms of the GNU Lesser General Public License as published
8  by the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  This program 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 Lesser General Public License for more details.
15 
16  You should have received a copy of the GNU Lesser General Public License
17  along with this program. If not, see <https://www.gnu.org/licenses/>.
18 */
19 
21 
22 
23 #include "world_builder/nan.h"
31 
32 
33 namespace WorldBuilder
34 {
35 
36  using namespace Utilities;
37 
38  namespace Features
39  {
40  namespace OceanicPlateModels
41  {
42  namespace Grains
43  {
44  Uniform::Uniform(WorldBuilder::World *world_)
45  :
46  min_depth(NaN::DSNAN),
47  max_depth(NaN::DSNAN)
48 
49  {
50  this->world = world_;
51  this->name = "uniform";
52  }
53 
55  = default;
56 
57  void
58  Uniform::declare_entries(Parameters &prm, const std::string & /*unused*/)
59  {
60  // Document plugin and require entries if needed.
61  // Add compositions the required parameters.
62  prm.declare_entry("", Types::Object({"compositions"}),
63  "Uniform grains model. All grains start exactly the same.");
64 
65  // Declare entries of this plugin
67  "The depth in meters from which the composition of this feature is present.");
68  prm.declare_entry("max depth", Types::OneOf(Types::Double(std::numeric_limits<double>::max()),Types::Array(Types::ValueAtPoints(std::numeric_limits<double>::max(), 2.))),
69  "The depth in meters to which the composition of this feature is present.");
70 
71  prm.declare_entry("compositions", Types::Array(Types::UnsignedInt(),0),
72  "A list with the integer labels of the composition which are present there.");
73 
74  prm.declare_entry("rotation matrices", Types::Array(Types::Array(Types::Array(Types::Double(0),3,3),3,3),0),
75  "A list with the labels of the grains which are present there for each compositions.");
76 
77  prm.declare_entry("Euler angles z-x-z", Types::Array(Types::Array(Types::Double(0),3,3),0),
78  "A list with the z-x-z Euler angles of the grains which are present there for each compositions.");
79 
80  prm.declare_entry("orientation operation", Types::String("replace", std::vector<std::string> {"replace"}),
81  "Whether the value should replace any value previously defined at this location (replace) or "
82  "add the value to the previously define value (add, not implemented). Replacing implies that all values not "
83  "explicitly defined are set to zero.");
84 
85  prm.declare_entry("grain sizes",
87  "A list of the size of all of the grains in each composition. If set to <0, the size will be set so that the total is equal to 1.");
88 
89 
90  }
91 
92  void
93  Uniform::parse_entries(Parameters &prm, const std::vector<Point<2>> &coordinates)
94  {
95  min_depth_surface = Objects::Surface(prm.get("min depth",coordinates));
97  max_depth_surface = Objects::Surface(prm.get("max depth",coordinates));
99  compositions = prm.get_vector<unsigned int>("compositions");
100 
101  const bool set_euler_angles = prm.check_entry("Euler angles z-x-z");
102  const bool set_rotation_matrices = prm.check_entry("rotation matrices");
103 
104  WBAssertThrow(!(set_euler_angles == true && set_rotation_matrices == true),
105  "Only Euler angles or Rotation matrices may be set, but both are set for " << prm.get_full_json_path());
106 
107 
108  WBAssertThrow(!(set_euler_angles == false && set_rotation_matrices == false),
109  "Euler angles or Rotation matrices have to be set, but neither are set for " << prm.get_full_json_path());
110 
111  if (set_euler_angles)
112  {
113  std::vector<std::array<double,3> > euler_angles_vector = prm.get_vector<std::array<double,3> >("Euler angles z-x-z");
114  rotation_matrices.resize(euler_angles_vector.size());
115  for (size_t i = 0; i<euler_angles_vector.size(); ++i)
116  {
117  rotation_matrices[i] = WorldBuilder::Utilities::euler_angles_to_rotation_matrix(euler_angles_vector[i][0],euler_angles_vector[i][1],euler_angles_vector[i][2]);
118  }
119 
120  }
121  else
122  {
123  rotation_matrices = prm.get_vector<std::array<std::array<double,3>,3> >("rotation matrices");
124  }
125  operation = prm.get<std::string>("orientation operation");
126  grain_sizes = prm.get_vector<double>("grain sizes");
127 
128 
130  "There are not the same amount of compositions (" << compositions.size()
131  << ") and rotation_matrices (" << rotation_matrices.size() << ").");
132  WBAssertThrow(compositions.size() == grain_sizes.size(),
133  "There are not the same amount of compositions (" << compositions.size()
134  << ") and grain_sizes (" << grain_sizes.size() << ").");
135  }
136 
137 
139  Uniform::get_grains(const Point<3> & /*position_in_cartesian_coordinates*/,
140  const Objects::NaturalCoordinate &position_in_natural_coordinates,
141  const double depth,
142  const unsigned int composition_number,
143  WorldBuilder::grains grains_,
144  const double /*feature_min_depth*/,
145  const double /*feature_max_depth*/) const
146  {
147  WorldBuilder::grains grains_local = grains_;
148  if (depth <= max_depth && depth >= min_depth)
149  {
150  const double min_depth_local = min_depth_surface.constant_value ? min_depth : min_depth_surface.local_value(position_in_natural_coordinates.get_surface_point()).interpolated_value;
151  const double max_depth_local = max_depth_surface.constant_value ? max_depth : max_depth_surface.local_value(position_in_natural_coordinates.get_surface_point()).interpolated_value;
152  if (depth <= max_depth_local && depth >= min_depth_local)
153  {
154  for (unsigned int i =0; i < compositions.size(); ++i)
155  {
156  if (compositions[i] == composition_number)
157  {
158  std::fill(grains_local.rotation_matrices.begin(),grains_local.rotation_matrices.end(),rotation_matrices[i]);
159 
160  const double size = grain_sizes[i] < 0 ? 1.0/static_cast<double>(grains_local.sizes.size()) : grain_sizes[i];
161  std::fill(grains_local.sizes.begin(),grains_local.sizes.end(),size);
162 
163 
164  return grains_local;
165  }
166  }
167  }
168  }
169  return grains_local;
170  }
172  } // namespace Grains
173  } // namespace OceanicPlateModels
174  } // namespace Features
175 } // namespace WorldBuilder
176 
177 
bool check_entry(const std::string &name) const
Definition: parameters.cc:205
#define WB_REGISTER_FEATURE_OCEANIC_PLATE_GRAINS_MODEL(classname, name)
Definition: interface.h:148
const double DSNAN
Definition: nan.h:41
void parse_entries(Parameters &prm, const std::vector< Point< 2 >> &coordinates) override final
Definition: uniform.cc:93
std::vector< double > sizes
Definition: grains.h:47
std::vector< std::array< std::array< double, 3 >, 3 > > rotation_matrices
Definition: uniform.h:106
WorldBuilder::grains get_grains(const Point< 3 > &position, const Objects::NaturalCoordinate &position_in_natural_coordinates, const double depth, const unsigned int composition_number, WorldBuilder::grains grains, const double feature_min_depth, const double feature_max_depth) const override final
Definition: uniform.cc:139
std::string get_full_json_path(size_t max_size=std::numeric_limits< size_t >::max()) const
Definition: parameters.cc:1933
#define WBAssertThrow(condition, message)
Definition: assert.h:40
SurfaceValueInfo local_value(const Point< 2 > &check_point) const
Definition: surface.cc:149
std::array< std::array< double, 3 >, 3 > euler_angles_to_rotation_matrix(double phi1_d, double theta_d, double phi2_d)
Definition: utilities.cc:1147
std::vector< std::array< std::array< double, 3 >, 3 > > rotation_matrices
Definition: grains.h:51
void declare_entry(const std::string &name, const Types::Interface &type, const std::string &documentation)
Definition: parameters.cc:197
std::vector< T > get_vector(const std::string &name)
T get(const std::string &name)
static void declare_entries(Parameters &prm, const std::string &parent_name="")
Definition: uniform.cc:58