Continental geotherm should support new syntax for Densities from visco plastic material model

[mitchellmcm27]

The visco plastic material model now supports multiple phases using a key-value syntax, e.g., set Densities = background:3000, lithosphere:3000|3300. The continental geotherm plugin reads the densities from the visco plastic model, but does not support this new syntax. It throws from dealii with "Can't convert <background:3000> to a double."

A quick github search for prm.get("Densities") didn't reveal any other code that depends on the densities from the visco plastic model, but I could've missed it.

[gassmoeller]

Hi @mitchellmcm27, thanks for letting us know this is indeed not supported at the moment and we should provide a clearer message about that. Just to clarify: if the densities for the visco plastic model are provided using the old syntax (a list of numbers) the continental geotherm plugin still works as expected, is this correct?

In order to support the new syntax the continental geotherm plugin would need to read the density as is now done here: https://github.com/geodynamics/aspect/blob/5fee496288104da75bc9443f8aa0de523ec0708b/source/material_model/equation_of_state/multicomponent_incompressible.cc#L104 and would need to support multiple phases (i.e. computing phase function values for the current conditions and computing the effective density at this depth. Not impossible, but not as straightforward as the current implementation.

Is this an important feature for your model? Maybe @naliboff or @anne-glerum have a better idea about how to implement this?

[mitchellmcm27]

Thanks for the quick response @gassmoeller. AFAIK the plugin works correctly using the old syntax. I didn't test it rigorously, but my model runs for several time steps after I make the change.

I'd thought about doing a PR but realized I didn't know the best way to handle the phase calculations. I just wanted to make sure the team was aware of the incompatibility.

From my perspective, one solution would be for the continental geotherm plugin to use its own Densities parameter. The user is already duplicating information on layer thicknesses in the prm file to use the continental geotherm, so it wouldn't add much mental overhead to have a bit of duplication in density values. For now, I can work around the incompatibility by patching the plugin to do this on my end.

[naliboff]

@mitchellmcm27 - Thanks for catching this and posting the issue. The plugin should be working correctly with the old syntax, as there is a currently a test that uses the continental geotherm plugin.

However, that still leaves the question about how to deal with cases where there are phase changes. The densities in the plugin are used to calculate thermal conductivity, which also requires the thermal diffusivity and heat capacities.

Long term, we actually decided to make the switch to directly specifying thermal conductivity instead of thermal diffusivities in Visco Plastic. Now that 2.2 has been released I think we can go ahead with that, which would also require updating the continental geotherm plugin.

Does this seem like a reasonable long-term solution?

[mitchellmcm27]

@naliboff that sounds like a great solution. I was aware of the switch to thermal conductivity, but didn't connect the dots to the continental geotherm plugin. For now, I was able to successfully patch the plugin to use its own Densities parameter. Thanks!

[naliboff]

@mitchellmcm27 - Great, we'll go ahead with that solution unless a better idea surfaces.

For now, I was able to successfully patch the plugin to use its own Densities parameter. Thanks!

Just to be sure, the patch is only required if one has multiple phases? If the patch is required for the case with no phases (i.e., prior syntax) then we should a fix.

[mitchellmcm27]

Correct, only required for multiple phases using the new key:value syntax.