[aspect-devel] Thermodynamic consistency of Aspect's temperature and momentum equations

Wed Feb 13 10:18:39 PST 2013

Hi,

  I haven't worked much with Aspect but I do know about the governing
equation approximations.

  The viscous dissipation and adiabatic heating terms will only
perfectly cancel each out in the ALA governing equation. The TALA
equations cause the viscous dissipation term to be a little bigger(?)
than the adiabatic term. This is explained in the Leng and Zhong, 2008
paper I've attached.

  In CitcomS, the incompressible Boussinesq has no viscous dissipation
and adiabatic heating terms. The extended Boussinesq approximation has
those two terms, with an adiabat but is still incompressible. TALA is
compressible but assumes something about the pressure (can't
remember), which separates it from ALA. ALA and the full pressure term
is too hard to solve for numerically so we tend to stick with those
other approximations.

 - Katrina Arredondo

On Wed, Feb 13, 2013 at 4:27 AM, Rene Gassmoeller <rengas at gfz-potsdam.de> wrote:
> Hi all,
> Juliane and myself had a quite lengthy discussion about this issue (it
> started yesterday and it is actually not finished yet ;-)) but I just wanted
> you to know what would be our opinion up to now (to save double effort in
> thinking).
> We think Ian had an interesting point here, although Thomas is right about
> the correctness of the currently used term in a way as well. Ian's term is
> derived from the momentum equation and exactly balances out the viscous
> dissipation (since it is derived that way). However, there is no guarantee
> that this term actually is(!) the adiabatic heating in all material models /
> approximations, except that it seems to be the right term for the boussinesq
> / anelastic liquid approximation with a simple material model.
> The term that is currently used in Aspect on the other hand is derived from
> the actual thermodynamic equations (dS/dt = ... the root of our equation of
> conservation of energy as stated in "Mantle convection in the Earth and
> Planets") and is therefore consistent, but it uses at least one
> approximation, namely that the pressure is only lithostatic and there is no
> dynamic pressure. The right term (in a thermodynamic sense) would not be:
>
>    Q_a = ( velocity * gravity ) * alpha * density * temperature
> (1)
>
> but rather:
>    Q_a = ( velocity * pressure_gradient ) * alpha * temperature
> (2)
>
> The gravity and density in the currently used term (1) just came out of the
> assumption that p = rho * g * h, and by the way, this justifies the use of
> the real rho instead of the reference_rho in the equation (1) although it
> does not seem to fit with the viscous heating. The use of the above
> mentioned formula (2) reduces the error in a test case (attached) from
> around 3% to 0.2% ... that is not perfect considering the fit of Ian's term
> of 1e-13 or so. Still I would currently rather use term (2) since it also
> makes us independent of the choice of a right reference density and it is
> also consistent in the compressible case (because it makes no assumptions on
> the material model at all, it just needs the real alpha und the local
> solution variables). Additionally with this term we also capture the heating
> of material that moves laterally along a pressure gradient and therefore
> also cools/heats up adiabatically.
> Nevertheless, the misfit of 0.2% seems too large to be satisfied. I first
> thought that the inconsistency between the definition of alpha = - (1 / rho)
> * drho / dT = const and the definition of rho = rho_ref (1 - alpha * (T -
> T_ref)) leads to the misfit and tried a definition of rho = rho_ref * exp (-
> alpha * (T - T_ref)), which is the solution of the definition of alpha
> above. But this does not seem to make a large difference in the error. So
> currently we are thinking in two ways:
> 1. Find the reason for the misfit between thermodynamically derived
> adiabatic heating term and viscous dissipation that is computed out of the
> stokes velocity field.
> 2. If we can derive that Ian's term actually is the adiabatic heating term
> in a general sense (perhaps with an additional or without an currently used
> approximation in term (2)), then of course that term would be more precise,
> although then we end up asking how to use this in a compressible model.
>
> Currently, I am a bit sceptical about the possibility of point 2. therefore
> I am thinking of point 1.. The next part is speculation but I find it a bit
> intriguing that Leng and Zhong (2008) also end up with a difference of the
> order of 0.2 % although they do not consider the formula (2) but rather
> Ian's term extended by an compressible correction. Could it be that the
> remaining 0.2% originate from an approximation that is not bound to the
> formulation of the adiabatic heating term but rather something influencing
> the velocity field (and therefore the viscous dissipation)?
>
>
> By the way during this discussion we came across the formulation of
> compressible stokes flow in aspect and think there is an inconsistency as
> well. In the right-hand side of the stokes equation we have the term rho*g
> as well as an additional compressible term. Considering a compressible term
> in the compressible case seems reasonable if the first term just covers the
> temperature dependency of rho, but since in Aspect the material model takes
> care of all the influences on the density, this term should be moved to any
> material model that wants to be capable of doing compressible convection.
> For example a material model could take the density values directly from a
> thermodynamic calculation rather than a simplified equation of state. In
> that case the current formulation would double-count the pressure effect
> (without any chance for the user to know this, since he thinks the material
> model need to take care of all the influences).
>
> We are continuing working on this, but please let us know, if you have an
> opinion or suggestion on this (and also in case something in this is wrong
> ;-)),
>
> Cheers,
> Rene
>
>
> --
> Rene Gassmoeller
>
> 2.5/Geodynamic Modelling
> Tel.: +49 (0)331/288-28744
> Fax:  +49 (0)331/288-1938
> Email: rengas at gfz-potsdam.de
> ___________________________________
>
> Helmholtz Centre Potsdam
> German Research Centre for Geosciences GFZ
> Telegrafenberg, 14473 Potsdam
>
>
>
> On 02/13/2013 01:21 AM, Magali Billen wrote:
>
> Hello All,
>
> Each of the different approximations of the equations have specific sets of
> terms that drop in or out together,
> so you need to be careful adding back in just one term or another without
> taking into account which approximation
> has leads to that term being assumed small. When people state Boussinesq
> approximation in mantle convection
> calculations, this has historically meant an approximation that does not
> include adiabatic heating. Also, I agree with
> Thomas that it doesn't really make sense to have such a term in an
> incompressible convection case, since without compression there is no
> physical cause for an adiabatic gradient.
>
> Speaking for myself, I don't have the specific approximations memorized,
> however, I've found Chapter 6 of
> "Mantle Convection in the Earth and Planets" by Schubert, Turcotte and
> Olson, very helpful when trying to understand
> the origin and loss of terms for boussinesq, extended-boussinesq and TALA
> approximations - it goes through
> each of these approximation in detail and explains the specific assumptions
> for each approximation and which
> terms drop out. I worked through all these once about a year ago and I think
> it might help with this particular question.
>
> Magali
>
> On Feb 12, 2013, at 4:09 PM, Ian Rose wrote:
>
> Hmm, I am not sure I agree.  Di is frequently assumed to be zero in mantle
> convection problems, but that is not a result of the Boussinesq
> approximation.  That is to say, the "work done by gravity" term in the
> kinetic energy equation arises just fine with Boussinesq (there are just
> more terms that come from the div velocity terms in the compressible case).
>
> Even though this term (along with viscous dissipation) are likely to be
> smallish, I see no reason not to allow them to be turned on and off with
> flags as they are now.  But if it is turned on, it should be consistent with
> what you get from integrating the momentum equation.
>
> Cheers,
> Ian
>
>
>
> On Tue, Feb 12, 2013 at 3:28 AM, Thomas Geenen <geenen at gmail.com> wrote:
>>
>> he Timo,
>>
>> there is no such thing as adiabatic heating in the incompressible
>> Boussinesq case Di (alpha*g/cp) is assumed zero .
>> for extended Boussinesq there should also be no problem since there is no
>> density in the net adiabatic heating term.
>>
>> setting  thermal diffusion, viscous dissipation and internal heating to
>> zero (dS/dt=0) we end up with
>> rhocp(dT/dt) - alphaTdP/dt=0
>> or
>> rho*cp*(dT/dt) - alpha*rho*g*u_r*T=0
>>
>> this will give for an adiabatic temperature profile
>> T(r) = T_0*exp(alpha*g*r/cp)
>>
>> iow the density does not play a role since its devided out of the
>> equation.
>>
>> this also holds for the compressible case i would say.
>>
>> cheers
>> Thomas
>>
>>
>>
>> On Tue, Feb 12, 2013 at 5:57 AM, Timo Heister <heister at math.tamu.edu>
>> wrote:
>>>
>>> Hey everyone,
>>>
>>> Ian approached me about this and I asked him to write it down here.
>>> Does anyone have any feedback about this, especially (assuming this is
>>> correct), what to do in the compressible case?
>>>
>>> On Wed, Feb 6, 2013 at 6:33 PM, Ian Rose <ian.rose at berkeley.edu> wrote:
>>> > Hi Aspect folks,
>>> >
>>> > I was working through some tests with Aspect and came across what I
>>> > believe
>>> > is an inconsistency in the governing equations.
>>> >
>>> > For incompressible Boussinesq flow, the global viscous dissipation
>>> > should
>>> > exactly cancel the global adiabatic heating.  This can be seen by
>>> > multiplying the momentum equation by velocity and integrating over the
>>> > domain.
>>> >
>>> > As it stands in assembly.cc, the formula used for calculating adiabatic
>>> > heating is different from that you would get by integrating the
>>> > momentum
>>> > equation.  I wrote a simple postprocessor that compares the two
>>> > integrated
>>> > quantities which I am attaching.  The difference is quite a lot for the
>>> > current formula.
>>> >
>>> > Put another way, this is the formula that is currently used:
>>> >
>>> >    Q_a = ( velocity * gravity ) * alpha * density * temperature
>>> >
>>> > The density at this point however, has already been adjusted for
>>> > temperature, so we are in effect double counting the thermal expansion.
>>> > Instead, I believe it should be
>>> >
>>> >   Q_a = ( velocity * gravity ) * ( density - reference_density )
>>> >
>>> >
>>> > The compressible case, too, should require some thought, though I have
>>> > not
>>> > gone through the paces there.
>>> >
>>> > Thoughts?
>>> >
>>> > Best,
>>> > Ian
>>> >
>>> > PS, for some details on the derivations, I refer you to Leng and Zhong
>>> > (2008)
>>> >
>>> >
>>> > _______________________________________________
>>> > Aspect-devel mailing list
>>> > Aspect-devel at geodynamics.org
>>> > http://geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
>>>
>>>
>>>
>>> --
>>> Timo Heister
>>> http://www.math.tamu.edu/~heister/
>>> _______________________________________________
>>> Aspect-devel mailing list
>>> Aspect-devel at geodynamics.org
>>> http://geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
>>
>>
>>
>> _______________________________________________
>> Aspect-devel mailing list
>> Aspect-devel at geodynamics.org
>> http://geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
>
>
> _______________________________________________
> Aspect-devel mailing list
> Aspect-devel at geodynamics.org
> http://geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
>
>
> --------------------------------------------------
> Associate Professor & Chancellor Fellow, U.C. Davis
> Chair, Geology Graduate Program
> Department of Geology & KeckCAVES
> 2129 Earth & Physical Sciences Bldg
> Davis, CA 95616
>
> E-mail: mibillen at ucdavis.edu
> Phone: no extension, please e-mail
> --------------------------------------------------
>
>
>
>
>
>
>
>
>
> _______________________________________________
> Aspect-devel mailing list
> Aspect-devel at geodynamics.org
> http://geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
>
>
>
> _______________________________________________
> Aspect-devel mailing list
> Aspect-devel at geodynamics.org
> http://geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
-------------- next part --------------
A non-text attachment was scrubbed...
Name: Leng_ZhongGJI2008-1.pdf
Type: application/pdf
Size: 161020 bytes
Desc: not available
Url : http://geodynamics.org/pipermail/aspect-devel/attachments/20130213/c2b70256/attachment-0001.pdf 