along interface velocity for subducting slab

[lhy11009]

Add an along surface velocity option for the subducting slab. This will trace the surface of the subducting slab and rotate the velocity to align with the surface. The free parameter is the magnitude of the velocity.

Here is the layout of the test, which connects a horizontal velocity in the plate to the along-surface velocity in the slab.

This is tentative. Eventually, we might want to tweak the velocity in each of the slab segments, but this would be an easy start.

[coveralls]

Pull Request Test Coverage Report for Build 23763068587

Warning: This coverage report may be inaccurate.

This pull request's base commit is no longer the HEAD commit of its target branch. This means it includes changes from outside the original pull request, including, potentially, unrelated coverage changes.

• For more information on this, see Tracking coverage changes with pull request builds.
• To avoid this issue with future PRs, see these Recommended CI Configurations.
• For a quick fix, rebase this PR at GitHub. Your next report should be accurate.

Details

• 15 of 42 (35.71%) changed or added relevant lines in 2 files are covered.
• 59 unchanged lines in 7 files lost coverage.
• Overall coverage decreased (-0.4%) to 98.079%

---

Changes Missing Coverage	Covered Lines	Changed/Added Lines	%
source/world_builder/features/subducting_plate_models/velocity/along_surface.cc	7	34	20.59%

Files with Coverage Reduction	New Missed Lines	%
source/gwb-dat/main.cc	1	99.38%
source/world_builder/features/fault.cc	1	99.74%
source/world_builder/features/subducting_plate.cc	1	99.75%
source/world_builder/features/plume.cc	3	98.45%
source/gwb-grid/main.cc	7	97.58%
source/world_builder/world.cc	20	92.51%
source/world_builder/parameters.cc	26	97.13%

Totals
Change from base Build 23666386453:	-0.4%
Covered Lines:	9904
Relevant Lines:	10098

---

💛 - Coveralls

[github-actions]

Benchmark	Main	Feature	Difference (99.9% CI)
Slab interpolation simple none	1.316 ± 0.010 (s=345)	1.316 ± 0.011 (s=341)	-0.2% .. +0.2%
Slab interpolation curved simple none	1.310 ± 0.007 (s=343)	1.311 ± 0.009 (s=346)	-0.0% .. +0.3%
Spherical slab interpolation simple none	1.279 ± 0.007 (s=361)	1.283 ± 0.010 (s=344)	+0.2% .. +0.5%
Slab interpolation simple curved CMS	1.369 ± 0.010 (s=331)	1.367 ± 0.005 (s=329)	-0.3% .. +0.0%
Spherical slab interpolation simple CMS	1.634 ± 0.013 (s=269)	1.639 ± 0.014 (s=283)	+0.1% .. +0.6%
Spherical fault interpolation simple none	1.287 ± 0.010 (s=373)	1.288 ± 0.010 (s=329)	-0.1% .. +0.3%
Cartesian min max surface	3.132 ± 0.014 (s=144)	3.146 ± 0.017 (s=145)	+0.3% .. +0.7%
Spherical min max surface	8.635 ± 0.091 (s=51)	8.634 ± 0.101 (s=56)	-0.7% .. +0.7%

[MFraters]

Cool! Yes, the next step would be have a velocity per segment something like the angles, either [1.5] or [1.5,2.2].

[MFraters]

Looking at it a bit more, I am not sure how to get the result you got, I get the following result.

with the following wb file:

{
  "version":"1.2",
  "coordinate system":{"model":"cartesian"},
  "cross section":[[0,450e3],[10e3,450e3]],
  "features":
  [
   { 
     "model":"subducting plate", "name":"Slab", "dip point":[0,0],
     "coordinates":[[1500e3,1000e3],[1600e3,350e3],[1500e3,0]],
     "segments":
       [
         {"length":300e3, "thickness":[300e3], "top truncation":[-100e3], "angle":[0,60],
          "composition models":[
            {"model":"uniform", "compositions":[3], "max distance slab top":50e3},
            {"model":"uniform", "compositions":[2], "min distance slab top":50e3, "max distance slab top":100e3}]},
         {"length":500e3, "thickness":[300e3], "top truncation":[-100e3], "angle":[60,20]}
       ],
     "sections":
       [
         {"coordinate":0,
          "segments":[{"length":300e3, "thickness":[300e3], "top truncation":[-100e3], "angle":[0,60]},
                      {"length":400e3, "thickness":[300e3], "top truncation":[-100e3], "angle":[60]}],
                       "composition models":[{"model":"uniform", "compositions":[1],"max distance slab top":100e3}]}
       ],
     "temperature models":[{"model":"mass conserving", "density":3300, "spreading velocity":0.02, "subducting velocity":0.02,
                            "ridge coordinates":[[[5e6,-1],[5e6,2000e3]]],"coupling depth":50e3,
                            "min distance slab top":-200e3, "max distance slab top":300e3}],
     "composition models":[{"model":"uniform", "compositions":[2], "max distance slab top":100e3}],
      "velocity models": [{"model":"along surface", "velocity": [1,2,3]}]
   }
  ]
}

and world builder file:

# output variables
grid_type = cartesian
dim = 2
compositions = 6
vtu_output_format = ASCII

# domain of the grid
x_min = -1000e3
x_max = 2000e3 
z_min = 0 
z_max = 600e3 

# low res grid properties
n_cell_x = 600 
n_cell_z = 100

# shown grid properties
# n_cell_x = 2400 
# n_cell_z = 400

I also thought the user would provide a scalar instead a vector, and maybe later on we could add an angle offset.

[lhy11009]

If we add a sign modification to Vx depending on the relative position of the dip point and the trench point, and if dip point have smaller x than the trench point, we change the sign of Vx.
With Menno's help, I also adjusted the angle I use to be the local angle; now the velocity perfectly aligns with the surface.
This angle is added to the distance_from_plane object

[lhy11009]

@MFraters. Now the question is how to figure out when sign_vx should be 1.0 versus -1.0. -1.0 is required to fix the Vx when the query point is left to the trench point with a smaller X coordinate. I wonder this could be figured out from the distance_from_plane object, or we have to pass the reference point to this function.

[MFraters]

Yes, I think reference_on_side_of_line is the variable you are looking for to add to the structure.