# "Shrinking earth" from dynamic topography plugin

**URL:** <https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572>\
**Category:** ASPECT\
**Created:** [September 12, 2022, 6:32pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572 "2022-09-12T18:32:53Z")\
**Posts on this page:** 14\
**Page:** 2

<div class="post-metadata">

**Author:** ![bobmyhill](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/bobmyhill/32/38_2.png) [@bobmyhill](https://community.geodynamics.org/u/bobmyhill)\
**Post date:** [September 28, 2022, 8:21pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/21 "2022-09-28T20:21:06Z")

</div>

Great news! Thanks for letting us know 🙂

---

<div class="post-metadata">

**Author:** ![goldberg-s](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/goldberg-s/32/982_2.png) [@goldberg-s](https://community.geodynamics.org/u/goldberg-s)\
**Post date:** [October 5, 2022, 7:34pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/22 "2022-10-05T19:34:11Z")

</div>

A puzzling twist: I used an input thermal model with a lithosphere that goes smoothly to 273K at the surface over a distance that is well-resolved by the model mesh, set the surface temperature boundary condition to 273K, and the problem came back. But when I set the surface temperature boundary condition to 1573, even though the input thermal structure has 273 at the surface, the problem went away again. So the problem is solved when I set the bottom and top boundary temperatures to 1573, regardless of what the input thermal structure has at the boundaries.

Other than the input .txt file for the thermal structure and the `Outer temperature` boundary condition, the two other places I specify 1573 are the `Adiabatic surface temperature` to compute the reference profile, and the `Reference temperature` in the `Visco Plastic` material model subsection. I tried changing the `Adiabatic surface temperature` to 273, but that didn’t solve the problem, and did make the solver take a lot longer since that’s far away from the ambient mantle temperature.

---

<div class="post-metadata">

**Author:** ![jbnaliboff](https://avatars.discourse-cdn.com/v4/letter/j/58956e/32.png) [@jbnaliboff](https://community.geodynamics.org/u/jbnaliboff)\
**Post date:** [October 6, 2022, 2:47am UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/23 "2022-10-06T02:47:35Z")

</div>

Hi Sam,

Odd indeed! Can you post images of the boundary layer viscosity and temperature with the mesh also visible?

As in previous tests, does decreasing the min/max viscosity (or using a constant viscosity) remove the issue?

Cheers,  
John

---

<div class="post-metadata">

**Author:** ![goldberg-s](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/goldberg-s/32/982_2.png) [@goldberg-s](https://community.geodynamics.org/u/goldberg-s)\
**Post date:** [October 6, 2022, 3:59pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/24 "2022-10-06T15:59:32Z")

</div>

Hi John,

Yes of course - here are some Paraview screenshots showing the lithosphere, for the model run that did not have the “shrinking earth.” You can see the smoothly decreasing temperature towards the surface within the lithosphere, but there’s a 1-element boundary layer that is created by the 1573 K `Outer temperature` boundary condition. The input thermal structure has 273 K at the surface. When I set the `Outer temperature` to 273, the thermal structure looks as it should, with the temperature smoothly going to 273 at the surface. But this makes the “shrinking earth” come back.

Best,  
Sam

![image1](https://canada1.discourse-cdn.com/flex036/uploads/geodynamics/original/2X/8/83f2ff33926a30510d57a2a030b47f09c0d3157b.png)

 ![image2](https://canada1.discourse-cdn.com/flex036/uploads/geodynamics/original/2X/d/d75834b347737484118c27f8ed7e34d86c4a69c7.jpeg)  
 ![image3](https://canada1.discourse-cdn.com/flex036/uploads/geodynamics/original/2X/1/1cc27c93bbeb614029466b686d9e18ab4aa5076f.png)

---

<div class="post-metadata">

**Author:** ![jbnaliboff](https://avatars.discourse-cdn.com/v4/letter/j/58956e/32.png) [@jbnaliboff](https://community.geodynamics.org/u/jbnaliboff)\
**Post date:** [October 6, 2022, 7:08pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/25 "2022-10-06T19:08:53Z")

</div>

Hi Sam,

Indeed, the boundary layer is certainly well resolved. In this case my guess is the issue lies within the postprocessor or how it is being used with this specific model.

That is to say unless there are other indicators (e.g., viscosity or velocity anomalies, poor solver convergence), I don’t see anything physically wrong with the initial conditions for the model with the surface temperature at 273 K.

My recollection is that many studies calculate dynamic topography using stresses from the base of the lithosphere, rather than at the surface. This would require changing where the dynamic topography post processor calculates dynamic topography, or manually calculating the values after extracting the stresses at a specific depth, isotherm, etc.

Is there a prior study that is similar in design to what you are trying to achieve here, which we could try benchmarking against?

Cheers,  
John

---

<div class="post-metadata">

**Author:** ![jaustermann](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/jaustermann/32/207_2.png) [@jaustermann](https://community.geodynamics.org/u/jaustermann)\
**Post date:** [October 10, 2022, 1:18pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/26 "2022-10-10T13:18:52Z")

</div>

Hi Sam,

Can you remind me whether this happens only if you use the viscoplastic rheology or also for the linear rheology? If it does happen with the linear rheology as well, we should have some runs that you could compare things to to hopefully figure our what’s going on.

Best,  
Jacky

---

<div class="post-metadata">

**Author:** ![goldberg-s](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/goldberg-s/32/982_2.png) [@goldberg-s](https://community.geodynamics.org/u/goldberg-s)\
**Post date:** [October 10, 2022, 1:39pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/27 "2022-10-10T13:39:18Z")

</div>

Hi Jacky,

It only happens with the `Visco Plastic` rheology; it works fine with the linear rheology. I’m working on making my own Python post-processor to extract the pressure at the base of the lithosphere to compute dynamic topography that way and see what results that gives.

Best,  
Sam

---

<div class="post-metadata">

**Author:** ![jaustermann](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/jaustermann/32/207_2.png) [@jaustermann](https://community.geodynamics.org/u/jaustermann)\
**Post date:** [October 10, 2022, 2:11pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/28 "2022-10-10T14:11:23Z")

</div>

Hi Sam,

Okay sounds good! One of my postdocs extracted the deviatoric stress tensors and the associated traction field at the base of the lithosphere using existing postprocessors to then use that as a boundary condition in a separate lithosphere deformation code (in [this](https://www.nature.com/articles/s41467-022-31903-2) paper). Happy to connect you to him if that is useful at any point.

Best,  
Jacky

---

<div class="post-metadata">

**Author:** ![goldberg-s](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/goldberg-s/32/982_2.png) [@goldberg-s](https://community.geodynamics.org/u/goldberg-s)\
**Post date:** [October 31, 2022, 1:52pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/29 "2022-10-31T13:52:54Z")

</div>

Well it took me a while to get it all figured out, but I extracted the pressure from the base of the lithosphere, subtracted the lithostatic pressure computed from the thermal/density profile, and divided the difference by rho\*g to get (air-loaded) dynamic topography. The dynamic topography field computed this way is pretty much identical to what comes out of the ASPECT post-processor (which uses surface stresses), so that’s good! I think the correspondence is because of the strong lithosphere, which transmits stresses from base to surface without much internal deformation. Unfortunately, the “shrinking earth” problem persists with this method of dynamic topography calculation as well, with the average rate of change is about the same as from the ASPECT post-processor. So this means that it’s not some quirk of the post-processor, but something more fundamental to the solution.

If I subtract the mean from the DT field at each timestep, to explicitly force a zero average rate-of-change, is that kosher? Or would I be missing something physical if I do that?

---

<div class="post-metadata">

**Author:** ![jbnaliboff](https://avatars.discourse-cdn.com/v4/letter/j/58956e/32.png) [@jbnaliboff](https://community.geodynamics.org/u/jbnaliboff)\
**Post date:** [November 2, 2022, 4:45am UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/30 "2022-11-02T04:45:27Z")

</div>

Hi Sam,

Great to hear that the new approach gives the same result at the ASPECT post processor.

Are you by chance available to join the ASPECT bi-weekly user meetings on Wednesdays (8-9 am Pacific)? The next meeting will be on November 9th. I think it would help to talk through all of this live and get input from others.

> If I subtract the mean from the DT field at each timestep, to explicitly force a zero average rate-of-change, is that kosher? Or would I be missing something physical if I do that?

My initial thought is that with sufficiently stricter solver tolerances, the global rate of change should be very close to zero. Again, I think we have reached the point where a live recap with the broader group would help 🙂

Cheers,  
John

---

<div class="post-metadata">

**Author:** ![Malley1](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/malley1/32/1366_2.png) [@Malley1](https://community.geodynamics.org/u/Malley1)\
**Post date:** [January 31, 2023, 12:20pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/31 "2023-01-31T12:20:28Z")

</div>

Dear All,

Was there any update on this issue from the meeting? I am encountering a similar phenomenon with TERRA model output.  
I was hopeful that taking stress at the base of the lithosphere would be the solution, and that the cooling lithosphere might be responsible for the net negative DT, so would like to hear if there was consensus on why that didn’t work!

Best,

Conor O’Malley

---

<div class="post-metadata">

**Author:** ![goldberg-s](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/goldberg-s/32/982_2.png) [@goldberg-s](https://community.geodynamics.org/u/goldberg-s)\
**Post date:** [January 31, 2023, 2:06pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/32 "2023-01-31T14:06:53Z")

</div>

Hi Conor,

I had a baby at the beginning of November, which put a pause on my scientific progress. So I have not done any further digging since the last messages on here. But I’m intrigued that you’re having a similar issue with a different code. What sort of rheology are you using in your models?

Best,  
Sam

---

<div class="post-metadata">

**Author:** ![jbnaliboff](https://avatars.discourse-cdn.com/v4/letter/j/58956e/32.png) [@jbnaliboff](https://community.geodynamics.org/u/jbnaliboff)\
**Post date:** [January 31, 2023, 3:47pm UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/33 "2023-01-31T15:47:17Z")

</div>

Hi all,

First, congratulations Sam 🙂

Second, Conor thank you for posting to the forum and welcome. Despite not knowing the underlying root cause for the issue, perhaps somewhat encouraging that the phenomenon is not just specific to ASPECT.

Perhaps a path forward would be to setup a common benchmark run with ASPECT, TERRA, and CitComS and see if there is a correlation between negative DT and code/numerical method?

Cheers,  
John

---

<div class="post-metadata">

**Author:** ![Malley1](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/malley1/32/1366_2.png) [@Malley1](https://community.geodynamics.org/u/Malley1)\
**Post date:** [April 12, 2023, 9:35am UTC](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572/34 "2023-04-12T09:35:42Z")

</div>

Hi John & Sam,

Congrats Sam!  
We were observing the offset for visco-elastic runs with & without temperature-dependent viscosity. I am not running the models myself; they’re being run by Huw Davies and James Panton in Cardiff. Having spoken to them, I think the issue was we weren’t using the appropriate reference density profile to calculate normal stresses; we have resolved that issue now and get ~0 mean DT. To calculate DT, the TERRA models re-run the requested timestep instantaneously, relax plate motion constraints and impose free slip surface boundary conditions, then calculate radial normal stresses. Those are then converted to DT via simple isostatic scaling. Is that similar to how the ASPECT post-processor calculates DT?

I think benchmarking has been done before to some degree by e.g. Liu & King (2019), who benchmarked DT kernels in CitComS, ASPECT and TERRA. However, they published min/max DT and not mean values. I may have a more careful read and see if we could add anything to their results, what do you think?

Best,

Conor

[Previous page](https://community.geodynamics.org/t/shrinking-earth-from-dynamic-topography-plugin/2572.md?page=1)
