Question about units and usage of global variable `Use years in output instead of seconds`

ASPECT
1

Dear all,

I am a beginner ASPECT user and I have a question about how units are handled across plugins. I would really appreciate any guidance people can provide. Specifically, I would like to understand the intended usage of the global variable called Use years in output instead of seconds. I understand from the documentation that this parameter does not just affect output, but also how inputs to certain plugins are interpreted.

I have run the cookbook called “Kinematically-driven 2d oceanic subduction” two different ways:

  1. With the global variable Use years in output instead of seconds set to true, and velocity boundary conditions set using units of m/yr.
  2. With the global variable Use years in output instead of seconds set to false, and velocity boundary conditions set using units of m/s.

The resulting velocity fields output by ASPECT are a bit different from each other, which I did not necessarily expect. Also, the output for both cases are on the same order of magnitude (m/yr), when I would have expected the output from the second case to be in m/s.

I wonder if anyone can help me understand:

  • Is there some other part of this cookbook (a different plugin) that is being affected by this flag? For example, should the units of Viscosity in the Multicomponent Material Model be converted also?
  • Should I expect the output from the two cases I specify above to be the same?
  • If I set boundary conditions in m/s, and the flag is set to false, why is the output on the order I’d expect from m/yr?
  • How would you recommended using this global variable when setting up models?

I have attached the two parameter files here, and the magnitude of velocity in each case is plotted below.

Thank you very much in advance for any help!

Best regards,
Gabrielle Hobson

velocity_magnitude_flag_true
velocity_magnitude_flag_true2124×456 75.2 KB

velocity_magnitude_flag_false
velocity_magnitude_flag_false2124×456 75 KB

kinematically_driven_subduction_2d_case1_flag_true.prm (5.8 KB)
kinematically_driven_subduction_2d_case1_flag_false.prm (5.8 KB)

2

@ghobson - Welcome and thank you for posting your questions to the forum!

I understand from the documentation that this parameter does not just affect output, but also how inputs to certain plugins are interpreted.

Indeed, this parameter affects both outputs and inputs. For others reading and unfamiliar with this parameter, here is a link to the documentation.

The resulting velocity fields output by ASPECT are a bit different from each other, which I did not necessarily expect. Also, the output for both cases are on the same order of magnitude (m/yr), when I would have expected the output from the second case to be in m/s.

Yes, I would expect them to be the same as well. After doing a diff on the parameter file I have a few ideas on what the underlying causes might be. In detail, the only differences are whether Use years in output instead of seconds is set to true or false and the user defined constant used for the prescribed velocity:

<     set Function constants  = cm=100.0 #*3600.0*24.0*365.25
---
>     set Function constants  = cm=100.0*3600.0*24.0*365.25

The origin of the slightly different (but same order of magnitude velocities) may be related to how the number of seconds per year is calculated in above versus in ASPECT.

The parameter years_in_seconds in ASPECT (see here) is defined as 60*60*24*365.2425, which gives a slightly different value (31556952.0) than the equation used in the function constant definition (31557600.0).

That may be enough to explain the slight difference in velocity values?

For the case where Use years in output instead of seconds = false, the units for the output velocity will indeed by m/s, which should give values lower than the other case by a factor of 31556952.0.

I think the reason the magnitudes stay approximately the same is that you modified your boundary velocity to be multiplied by 3600.0*24.0*365.25.

This effectively gets you back to units of m/yr for the boundary velocity, but again is slightly different than ASPECT’s internal seconds to years conversion.

So, if you just use set Function constants = cm=100.0 in both PRM files I think you should see velocities fields that are different from each other by exactly a factor of 31556952.0?

I hope this helps!

Cheers,
John

3

Dear John,

Thank you very much for your response! I appreciate you taking the time to do so.

That’s good to know regarding the slightly different definition of the parameter years_in_seconds. I have modified the unit conversion in my parameter file to match accordingly.

Your details about the factor being used for conversion raised a good point and prompted me to test a bit more, and I now understand what was causing my issue. It turns out I was expecting function constants to be computed/evaluated somehow before being passed into the expression (I suppose I’m too used to Python) and therefore I was missing many factors in my unit conversion.

In detail: I was using the function constant defined in the Boundary velocity model function for unit conversion, rather than as the velocity directly. For example, in the case where Use years in output instead of seconds = true, the function looks like so:

subsection Boundary velocity model
  set Tangential velocity boundary indicators = left, bottom, top
  set Prescribed velocity boundary indicators = right x:function

  subsection Function
    set Function constants  = cm=100.0 # convert cm to yr
    set Function expression = if(z<540000.0, (600.0/550.0)/cm, \
                                 if(z>560000.0, -5.0/cm, \
                                    ((((-600.0/550.0)-5.0)/-20.0)*((z/1000.0)-560.0)+5.0)*(-1.0/cm))); \
                              0
    set Variable names      = x,z
  end
end

And so cm=100.0 is just being used to divide the velocity to convert it from cm/yr to m/yr. So I would expect the velocity output to be (on the right half of the domain for example), 0.05 (m/yr). And this is what I see in the output, so that made sense.

In the case where Use years in output instead of seconds = false, if I write the function like so:

subsection Boundary velocity model
  set Tangential velocity boundary indicators = left, bottom, top
  set Prescribed velocity boundary indicators = right x:function

  subsection Function
    set Function constants  = cm=100.0*60*60*24*365.2425
    set Function expression = if(z<540000.0, (600.0/550.0)/cm, \
                                 if(z>560000.0, -5.0/cm, \
                                    ((((-600.0/550.0)-5.0)/-20.0)*((z/1000.0)-560.0)+5.0)*(-1.0/cm))); \
                              0
    set Variable names      = x,z
  end
end

It turns out this is not equivalent to defining the function constant directly, like this:

subsection Boundary velocity model
  set Tangential velocity boundary indicators = left, bottom, top
  set Prescribed velocity boundary indicators = right x:function

  subsection Function
    set Function constants  = cm=3155695200.0
    set Function expression = if(z<540000.0, (600.0/550.0)/cm, \
                                 if(z>560000.0, -5.0/cm, \
                                    ((((-600.0/550.0)-5.0)/-20.0)*((z/1000.0)-560.0)+5.0)*(-1.0/cm))); \
                              0
    set Variable names      = x,z
  end
end

Because in the former case, if cm=100.0*60*60*24*365.2425 and I divide -5.0/cm, it seems like it is actually only interpreting the first part and using cm=100.0, rather than passing it in and doing -5.0/(100.0*60*60*24*365.2425) which is what I had expected/intended it to do. This explains why I was getting the same output when Use years in output instead of seconds = false whether I used cm=100.0*60*60*24*365.2425 or cm=100.0.

I’m now getting the expected output in both cases (with Use years in output instead of seconds = true and Use years in output instead of seconds = false), when I correctly define the constant.

Long story short, I’ve learned that function constants must be exactly that, constants, and these constants don’t necessarily get evaluated/computed before being passed into the expressions. It seems obvious when I write it down like that - even though it’s just a misunderstanding on my part, I figured I’d write down this conclusion here in case other beginners find it useful to know in future.

Thanks again for your help, John.

Best regards,
Gabrielle Hobson

4

@ghobson This seems like the right conclusion. I don’t think anyone ever thought of defining a constant as an expression. I looked it up in the code, where the Function constants parameter is parsed, and would have hoped that we produce an error in a case like yours. It’s not clear to me why not, but I opened a report about it: Ensure 'Function constant' in ParsedFunction is a constant, not an expression. · Issue #17535 · dealii/dealii · GitHub

Either way, I’m glad you figured it out!

Best
W.

5

@ghobson - Indeed, agree this is the right conclusion and glad my response helped guide the next round of testing, even if my initial guess was off :slight_smile:

I think this is issue of the equations being used in function constant definitions has come up before, and given that an error message is not thrown it is to miss.

Thank you again for posting to the forum as I’m sure this discussion has been useful for a quite a few folks!

Cheers,
John