Add tests for variable spreading rate

GeodynamicWorldBuilder · Feb 15, 2024 · 842deb5 · 842deb5
1 parent 3e4b573
commit 842deb5
Show file tree

Hide file tree

Showing 7 changed files with 53 additions and 27 deletions.
diff --git a/include/world_builder/features/oceanic_plate_models/temperature/half_space_model.h b/include/world_builder/features/oceanic_plate_models/temperature/half_space_model.h
@@ -91,7 +91,7 @@ namespace WorldBuilder
             double bottom_temperature;
             std::pair<std::vector<double>,std::vector<double>> spreading_velocities;
             std::vector<std::vector<Point<2> > > mid_oceanic_ridges;
-            std::vector<std::vector<double>> reworked_velocities;
+            std::vector<std::vector<double>> spreading_velocities_at_each_ridge_point;
             Operations operation;
 
         };

diff --git a/include/world_builder/parameters.h b/include/world_builder/parameters.h
@@ -126,7 +126,7 @@ namespace WorldBuilder
        * A specialized version of get which can return a values at times type.
        * \param name The name of the entry to retrieved
        */
-      std::pair<std::vector<double>,std::vector<double>> get_ya(const std::string &name);
+      std::pair<std::vector<double>,std::vector<double>> get_value_at_array(const std::string &name);
 
       /**
        * A specialized verions of get which can return vectors/arrays.

diff --git a/source/world_builder/features/oceanic_plate_models/temperature/half_space_model.cc b/source/world_builder/features/oceanic_plate_models/temperature/half_space_model.cc
@@ -105,44 +105,29 @@ namespace WorldBuilder
           operation = string_operations_to_enum(prm.get<std::string>("operation"));
           top_temperature = prm.get<double>("top temperature");
           bottom_temperature = prm.get<double>("bottom temperature");
-          spreading_velocities = prm.get_ya("spreading velocity");
+          spreading_velocities = prm.get_value_at_array("spreading velocity");
 
           mid_oceanic_ridges = prm.get_vector<std::vector<Point<2>>>("ridge coordinates");
-          // std::vector<std::vector<double>> reworked_velocities; 
           const double dtr = prm.coordinate_system->natural_coordinate_system() == spherical ? Consts::PI / 180.0 : 1.0;
 
           unsigned int index_x = 0;
           unsigned int index_y = 0;
           unsigned int test_ind = 0;
           for (index_x = 0; index_x < mid_oceanic_ridges.size(); index_x++)
             {
-              std::vector<double> relevant_spreading;
+              std::vector<double> spreading_rates_for_ridge;
               for (index_y = 0; index_y < mid_oceanic_ridges[index_x].size(); index_y++)
                 {
-                  relevant_spreading.push_back(spreading_velocities.second[test_ind]);
+                  if (spreading_velocities.second.size() <= 1)
+                    spreading_rates_for_ridge.push_back(spreading_velocities.first[0]);
+                  else
+                    spreading_rates_for_ridge.push_back(spreading_velocities.second[test_ind]);
                   test_ind += 1;
                 }
-              reworked_velocities.push_back(relevant_spreading);
-              std::cout << "Heyyy";
+              spreading_velocities_at_each_ridge_point.push_back(spreading_rates_for_ridge);
             }
-
-
-          // for (auto &ridge_coordinates : mid_oceanic_ridges)
-          //   {
-          //     std::vector<double> relevant_spreading;
-          //     index_y = 0;
-          //     for (auto &ridge_coordinate : ridge_coordinates)
-          //       {
-          //         ridge_coordinate *= dtr;
-          //         relevant_spreading.push_back(spreading_velocities.second[index_x + index_y]);
-          //         index_y += 1;
-          //       }
-          //     index_x += 1;
-          //     reworked_velocities.push_back(relevant_spreading);
-          //   }
         }
 
-
         double
         HalfSpaceModel::get_temperature(const Point<3> &position,
                                         const Objects::NaturalCoordinate &position_in_natural_coordinates,

diff --git a/source/world_builder/parameters.cc b/source/world_builder/parameters.cc
@@ -652,7 +652,7 @@ namespace WorldBuilder
   }
 
 
-  std::pair<std::vector<double>,std::vector<double>> Parameters::get_ya(const std::string &name)
+  std::pair<std::vector<double>,std::vector<double>> Parameters::get_value_at_array(const std::string &name)
   {
     // There are four cases:
     // 1. No value provided: use the default value everywhere. Return first with one value and second with size 0.
@@ -717,7 +717,7 @@ namespace WorldBuilder
                     WBAssertThrow(false, "Could not convert values of " << strict_base << "/" << name << "/0/0 into a double. "
                                   << "The provided value was \"" <<  Pointer((strict_base + "/" + name + "/0/0").c_str()).Get(parameters)->GetString() << "\".");
                   }
-                // result.first.emplace_back(value);
+                result.first.emplace_back(value);
               }
             else
               {
@@ -773,7 +773,7 @@ namespace WorldBuilder
                                                   << "\".");
                                   }
 
-                            result.first.emplace_back(value);
+                                result.first.emplace_back(value);
                                 result.second.emplace_back(testing_value);
                               }
                           }

diff --git a/tests/gwb-dat/cartesian_variable_spreading.dat b/tests/gwb-dat/cartesian_variable_spreading.dat
@@ -0,0 +1,14 @@
+# This is a comment in the data
+# file. 
+# Now define parameters:
+# dim = 3
+# compositions = 0
+# x y z d T 
+250e3 -90e3 0 50e3
+250e3 -60e3 0 50e3
+250e3 -30e3 0 50e3
+250e3 -1e3 0 50e3
+250e3 1e3 0 50e3
+250e3 30e3 0 50e3
+250e3 60e3 0 50e3
+250e3 90e3 0 50e3
diff --git a/tests/gwb-dat/cartesian_variable_spreading.wb b/tests/gwb-dat/cartesian_variable_spreading.wb
@@ -0,0 +1,17 @@
+{
+  "version": "0.6",
+  "coordinate system":{"model":"cartesian"},
+  "gravity model":{"model":"uniform", "magnitude":10},
+  "surface temperature":273, "force surface temperature":true,
+  "potential mantle temperature":1673, "thermal expansion coefficient":3.1e-5, 
+  "specific heat":1250, "thermal diffusivity":1.0e-6,
+  "features":
+  [
+    { "model":"oceanic plate", "name":"Subducting", "max depth":250e3,"min depth":0,
+      "coordinates" :[[500e3, -100e3],[500e3, 100e3],[0, 100e3],[0, -100e3]],
+      "temperature models":[
+        {"model":"half space model", "min depth":0,  "max depth":100e3, "spreading velocity":[ [0,[[0.001, 0.0005, 0.0001]]], [1,[[0.000005, 0.0005]]] ],
+        "top temperature":273,
+         "ridge coordinates": [[[100e3,-110e3], [150e3, -60e3], [100e3, 0]], [[400e3, 0], [400e3, 110e3]]]}]}      
+  ]
+}
diff --git a/tests/gwb-dat/cartesian_variable_spreading/screen-output.log b/tests/gwb-dat/cartesian_variable_spreading/screen-output.log
@@ -0,0 +1,10 @@
+# x y z d g T 
+250e3 -90e3 0 50e3 751.596 
+250e3 -60e3 0 50e3 761.356 
+250e3 -30e3 0 50e3 751.596 
+250e3 -1e3 0 50e3 728.246 
+250e3 1e3 0 50e3 329.759 
+250e3 30e3 0 50e3 489.645 
+250e3 60e3 0 50e3 574.848 
+250e3 90e3 0 50e3 639.408 
+