adjust function Utilities::calculate_effective_trench_and_plate_ages

[lhy11009]

@danieldouglas92 , @MFraters Address what left in PR #654

mainly adust the function of Utilities::calculate_effective_trench_and_plate_ages

also added a few tests

[github-actions]

Benchmark	Main	Feature	Difference (99.9% CI)
Slab interpolation simple none	1.023 ± 0.003 (s=445)	1.023 ± 0.003 (s=437)	-0.1% .. +0.0%
Slab interpolation curved simple none	1.024 ± 0.004 (s=445)	1.024 ± 0.005 (s=436)	-0.1% .. +0.1%
Spherical slab interpolation simple none	1.170 ± 0.006 (s=393)	1.167 ± 0.007 (s=380)	-0.4% .. -0.1%
Slab interpolation simple curved CMS	1.062 ± 0.004 (s=423)	1.061 ± 0.004 (s=427)	-0.1% .. +0.0%
Spherical slab interpolation simple CMS	1.557 ± 0.010 (s=273)	1.552 ± 0.011 (s=308)	-0.5% .. -0.2%
Spherical fault interpolation simple none	1.176 ± 0.005 (s=392)	1.172 ± 0.006 (s=377)	-0.5% .. -0.2%
Cartesian min max surface	2.304 ± 0.019 (s=187)	2.305 ± 0.020 (s=206)	-0.2% .. +0.3%
Spherical min max surface	7.617 ± 0.060 (s=65)	7.620 ± 0.068 (s=56)	-0.5% .. +0.6%

[danieldouglas92]

@lhy11009 This looks good! So this is from the perspective of the trench, and if subducting_velocity > spreading_velocity the trench advances and age_at_trench increases, and if subducting_velocity < spreading_velocity, the trench retreats and age_at_trench decreases. I think this makes sense, let me know if I can help with anything but I like what you've done here!

[MFraters]

Thanks for working on this @lhy11009!

I liked that you used unit tests to test the functions. In this case, I would also like to see some test where is is used in an input file, so we can make sure the whole system works as intended. Can you also add some tests for this (preferably .dat tests)?

Your new assert seems to be triggered in some of the existing test, could you have a look at that as well?

[lhy11009]

have this issue I want to assert the subducting velocity as positive in the code

  WBAssertThrow(subducting_velocity >= 0, "The subducting velocity is less than 0. "
                "Subducting velocity: " << subducting_velocity);

And then some errors are triggered. Then I realized the default value is -1.0 indicating we want to use the spreading velocity instead.

Is there a workaround where we can still check if a incorrect negative value is entered versus if we want to use the spreading velocity instead?

Some trials: I tried to initiate the subducting velocity as nan value, but there is a check for the input values. So this would trigger another error.

      prm.declare_entry("subducting velocity", Types::OneOf(Types::Double(NaN::DQNAN), Types::Array(Types::Array(Types::Double(NaN::DQNAN), 1), 1)),
                        "The velocity with which the ridge is moving through time, and how long the ridge "
                        "has been moving. First value is the velocity, second is the time. Default is [0 cm/yr, 0 yr]");

@MFraters ? Any idea.

[MFraters]

Is there a workaround where we can still check if a incorrect negative value is entered versus if we want to use the spreading velocity instead?

Thinking about this, the ideal case would be that you use an WBAssert, not WBAssertThrow here. The real check (WBAssertThrow) should be in parse parameters. In there you want to convert the value to the actual value over there. Do you think it is feasible to do something like that here?

[lhy11009]

@MFraters, thanks for your anwser but not sure if I get the idea. Are you saying inserting an assert into the parse parameters to check if subducting velocity is positive instead of checking for it in the function?
Other than this, if we want to use the new parameter of "subducting velocity" as required in the file optional. The first way would require modification of previous prm files.
Maybe at this point, a short discussion would help?
I can meet Friday 11:00 am - 3 pm. When would be a good timing?

[danieldouglas92]

What if you change the default value to 0? Require that if subducting velocity = 0 to use the spreading velocity, and throw an error if it’s negative?

[MFraters]

I can meet Friday 11:00 am - 3 pm. When would be a good timing?
Sure, 11:00 am - 2 pm PST work for me. Daniel, can you make it then?

[danieldouglas92]

@lhy11009 @MFraters Yes I can make it!

[lhy11009]

I have no problme running the gwb-dat command with

/home/lochy/Softwares/WorldBuilder/build_effective_trench_and_plate_ages/bin/gwb-dat /home/lochy/Softwares/WorldBuilder/tests/gwb-dat/slab_temperature_trench_retreat.wb /home/lochy/Softwares/WorldBuilder/tests/gwb-dat/slab_temperature_trench_retreat.dat

but running gwd-grid with command would produce an error

/home/lochy/Softwares/WorldBuilder/build_effective_trench_and_plate_ages/bin/gwb-grid /home/lochy/Softwares/WorldBuilder/tests/gwb-dat/slab_temperature_trench_retreat.wb /home/lochy/Softwares/WorldBuilder/tests/gwb-dat/slab_temperature_trench_retreat.dat

[lochy-Inspiron-16-Plus-7620:24091] *** Process received signal *** [lochy-Inspiron-16-Plus-7620:24091] Signal: Floating point exception (8) [lochy-Inspiron-16-Plus-7620:24091] Signal code: Invalid floating point operation (7) [lochy-Inspiron-16-Plus-7620:24091] Failing at address: 0x55bf02c6d9af [lochy-Inspiron-16-Plus-7620:24091] [ 0] /lib/x86_64-linux-gnu/libc.so.6(+0x42990)[0x7f25c0e42990] [lochy-Inspiron-16-Plus-7620:24091] [ 1] /home/lochy/Softwares/WorldBuilder/build_effective_trench_and_plate_ages/bin/gwb-grid(_ZSt5isnand+0x12)[0x55bf02c6d9af] [lochy-Inspiron-16-Plus-7620:24091] [ 2] /home/lochy/Softwares/WorldBuilder/build_effective_trench_and_plate_ages/bin/gwb-grid(main+0x1bdc)[0x55bf02c5f8f8] [lochy-Inspiron-16-Plus-7620:24091] [ 3] /lib/x86_64-linux-gnu/libc.so.6(+0x28150)[0x7f25c0e28150] [lochy-Inspiron-16-Plus-7620:24091] [ 4] /lib/x86_64-linux-gnu/libc.so.6(__libc_start_main+0x89)[0x7f25c0e28209] [lochy-Inspiron-16-Plus-7620:24091] [ 5] /home/lochy/Softwares/WorldBuilder/build_effective_trench_and_plate_ages/bin/gwb-grid(_start+0x25)[0x55bf02c5cb95] [lochy-Inspiron-16-Plus-7620:24091] *** End of error message *** Floating point exception (core dumped)

[lhy11009]

@MFraters do you have an idea why this happens?

[MFraters]

I assume that /home/lochy/Softwares/WorldBuilder/tests/gwb-dat/slab_temperature_trench_retreat.dat in the last command was /home/lochy/Softwares/WorldBuilder/tests/gwb-dat/slab_temperature_trench_retreat.grid?

can you post the gdb backtrace and the valgrind output?

[lhy11009]

Oh, I happened to use the dat file instead of the grid file, my bad.

[lhy11009]

The test with a retreat trench.

Temperature profiles
150ma slab with a spreading rate and subducting rate of 5 cm / yr.

subducting rate of 10 cm / yr (retreat with 5 cm/yr)

difference (second - first)

The differences lie in the mantle temperature is hotter and the slab temperature is colder.
The reason is that the slab is heating the mantle above and a higher subducting rate results in less time for heat conduction and lead to the observed structure.

[lhy11009]

Results with an advancing trench:

temperature:

difference to the previous profile (advance - no advance).

The temeprature profile with trench advance motion is hotter overall. I attribute this to the slab being younger back in time.

[lhy11009]

Another note here, I mean to change the variable "plate velocity" to "spreading_velocity" or "ridge_spreading_velocity". I saw Daniel already call it this way in pieces of comments and inside other functions. I'd like to know what do you think?

[danieldouglas92]

@lhy11009 I think it makes sense to distinguish the velocities more clearly by using "spreading velocity", this would also break all other previous input files, but since you're making "subducting velocity" a required parameter I guess it better to break multiple things at once! Thoughts on this @MFraters?

[MFraters]

Yes, this would break input files even more, but in this case I think it makes sense. It is best to get the names sorted out now properly, before the 1.0 release. Feel free to make this change.

[MFraters]

@lhy11009 Thanks for posting the pictures, they look cool!

The bottom seems to be quite sharp, is that because it is cut off? I am not sure I follow what is happening with the slab which is both warmer and cooler, but I will have to think a bit more about that :)

[lhy11009]

@MFraters , that's a good point and I mean to write more about the temperature differences.

In the case of trench retreat, an older age in the past addressed the colder temperature in the slab's internal. For the mantle above the slab, the warmer temperature is from the heat conduction. Specifically, the subducting velocity increases from 5 cm/yr in the reference case to 10 cm/yr here. So the effect is, that the transit time of the slab is decreased by half everywhere. Thus the slab has only half the time to cool the mantle above resulting in a warmer temperature in the overlying mantle.

I actually have a question myself regarding the advance case, I think it should be opposite to the retreat case.
But, as mentioned before, there are two variables: the slab age and the transit time. So it's not very straightforward to explain what happens to the overlying mantle. In the case of advanced motion, the subducting velocity is 2 cm/yr. So there might be a slow-subduction effect that conduction is more complete. This might address the temperature is warmer everywhere.

[lhy11009]

Addressed tasks:

• Change all the "plate velocity" to "spreading velocity" in mass_converving.cc and in all .wb files.
• Explicitly add "subducting velocity" to every .wb file that uses the "mass conserving" model.
• Add change log.

Any other thoughts?

[danieldouglas92]

@lhy11009 I think to isolate the effect of trench advancement/retreat you could make the subducting velocity the same in both cases, but make the spreading velocity different. This will make the heat conduction the same, and if you make it so that the age of the plate at the trench is the same in both cases (by moving the location of the ridge based on the spreading rate), this would effectively create a case where you are only capturing the temperature difference based on trench retreat/trench advancement.

[lhy11009]

@lhy11009 I think to isolate the effect of trench advancement/retreat you could make the subducting velocity the same in both cases, but make the spreading velocity different. This will make the heat conduction the same, and if you make it so that the age of the plate at the trench is the same in both cases (by moving the location of the ridge based on the spreading rate), this would effectively create a case where you are only capturing the temperature difference based on trench retreat/trench advancement.

Good point, I'll generate these new tests.

[lhy11009]

In these tests, I prescribe the same subducting velocity of 10 cm/yr and different spreading velocities (5cm/yr - retreat, 10 cm/yr - static, 15 cm/yr - advance).
Meanwhile, the same age is maintained at the trench by varying the length of the plate (4000 km, 8000 km, 12000 km).
Here are the results:

Case with a static trench (slab_temperature_trench_static):

Case with a retreating trench (slab_temperature_trench_retreat_1):

Case with an advancing trench (slab_temperature_trench_advance_1):

Trending in the temperature field: retreat < static < advance, which matches the change of the age on the plate (retreat has the oldest age on the plate from a small spreading rate).

[codecov]

Codecov Report

All modified and coverable lines are covered by tests ✅

Project coverage is 98.21%. Comparing base (c547b29) to head (33c703e).
Report is 8 commits behind head on main.

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #694      +/-   ##
==========================================
+ Coverage   98.19%   98.21%   +0.02%     
==========================================
  Files         107      107              
  Lines        7464     7462       -2     
==========================================
  Hits         7329     7329              
+ Misses        135      133       -2     

Files	Coverage Δ
...ucting_plate_models/temperature/mass_conserving.cc	99.59% <100.00%> (-0.01%)	⬇️
source/world_builder/utilities.cc	97.22% <100.00%> (ø)

... and 3 files with indirect coverage changes

---

Continue to review full report in Codecov by Sentry.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update c547b29...33c703e. Read the comment docs.

[danieldouglas92]

Very cool!!!

[lhy11009]

I see there is some coverage issue. I will fix this shortly

[lhy11009]

@MFraters, I have fixed two spots of the code coverage issue. But there is still one in parameters.cc. I haven't changed any code there, can you take a look at why that code coverage decreases?

[lhy11009]

@MFraters Can you take a look at this again ?

[MFraters]

I am sorry it took so long for me to take a look at it. The code looks great to me, just a few small comments and questions for clarification. @danieldouglas92, could you approve this pull request as well if you think it it ready to merge?

The piece of code in paramters.cc you are refering to is filling in the default value if needed. It might be that you added a value in one of the tests, which first would have taking a default value. Don't worry to much about it for now.

[MFraters]

can you align this?

[lhy11009]

@MFraters I addressed all your comments.

[MFraters]

Why is the subducting velocity sticking out so strangely in the .wb file compared to the other lines?

[lhy11009]

It's fine on my computer. I think there is a problem converting \t to " ". For some reason these files are not covered by astyle.
And this is apparant when I am viewing changes in files, but it aligns well when I click "edit" and enter a page of editor.
I think I have fixed it by searching for all the tabs in the file and replacing them with 2 blanks. See if this makes sense.

[danieldouglas92]

Looks good to me! Super exciting that this is a feature now

[MFraters]

It's fine on my computer. I think there is a problem converting \t to " ". For some reason these files are not covered by astyle. And this is apparant when I am viewing changes in files, but it aligns well when I click "edit" and enter a page of editor. I think I have fixed it by searching for all the tabs in the file and replacing them with 2 blanks. See if this makes sense.

yes, the .wb files are not processed by astyle. I didn't think it would make sense to do that, but maybe it would be useful. Anyway thanks for your work on this!

[lhy11009]

@mibillen @rfildes We have completed adding ridge migration and trench retreat to the slab temperature in WB

[lhy11009]

I have added a few tests in this PR. You can scroll upward and figure it out from the prm and picture.