An issue about phase transition.

[lhy11009]

Right now, our phase transition implementation could breaks when we have a phase transition in a material that starts from the surface of the model. We may only want one phase transitions in that material but the visual result will be we have two phase transitions. Here is a test case. We want at phase transition at 80 km depth here, but we ended up also having the phase transition when we go to the shallower region and the first phase only presents in between. Summary: Parameter file: And there is a clear problem in the implementation: See if we take (80km, 1373K) as our start point and 2.58e6 as our Claypeyron slope, we may cross the 'transition' at (0km, 273k), which is not what we want. I am not sure what would be the solution to that. I think the final solution will be to have a phase diagram. But for now, I think maybe to put a lower limit of depth to the transition. But in that case, we won't also fit the need of people where their transition is going bottom-top instead of top-bottom.

[lhy11009]

Two missing links: https://lhy11009.github.io/LearnAspect/phase_transition_issue/summary/ https://lhy11009.github.io/LearnAspect/phase_transition_issue/case.prm

[gassmoeller]

I think I do not quite understand all of your description, it seems you are describing two separate problems? Let me start with the second one: Yes, if one has a phase transition that is defined in a shallow depth and with a steep clapeyron slope it might intersect with the surface, because of the large temperature jump there. This is a problem with the way we (and most everyone else) defines these transitions. A linear transition in pressure-temperature space is simply not what happens in nature, the Clapeyron slope is technically only a linearized approximation at a certain pressure temperature condition. Additionally our approximation does not allow for phase transitions to disappear (e.g. because one of the two phases that border the transition become unstable). Having a depth or pressure cutoff seems reasonable conceptually, but may be a bit tricky to implement (although maybe a single new parameter for the PhaseFunction class may be enough?).

I do not completely understand the second part of your description. You have a model in which you set one phase transition, but you observe two phase transitions? I ran the model you attached, and it seems that the viscosity correctly shows only one jump. The density however has some weird behavior. Maybe the difference between density and viscosity already tells us something about what goes wrong? Have you checked which profile you expect for density and which one you get from the model? That may be another clue what is happening.

[lhy11009]

"it seems that the viscosity correctly shows only one jump."

For sure, because one phase's viscosity is dependent on temperature while another is not. I take from Magali's setting. I guess I just did too much copy and paste. I'll change that a little and let's see again.

[lhy11009]

I have altered the dependence of viscosity on temperature. Now it looks more clear. And it is very near to what happens in @mibillen 's model with Oceanic crust where she wants to apply a phase transition to eclogite but it turns out the transition also happens at a shallower depth, leaving the basaltic material only presented in the slab but not actually in the Oceanic seafloor. https://lhy11009.github.io/LearnAspect/phase_transition_issue1/summary/ https://lhy11009.github.io/LearnAspect/phase_transition_issue1/case.prm When you have free time, please take a look again. We are not urgent, but once we have discussed it, I'll go ahead and change the implementation. I believe other people who use the same implementation will be bewildered by this problem. You are right that the Claypayron slope is just an approximation of the real thing, although it is what my textbook taught me years before. But when I took a look at the actual PT of the basaltic oceanic sea, it is really complicated. Like this one here. I wonder if we could do the same thing in our implementation as for putting in this kind of phase diagram. And for now, let's discuss some quick amend to the problem itself.

[lhy11009]

@gassmoeller any update on this? Take your time, it would be great if you can reply me next week.

[gassmoeller]

Hey, got sidetracked. So here are two options of solving this:

1. Introduce two new parameters to the PhaseFunction class in material_model/utilities.cc named Phase transition minimum pressures, and Phase transition minimum depths and let the phase function and derivative always be zero below this pressure / depth. You will need to check if this fixes the problem (easiest if you use a simple example case like the tests we introduced when we wrote the phase function initially).
2. Alternatively you could look into using the functionality that @bobmyhill introduced in #3733. This would allow you to read in the whole complex phase diagram above and let the material model figure out what to do at different phases. This would be much more work though and would essentially replace the phase functions with a lookup for stable phases. Probably not what you want.

[lhy11009]

_Introduce two new parameters to the PhaseFunction class in material_model/utilities.cc named Phase transition minimum pressures, and Phase transition minimum depths_ Yes, that seems plausible. I am not going to do it right away though and I'll probably try to fix it right away, but probably very soon. Alternatively you could look into using the functionality that @bobmyhill introduced in #3733. I think I can manage it perhaps, but you are right that this is not what we want for now. We'll wait until needs pop-up.

[lhy11009]

Hi Rene How is it going lately? It has been a few months after our last message on this issue. I am distracted in the process of setting up and testing my model. Now I come to the point to make some changes to the phase function and I am ready to submit a pull request to fix this issue. The approach follows your first point: introduce two new parameters to the PhaseFunction class in material_model/utilities.cc named Phase transition minimum pressures, and Phase transition minimum depths Contact me back if you have the time to do a review. Thanks

[lhy11009]

Hi I found another issue in Material/utilities.cc. Now we use a "use_depth_instead_of_pressure" parameter to determine the way to define phase functions. This is from the "Define transition by depth instead of pressure" in the prm file. However, a default value is not give in the code so that this is always 1 (true). This is a little bit confusing. I could create a pull request to deal with this. Contact me back if you have the time to do a review. Thanks @gassmoeller

[lhy11009]

Hi How is it going lately? I think in terms of the phase transition minimum and maximum, we can also define a minimum temperature and a maximum temperature. Tell me what you think of it and when you have time to do a review? @gassmoeller

[gassmoeller]

Hi Haoyuan, replies to your individual points below:

Contact me back if you have the time to do a review. Thanks Feel free to open a PR about this. Someone will take a look as soon as possible.

However, a default value is not give in the code so that this is always 1 (true) This is not correct as far as I can see, the default value is given in this line: https://github.com/geodynamics/aspect/blob/d254afa1d96f8c63036cdc9a1c87be022a9f336e/source/material_model/utilities.cc#L1053 as "true". As long as the parameters of the phase function are correctly declared and parsed the parameter is set to true by default if it is not set to false in the input file.

I think in terms of the phase transition minimum and maximum, we can also define a minimum temperature and a maximum temperature. Tell me what you think of it and when you have time to do a review?

We could, but the important question is do you need it at the moment? I would take this one step at a time and only add functionality that is immediately useful. So if you need this feel free to add it at the same time as the pressure boundaries, but otherwise I would suggest waiting until it is necessary.

[lhy11009]

Yeah, that's good advice. So I have already done that in my own branch, I would submit a new pull request shortly after I finish this stage of model setting. Thanks for your reply. I'll tell you when I do that.

Best

On Mon, Jan 18, 2021 at 12:07 PM Rene Gassmoeller notifications@github.com wrote:

Hi Haoyuan, replies to your individual points below:

Contact me back if you have the time to do a review. Thanks Feel free to open a PR about this. Someone will take a look as soon as possible.

However, a default value is not give in the code so that this is always 1 (true) This is not correct as far as I can see, the default value is given in this line: https://github.com/geodynamics/aspect/blob/d254afa1d96f8c63036cdc9a1c87be022a9f336e/source/material_model/utilities.cc#L1053 as "true". As long as the parameters of the phase function are correctly declared and parsed the parameter is set to true by default if it is not set to false in the input file.

I think in terms of the phase transition minimum and maximum, we can also define a minimum temperature and a maximum temperature. Tell me what you think of it and when you have time to do a review?

We could, but the important question is do you need it at the moment? I would take this one step at a time and only add functionality that is immediately useful. So if you need this feel free to add it at the same time as the pressure boundaries, but otherwise I would suggest waiting until it is necessary.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/geodynamics/aspect/issues/3785#issuecomment-762448817, or unsubscribe https://github.com/notifications/unsubscribe-auth/AG2WMFOMVHYD5P3P3TYBYQLS2SIJJANCNFSM4P6WHK3Q .

--

Ph.D student of Geophysics Earth & Planetary Sciences Dept., UC Davis Davis, CA 95616

[lhy11009]

I believe this issue is addressed with PR #4591. A pair of temperature values are included in that PR for each of the phase transitions, where these values could be used to restrict the phase transitions to their respective range.