[aspect-devel] Internal heating in aspect (Ludovic Jeanniot)

Jeanniot, L. (Ludovic) l.jeanniot at uu.nl
Sun Aug 26 08:03:58 PDT 2018


Rene, Scott, Ian and Max,

I am currently benchmarking thermal convection in a 3D hollow sphere following Arrial et al. 2014 setup to compare the obtained results (average temperature, V_rms, Nusselt number and thermal shape at steady-state) in ASPECT with CitcomS and a RBF-PS model (Arrial et al. 2014). The whole setup is made to be nondimensional. 
I did encountered similar issues of lowered average temperature than expected at low resolution (GMR = 3 (level 4)) and found a certain solution although I still need to understand why now it works. 

     The first issue is about the default mesh:
1.- The default mesh is (I think) depth-dependent as the mesh is finer beneath a certain depth and coarser above. What is happening is that at low resolution, the temperature solution at the boundary between the two finer/coarser meshes is misbehaving (under-solved?). In CitcomS, the mesh is composed of about 12 'caps' of identical number of cells and I therefore changed the parameter "Cells along circumference" to 12 in the geometry model subsection to reproduce a similar and more uniform grid. But I must admit that what Cells along circumference does remain obscure to me. 

     The second issue may be what you are looking for:
2.- I attached some experiments with changing the sensitivity of the beta factor in the entropy viscosity stabilization and the cR factor in the artificial viscosity stabilization (in the manual: A.41 Parameters in section Discretization/Stabilization parameters). By default in ASPECT, the values of the two factors are 0.078 and 0.33 respectively. In the figure, the 3 blue lines, which were obtained when turning the two factor values to 0, or one of the two to 0, are identical and they correspond at what I should get (note the level 4 resolution). Keeping at low resolution, by default with two factor values beta = 0.078 and cR = 0.33 results in dropping the average temperature from 0.218 to 0.19. Running the models with higher resolution then results in the average temperature I was expecting (similar to CitcomS and RBF-PS model) 

Maybe this can help too? 

Ludovic


-----Original Message-----
From: Aspect-devel [mailto:aspect-devel-bounces at geodynamics.org] On Behalf Of aspect-devel-request at geodynamics.org
Sent: 25 August 2018 22:17
To: aspect-devel at geodynamics.org
Subject: Aspect-devel Digest, Vol 81, Issue 7

Send Aspect-devel mailing list submissions to
	aspect-devel at geodynamics.org

To subscribe or unsubscribe via the World Wide Web, visit
	http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
or, via email, send a message with subject or body 'help' to
	aspect-devel-request at geodynamics.org

You can reach the person managing the list at
	aspect-devel-owner at geodynamics.org

When replying, please edit your Subject line so it is more specific than "Re: Contents of Aspect-devel digest..."


Today's Topics:

   1. Re: Internal heating in aspect (Scott King)
   2. Re: Internal heating in aspect (Ian Rose)


----------------------------------------------------------------------

Message: 1
Date: Sat, 25 Aug 2018 15:50:25 -0400
From: Scott King <sdk at vt.edu>
To: aspect-devel <aspect-devel at geodynamics.org>
Subject: Re: [aspect-devel] Internal heating in aspect
Message-ID: <1E881A37-48D1-46FA-A2F1-9FEF1FAC8308 at vt.edu>
Content-Type: text/plain; charset="utf-8"

Hi Rene, Max, Wolfgang, and Diogo;

Grant found when looking at the Zhong 2008 “Benchmark” calculations that he gets better agreement with mean temperature and mean RMS velocity than he does with surf/cmb Nusselt number.  We find generally better results with a nearly uniform radial mesh than the aspect default mesh pretty much across the board.  As we go to higher and higher resolution (5 or 6 refinements), we get convergence of the surf and cmb heatflow and good agreement with other codes.  We are puzzling about it because the mean temperature and mean velocity profiles (fn of r) are spot on with the ones from CitcomS, so the fluxes should really be correct.  The numerical diffusion still being included sounds like a possibility. On the other hand, that explanation puzzles me because the 2D Cartesian cases are spot on.  It made me wonder if there is a scaling issue that only shows up in the spherical case?  On my list of things to do is to check into that but life interrupted my to do list.

Interestingly we get a significantly higher mean temperatures (essentially the profile shifted by 10%) for the C cases (Ra=10^5) and we’ve never quite figured out why.  That sort of concerns me because all the recent spherical papers I’ve seen only address the Ra=7000 cases and assume this tell us enough about how these codes behave at higher Ra, which seems to me a bit optimistic.   We’ve also never seen any codes try to reproduce the 10^5 results or to use more than a single grid either :)   Of course that makes me wonder which code is correct. (LOL)  Convergence studies seem to be out of favor.

Scott

> On Aug 25, 2018, at 1:58 PM, Rene Gassmoeller <rene.gassmoeller at mailbox.org> wrote:
> 
> Hi Max,
> 
> let me cc this to the mailing list, as it might be of interest to others.
> I found this old email thread on the mailing list that seems to discuss the same issue: http://lists.geodynamics.org/pipermail/aspect-devel/2012-November/000138.html <http://lists.geodynamics.org/pipermail/aspect-devel/2012-November/000138.html> (there are some additional messages when you scroll through the thread). In short: At the time Ian found that this is a problem when the boundary layer is underresolved. Since the element at the surface are larger than the elements at the CMB they do resolve a different temperature gradient for the boundary layer if the boundary layer is not sufficiently resolved, and we compute the heat flux based on the temperature gradient.
> However, thinking about it again, I do not think this can be the only reason. If the difference were caused by numerical resolution alone, the average temperature should change, until the two heat fluxes are in balance as you mentioned. But another idea I had was that the computed heat flux is a postprocessor that only considers physical diffusion (i.e. diffusion due to thermal conductivity), while in the assembly we also include the numerical diffusion to stabilize the equation. So maybe the difference between the heat fluxes is caused by numerical diffusion (and it is so large because the cells are relatively coarse)? If that were the case the difference should increase/decrease with higher/lower resolution, as the stabilization scales with the element size. Unfortunately that is not a criterion to distinguish between the two possible causes (obviously resolution of the temperature gradient across the boundary layer also scales with element size). You could try to re-run the model with discontinuous temperature elements, that should disable the numerical diffusion. Let me know if that makes sense.
> Best,
> 
> Rene
> 
> On 08/24/2018 06:46 PM, Max Rudolph wrote:
>> Rene and Wolfgang,
>> We're running some models in 3D spherical geometry in aspect with internal heating and are unable to reproduce thermal histories consistent with what we saw in identically configured calculations in CitcomS. It is possible that there is some aspect of our model setup that is inconsistent, but we are using the same depth profile of viscosity, same boundary conditions, and same initial conditions. The ASPECT models cool much faster (which also drives up viscosity). We are trying to get understand these results, and I was looking at the heat transport in the 3D spherical model discussed in section 5.3.2 of the aspect manual 'simple convection in the 3D spherical shell'. In figure 40, it appears that this calculation has been run to nearly statistically steady state. Is this correct? I am wondering, basically, how to interpret the imbalance between reported surface and CMB heat flow. To follow up on this, I ran the cookbook .prm file from 5.3.2 (shell_simple_3d) for 5 Gyr, which appears to be very close to or at statistically steady state at the end of the calculation. Due to resource constraints, I only ran at global refinement level 3 (definitely under-resolved, mea culpa). However, at the end of the calculation, there is still a significant imbalance between the surface and core mantle boundary heat flows. This is a model with no internal heating, so at statistically steady state, the surface and cmb heat flows should balance to conserve energy. Do you have any insight into why we might be seeing this behavior? I've attached the input .prm file, the statistics file, and a plot showing the average temperature and surface and cmb heat flows. Thank you for thinking about this!
>> 
>> Best
>> Max
>> 
>> p.s. I sent this to both of you because it looked to me like Wolfgang contributed the cookbook. Please let me know if there is someone else that might be better positioned to help with this.
>> 
>> <image.png>
>> 
>> 
>> 
> 
> --
> Rene Gassmoeller
> https://gassmoeller.github.io/ 
> <https://gassmoeller.github.io/>______________________________________
> _________
> Aspect-devel mailing list
> Aspect-devel at geodynamics.org
> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.geodynamics.org/pipermail/aspect-devel/attachments/20180825/b72f5484/attachment-0001.html>

------------------------------

Message: 2
Date: Sat, 25 Aug 2018 16:17:13 -0400
From: Ian Rose <ian.rose at berkeley.edu>
To: aspect-devel <aspect-devel at geodynamics.org>
Subject: Re: [aspect-devel] Internal heating in aspect
Message-ID:
	<CAEm=8TR0FX7wiMDYXGh7zWtVouQ0cgH1GYnJgNQzz7JPnjQfmg at mail.gmail.com>
Content-Type: text/plain; charset="utf-8"

It is also possible that the accuracy could be improved by using the consistent-boundary-flux method (CBF) for calculating the fluxes in postprocessing. The CBF can be seen as a form of projection for the solution into the flux space. It essentially constructs a new linear system, the solution to which is the fluxes at the boundaries. That system can be made to be diagonal, so it is cheap to solve.

I found that using CBF in calculating the dynamic topography significantly improved the convergence of the solution (quadratic rather than linear, if memory serves). The mathematics for computing the heat flux is essentially identical, so a lot of the dynamic topography code could be copied whole-cloth into the heat flux postprocessor, in case anyone wants to give it a shot.

Best,
Ian

On Sat, Aug 25, 2018 at 4:01 PM Scott King <sdk at vt.edu> wrote:

> Hi Rene, Max, Wolfgang, and Diogo;
>
> Grant found when looking at the Zhong 2008 “Benchmark” calculations 
> that he gets better agreement with mean temperature and mean RMS 
> velocity than he does with surf/cmb Nusselt number.  We find generally 
> better results with a nearly uniform radial mesh than the aspect 
> default mesh pretty much across the board.  As we go to higher and 
> higher resolution (5 or 6 refinements), we get convergence of the surf 
> and cmb heatflow and good agreement with other codes.  We are puzzling 
> about it because the mean temperature and mean velocity profiles (fn 
> of r) are spot on with the ones from CitcomS, so the fluxes should 
> really be correct.  The numerical diffusion still being included 
> sounds like a possibility. On the other hand, that explanation puzzles 
> me because the 2D Cartesian cases are spot on.  It made me wonder if 
> there is a scaling issue that only shows up in the spherical case?  On 
> my list of things to do is to check into that but life interrupted my to do list.
>
> Interestingly we get a significantly higher mean temperatures 
> (essentially the profile shifted by 10%) for the C cases (Ra=10^5) and 
> we’ve never quite figured out why.  That sort of concerns me because 
> all the recent spherical papers I’ve seen only address the Ra=7000 
> cases and assume this tell us enough about how these codes behave at higher Ra, which seems to me a bit
> optimistic.   We’ve also never seen any codes try to reproduce the 10^5
> results or to use more than a single grid either :)   Of course that makes
> me wonder which code is correct. (LOL)  Convergence studies seem to be 
> out of favor.
>
> Scott
>
> On Aug 25, 2018, at 1:58 PM, Rene Gassmoeller < 
> rene.gassmoeller at mailbox.org> wrote:
>
> Hi Max,
>
> let me cc this to the mailing list, as it might be of interest to others.
>
> I found this old email thread on the mailing list that seems to 
> discuss the same issue:
> http://lists.geodynamics.org/pipermail/aspect-devel/2012-November/0001
> 38.html (there are some additional messages when you scroll through 
> the thread). In
> short: At the time Ian found that this is a problem when the boundary 
> layer is underresolved. Since the element at the surface are larger 
> than the elements at the CMB they do resolve a different temperature 
> gradient for the boundary layer if the boundary layer is not 
> sufficiently resolved, and we compute the heat flux based on the temperature gradient.
>
> However, thinking about it again, I do not think this can be the only 
> reason. If the difference were caused by numerical resolution alone, 
> the average temperature should change, until the two heat fluxes are 
> in balance as you mentioned. But another idea I had was that the 
> computed heat flux is a postprocessor that only considers physical 
> diffusion (i.e. diffusion due to thermal conductivity), while in the 
> assembly we also include the numerical diffusion to stabilize the 
> equation. So maybe the difference between the heat fluxes is caused by 
> numerical diffusion (and it is so large because the cells are 
> relatively coarse)? If that were the case the difference should 
> increase/decrease with higher/lower resolution, as the stabilization 
> scales with the element size. Unfortunately that is not a criterion to 
> distinguish between the two possible causes (obviously resolution of 
> the temperature gradient across the boundary layer also scales with 
> element size). You could try to re-run the model with discontinuous 
> temperature elements, that should disable the numerical diffusion. Let me know if that makes sense.
>
> Best,
>
> Rene
>
> On 08/24/2018 06:46 PM, Max Rudolph wrote:
>
> Rene and Wolfgang,
> We're running some models in 3D spherical geometry in aspect with 
> internal heating and are unable to reproduce thermal histories 
> consistent with what we saw in identically configured calculations in 
> CitcomS. It is possible that there is some aspect of our model setup 
> that is inconsistent, but we are using the same depth profile of 
> viscosity, same boundary conditions, and same initial conditions. The 
> ASPECT models cool much faster (which also drives up viscosity). We 
> are trying to get understand these results, and I was looking at the 
> heat transport in the 3D spherical model discussed in section 5.3.2 of 
> the aspect manual 'simple convection in the 3D spherical shell'. In 
> figure 40, it appears that this calculation has been run to nearly 
> statistically steady state. Is this correct? I am wondering, 
> basically, how to interpret the imbalance between reported surface and 
> CMB heat flow. To follow up on this, I ran the cookbook .prm file from 
> 5.3.2
> (shell_simple_3d) for 5 Gyr, which appears to be very close to or at 
> statistically steady state at the end of the calculation. Due to 
> resource constraints, I only ran at global refinement level 3 
> (definitely under-resolved, mea culpa). However, at the end of the 
> calculation, there is still a significant imbalance between the 
> surface and core mantle boundary heat flows. This is a model with no 
> internal heating, so at statistically steady state, the surface and 
> cmb heat flows should balance to conserve energy. Do you have any 
> insight into why we might be seeing this behavior? I've attached the 
> input .prm file, the statistics file, and a plot showing the average temperature and surface and cmb heat flows.
> Thank you for thinking about this!
>
> Best
> Max
>
> p.s. I sent this to both of you because it looked to me like Wolfgang 
> contributed the cookbook. Please let me know if there is someone else 
> that might be better positioned to help with this.
>
> <image.png>
>
>
>
>
> --
> Rene Gassmoellerhttps://gassmoeller.github.io/
>
> _______________________________________________
> Aspect-devel mailing list
> Aspect-devel at geodynamics.org
> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
>
>
> _______________________________________________
> Aspect-devel mailing list
> Aspect-devel at geodynamics.org
> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.geodynamics.org/pipermail/aspect-devel/attachments/20180825/41f8abc5/attachment.html>

------------------------------

Subject: Digest Footer

_______________________________________________
Aspect-devel mailing list
Aspect-devel at geodynamics.org
http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel

------------------------------

End of Aspect-devel Digest, Vol 81, Issue 7
*******************************************
-------------- next part --------------
A non-text attachment was scrubbed...
Name: test stabilization.png
Type: image/png
Size: 6460 bytes
Desc: test stabilization.png
URL: <http://lists.geodynamics.org/pipermail/aspect-devel/attachments/20180826/21f40642/attachment-0001.png>


More information about the Aspect-devel mailing list