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

Elbridge Gerry Puckett egp at math.ucdavis.edu
Wed Aug 29 14:10:09 PDT 2018


Hi Everyone,

I'm sorry to come to this discussion so late.  I have been meaning to 
read the entire sequence of posts before posting myself.

Ying He implemented DG with a Bound Preserving (BP) Limiter for the 
temperature equation (with incomprehensible flow).  We (Ying, Magali, 
Louise, and myself) submitted a poster to the 2016 AGU fall meeting and 
it was chosen for a talk, which Ying gave.

I have the talk slides and they appear to be relevant to to this 
discussion, even though we did not study internal heating.  Again, I 
would have liked to read the entire discussion before posting, but Max 
just mentioned DG and hence it seems appropriate for me to offer Ying's 
talk slides to the [aspect-devel] mailing list as well as the rest of 
the community.

A few words about what Ying / we did.  Ying was interested in studying 
the effectiveness of this new methodology as a function of the  Péclet 
number (Pe) number.  For many of the problems we were computing as a 
group, such as our recent LLSVP paper in the special issue of PEPI to 
associated with the 15th SEDI conference, we (Ying) determined that at 
that  Péclet number the changes to our results were negligible when one 
used the original algorithm in ASPECT with EV and the new DGBP algorithm.

However, for other problems, Magali had an argument that for some 
problems one might want to compute in ASPECT there would be an 
/effective local Péclet number, /for which there are decided differences 
between the two algorithms.  This difference is shown in the talk 
slides.  (There is a movie of this, which I can probably find.)

I don't see this material in the manual, which is unfortunate. I can put 
this on my to-do list if others think it would be helpful.

Two things to note:

     1. Part of Ying's to-do list was to add the ability to use AMR to
        algorithm with this DGBP temperature advection-diffusion
        algorithm, which I think remains to be done.

     2. Another issue to note is that that Ying found it necessary to
        use a second-order accurate Strong Stability Preserving time
        stepping algorithm (SSP2) instead of the usual BDF2 for the DGBP
        algorithm. This is from Ying's 2017 PEPI paper with Magali and me.

            "For the case when the FEM is used for the spatial
            discretization, we use the same semi-implicit BDF2 scheme as
            shown in Kronbichler et al. (2012) except we fix the
            time-step size dt for simplicity. However, if the DG method
            is used, then a Strong Stability-Preserving (SSP) high order
            time discretization method (Gottlieb, 2005; Gottlieb et al.,
            2001) is applied in order to maintain the bound preserving
            property from the limiter."

This does not seem to be mentioned in the talk slides, probably due to 
lack of time.  A search on SSP2 in the ASPECT manual generated June 22, 
2018 yields zero hits and and even a search on DGBP yields only three 
hits in thre parameter section.  There is probably room for more 
documentation there. : -(

Anyway, the point is that one might keep (1) and (2) in mind while using 
DGBP for the temperature equation.

In summary, I didn't intend to write such a long post, I hope some of 
you find it is worthwhile.

Shall I post the talk slides and if so, where?  I can put them on 
math.ucdavis.edu/~egp/... or post to this discussion with the slides 
appended, or ..

Sincerely,

- Gerry


On 08/29/2018 12:49 PM, Max Rudolph wrote:
> Rene,
> Thank you for summarizing the current state of affairs. I would 
> suggest that when you make the changes in (1)-(2), you also make 
> ASPECT print at least a warning when the entropy viscosity exceeds, 
> say, 1% of the thermal conductivity. It might be even better for the 
> code to exit in with an error unless a parameter is explicitly set to 
> ignore such a condition.
>
> It also seems that one way around (1) is to use DG elements for 
> composition. Is this what most users are doing anyways?
>
> Thank you for addressing these issues so quickly. I am happy to help 
> with testing as needed.
>
> Best,
> Max
>
> On Wed, Aug 29, 2018 at 12:16 PM Rene Gassmoeller 
> <rene.gassmoeller at mailbox.org <mailto:rene.gassmoeller at mailbox.org>> 
> wrote:
>
>     Hi all,
>
>     I think by now we have a pretty good understanding of the problem,
>     however there is no clear path forward yet. I have done some
>     further testing with the shell_simple_3d cookbook:
>
>      - As Max pointed out the artificial viscosity at resolution <= 2
>     global refinement in spherical models is bigger than the natural
>     diffusion, although it reduces drastically with resolution
>     (compare to a natural conductivity of ~4 W/m*K).
>
>     Global refinement / Maximum artificial viscosity timestep 0 /
>     timestep 2:
>
>     1 (4 radial elements): 94.1 W/m*K  /  59.9 W/m*K
>
>     2 (8 radial elements): 48.9 W/m*K / 6.12 W/m*K
>
>     3 (16 radial elements): 24.8 W/m*K / 0.92 W/m*K
>
>     There are some pitfalls here, in particular the artificial
>     viscosity in the output of timestep 0 is in general much bigger
>     than in later timesteps (and it only scales linearly with cell
>     size), because we do not have a velocity solution yet. Timestep 1
>     is somewhat unreliable as well, because the BDF2 scheme is not
>     fully initialized. Nevertheless, for later timesteps the
>     artificial viscosity scales at least with cell size squared (h^2)
>     as expected (actually slightly better, because it is not only h^2,
>     but also based on the residual of the equation, which reduces as
>     well, ideally with h^3). This means higher resolutions should
>     significantly reduce the total amount of artificial diffusion,
>     although it might still be significant (>1% of total diffusion) up
>     to resolution 4 or so, which seems unacceptably expensive.
>     Moreover, as Scott and Cedric mentioned there will always be an
>     imbalance between the heat fluxes at the top and bottom boundary
>     with this method, unless we use radially uniformly spaced
>     elements, or disable the stabilization, although the magnitude of
>     the imbalance should reduce with increased resolution. Also the
>     thermal conductivity within the domain will be unequally
>     distributed (artificial conductivity in the upper mantle will be
>     higher than in the lower mantle).
>
>     While the above point suggests that the artificial viscosity
>     scheme is correctly implemented, there still exist three points
>     for possible improvements:
>
>      - As Wolfgang mentioned the parameters that were originally
>     chosen for the artificial viscosity scheme were different from
>     what they are today, because we noticed that the original
>     smoothing was too weak to stabilize oscillations in compositional
>     fields. It is completely possible that the original values were
>     sufficient for the temperature (which already has natural
>     diffusion), and in hindsight I should have thought of just using
>     different parameters for composition and temperature (though that
>     was 6 years ago, and I was just a 2nd year PhD student at the time
>     I worked on modifying the parameters for stabilizing the
>     composition equation). I can easily make a change that allows for
>     different parameters for temperature and composition, which would
>     allow everyone to test their favorite values, without risking
>     oscillations in compositional fields.
>
>     -  Even if we can reduce the parameters it is of course possible
>     that as Max pointed out the anisotropic nature of the diffusion in
>     SUPG is more appropriate / "less wrong" than the isotropic entropy
>     viscosity for our problem. As Juliane and Timo pointed out we
>     could reimplement the content of
>     https://github.com/geodynamics/aspect/pull/412 and see if SUPG
>     performs better for low resolutions. Does anyone know if the
>     "artificial viscosity" of SUPG also scales with h^2? Because if it
>     is linear, we might implement something that is better for low
>     resolutions, but worse for high resolutions.
>
>     - As Max mentioned: We should not need a stabilization for pure
>     conduction problems (where velocity is 0), and should modify the
>     algorithm accordingly.
>
>
>     So the next steps could be:
>
>     1. Allow for different stabilization parameters for temperature
>     and composition, and check which values are still stable.
>     2. Do not stabilize advection/diffusion solutions where the
>     velocity is zero (because it is only a diffusion equation).
>     3. Reimplement the SUPG based on
>     https://github.com/geodynamics/aspect/pull/412 and see how it
>     performs (at low and high resolutions).
>
>     Does that summarize our discussion appropriately? I can easily
>     make the code adjustments 1 and 2 (they are easy and useful in any
>     case), and could also look into creating an initial version of 3
>     (although it would take a bit of time), but I currently do not
>     have the time for much testing of the methods, so I would be
>     greatful if someone else could do the testing and benchmarking of
>     the methods.
>
>     Best,
>     Rene
>
>
>     On 08/29/2018 09:19 AM, Max Rudolph wrote:
>>     On Tue, Aug 28, 2018 at 8:01 PM Wolfgang Bangerth
>>     <bangerth at colostate.edu <mailto:bangerth at colostate.edu>> wrote:
>>
>>         On 08/28/2018 05:33 PM, Max Rudolph wrote:
>>         >  From this, it is very obvious why the solution to the
>>         convection problem at
>>         > low resolution is very diffusive and also why the interior
>>         temperature is much
>>         > closer to the surface temperature than to the CMB
>>         temperature because the
>>         > artificial viscosity is on the order of 20 times larger
>>         than the thermal
>>         > conductivity near the surface.
>>
>>         Would it be easy to verify whether the artificial viscosity
>>         ("artificial
>>         conductivity") decreases at the expected rate with mesh
>>         refinement?
>>
>>
>>     What is the most helpful way for me to show this? Visualization
>>     of a couple of slices from the 3D conduction model? I tried to
>>     get the depth average of artificial viscosity but the
>>     postprocessor is not implemented.
>>
>>         > For the conduction problem, the default values of the
>>         artificial viscosity are
>>         > also much larger than the thermal conductivity.
>>
>>         I think that's the point worth investigating. Since in this
>>         case the velocity
>>         is zero, one would expect the artificial viscosity to also be
>>         at least quite
>>         small. Why is it not?
>>
>>
>>     *Maybe the spherical and 2D annulus geometry models are returning
>>     an unhelpful length scale, like planetary radius instead of layer
>>     depth?*
>>
>>     aspect/source/simulator/entropy_viscosity.cc (starting line 191):
>>     // If the velocity is 0 we have to assume a sensible velocity to
>>     calculate
>>     // an artificial diffusion. We choose similar to nondimensional
>>     // formulations: v ~ thermal_diffusivity / length_scale, which
>>     cancels
>>         // the density and specific heat from the entropy
>>     formulation. It seems
>>         // surprising at first that only the conductivity remains,
>>     but remember
>>         // that this actually *is* an additional artificial diffusion.
>>         if (std::abs(global_u_infty) < 1e-50)
>>           return parameters.stabilization_beta *
>>                  max_conductivity / geometry_model->length_scale() *
>>                  cell_diameter;
>>
>>
>>         Best
>>           W.
>>
>>         -- 
>>         ------------------------------------------------------------------------
>>         Wolfgang Bangerth          email: bangerth at colostate.edu
>>         <mailto:bangerth at colostate.edu>
>>                                     www:
>>         http://www.math.colostate.edu/~bangerth/
>>         <http://www.math.colostate.edu/%7Ebangerth/>
>>
>>
>>
>>     _______________________________________________
>>     Aspect-devel mailing list
>>     Aspect-devel at geodynamics.org <mailto:Aspect-devel at geodynamics.org>
>>     http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
>
>     -- 
>     Rene Gassmoeller
>     https://gassmoeller.github.io/
>
>     _______________________________________________
>     Aspect-devel mailing list
>     Aspect-devel at geodynamics.org <mailto: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/20180829/dad798e3/attachment-0001.html>


More information about the Aspect-devel mailing list