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

fix latent heat viscosity prefactor #6214

Open
wants to merge 2 commits into
base: main
Choose a base branch
from

Conversation

jdannberg
Copy link
Contributor

The list of viscosity prefactors that would change the viscosity at phase transitions in the latent heat material model did not work correctly if more than one compositional field was used. Specifically, individual layers would have the wrong viscosities, and the length of the list of prefactors that had to be specified in the input file made no sense. This specific setup had no test case, which is likely why we did not notice when we broke this.
The problem was that the prefactors are specified relative to each other, but when they are applied, they are multiplied on top of each other (due to the way the phase function works). This means they have to be scaled after they are read in, and this scaling depends on how many phase transitions are occurring above a given phase, and the old interface just had no good way of figuring that out. Therefore, there was no good way for fixing this without also changing the interface for this input parameter. I changed the interface to be consistent with what we use for all other phase transition parameters (the map parsing from the material model utilities), which has the advantage that there's a good default for models with and without chemical heterogeneities. I also fixed the computation of the viscosities in the individual phases, and changed the averaging between phases to a log average (which is consistent with the rest of the material model, which already uses log averaging for viscosities).

For new features/models or changes of existing features:

  • I have tested my new feature locally to ensure it is correct.
  • I have created a testcase for the new feature/benchmark in the tests/ directory.
  • I have added a changelog entry in the doc/modules/changes directory that will inform other users of my change.

@@ -131,18 +136,27 @@ namespace aspect
// for the first layer, so we have to use phase+1 as the index.
if (composition.size()==0)
{
phase_dependence += phaseFunction * density_jumps[phase];
viscosity_phase_dependence *= 1. + phaseFunction * (phase_prefactors[phase+1]-1.);
phase_dependence += phaseFunction * density_jumps[transition_index];
Copy link
Member

Choose a reason for hiding this comment

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

can you fix/delete the comment 3 lines above?


// Viscosities
if (transition_index+1 < phase_function.n_phases_for_each_composition()[0]) // 1st compositional field
viscosity_phase_prefactors[0] = std::pow(10,
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
viscosity_phase_prefactors[0] = std::pow(10,
viscosity_phase_prefactors[0] = std::pow(10.,

+ (1. - phaseFunction) * std::log10(viscosity_phase_prefactors[0]));
else if (transition_index+2 < phase_function.n_phases_for_each_composition()[0]
+ phase_function.n_phases_for_each_composition()[1]) // 2nd compositional field
viscosity_phase_prefactors[1] = std::pow(10,
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
viscosity_phase_prefactors[1] = std::pow(10,
viscosity_phase_prefactors[1] = std::pow(10.,

if (composition.size()==0)
out.viscosities[i] = out.viscosities[i]*viscosity_phase_prefactors[0];
else
out.viscosities[i] = (std::pow(10, ((1-composition[0]) * std::log10(out.viscosities[i]*viscosity_phase_prefactors[0])
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
out.viscosities[i] = (std::pow(10, ((1-composition[0]) * std::log10(out.viscosities[i]*viscosity_phase_prefactors[0])
out.viscosities[i] = (std::pow(10., ((1-composition[0]) * std::log10(out.viscosities[i]*viscosity_phase_prefactors[0])

"viscosity (modified by any compositional prefactors) will be multiplied "
"by this factor to get the corresponding viscosity for each phase. "
"List must have the same number of entries as there are phases, that is one "
"more than Phase transition depths for each composition that is used in the "
Copy link
Member

Choose a reason for hiding this comment

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

do you want to quote this parameter?

@jdannberg jdannberg force-pushed the latent_heat_viscosity branch from 1299048 to 38227d4 Compare January 29, 2025 13:26
@jdannberg jdannberg force-pushed the latent_heat_viscosity branch from 38227d4 to 2b2b703 Compare January 29, 2025 13:28
@jdannberg
Copy link
Contributor Author

@tjhei Thank you for the review!
I've addressed your comments and updated the one test that was failing. I checked that the results are still correct (the small change likely comes from fixing the averaging across phase transitions to use a log average rather than an arithmetic average).

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants