Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Be safer and more correct when using strain rate #5239

Merged
merged 1 commit into from Jul 15, 2023

Conversation

gassmoeller
Copy link
Member

This PR contains three changes that make using the strain rate in material models safer. All of these changes were necessary for #5211.

  1. The Stokes assembly currently only provides the strain rate to the material model if the viscosity is requested. However, with complex material models we dont know what other material properties the strain rate may be needed for. For example the visco_plastic and viscoelastic material models also need the viscosity to compute the elastic force terms.
  2. The visco-plastic rheology sometimes has to use a reference strain rate, but the condition for when to use the reference strain rate seems wrong. It needs to use the reference strain rate during the initialization (when timestep_number is undefined), during the first timestep (but only the first nonlinear iteration when we dont have the Stokes solution yet), and if the strain rate is very small (because we need to be able to divide through it).
  3. The visco-plastic material model only computes the viscosity when it is requested, but it also uses the viscosity for the additional outputs (the elastic force terms). Make sure it computes the viscosity if the additional_outputs are requested.

Copy link
Member

@bobmyhill bobmyhill left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks sensible. One grammatical suggestion.

@@ -218,9 +218,11 @@ namespace aspect
out.entropy_derivative_temperature[i] = MaterialUtilities::average_value (volume_fractions, eos_outputs.entropy_derivative_temperature, MaterialUtilities::arithmetic);

// Compute the effective viscosity if requested and retrieve whether the material is plastically yielding
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
// Compute the effective viscosity if requested and retrieve whether the material is plastically yielding
// Compute the effective viscosity if requested and retrieve whether the material is plastically yielding.

@gassmoeller
Copy link
Member Author

I think this broke some tests, but I am unsure if the new version is better or worse. @naliboff and @MFraters I think this affects your models most. Could you take a look at this PR and tell me what you think about the changed test results?

@bangerth
Copy link
Contributor

Specifically, the single failing test is:

	823 - solubility (Failed)

Isn't that test new? Did that test succeed on your branch before the commits?

Comment on lines -223 to +228
if (in.requests_property(MaterialProperties::viscosity))
if (in.requests_property(MaterialProperties::viscosity) || in.requests_property(MaterialProperties::additional_outputs))
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This seems right, but points out the fact that we should add an else to this if that poisons the output fields then not computed with signaling_nans.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@bangerth - are the material model outputs by default not already set to signaling nans?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Perhaps -- we weren't sure whether it's possible that the same object was reused from a state where information was requested to where it is no longer requested (but at that point no longer contains the signaling nans).

need_viscosity);
/*compute_strain_rate=*/ true);
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Doesn't this negate the point of the flag? Do you even still need need_viscosity in this case?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We still need the flag, because need_viscosity is used in the next line to tell the material model which outputs to compute. The wrong assumption here was that if we do not explicitly request the viscosity we also never need the strain rate input. This is not true anymore for complex material model (strain rate may be used for other properties than just viscosity).

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ah, that's a bummer. I recall that computing the strain rate was expensive, and it's a shame that we now compute it unconditionally because the material model does not have a way to indicate whether or not it actually wants it...

@bangerth
Copy link
Contributor

Since you observe that this patch changes output, the conclusion is that apparently we forgot to compute something previously and used invalid values (perhaps zeros, perhaps computed previously on the same object). This suggests what I mention in a comment above: We should have poisoned these values in an else branch. It would be interesting to write such an else branch for a separate PR and see what fails.

@bangerth
Copy link
Contributor

See #5247 for an attempt.

Copy link
Contributor

@naliboff naliboff left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@gassmoeller - Aside from addressing comments from @bangerth, this looks good to go. Thanks for adding this.

One general comment - at some point it may be worth adding the checks for using the reference strain rate to the material model utilities, as other material models have to do the same checks.

@bangerth bangerth merged commit 38b3201 into geodynamics:main Jul 15, 2023
6 checks passed
@gassmoeller gassmoeller deleted the fix_some_uses_of_strain_rate branch July 15, 2023 03:26
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

5 participants