Phase diagram option

[NicolasRiel]

Presently we cannot prescribe a density phase diagram using the Julia wrapper.

[boriskaus]

Do you mean these options?

rho_ph     = TestPD # name of the phase diagram you want to use (still needs rho to be defined for the initial guess of pressure)
rho_ph_dir =        # in case the phase diagram has a different path provide the path (without the name of the actual PD) here

[NicolasRiel]

yes!

[boriskaus]

this is actually there, so I am wondering why it doesn't work

[NicolasRiel]

Oh ok, maybe I missed it I will give it a try.

[NicolasRiel]

After checking, the rho_ph option works.

However the densities that are computed seem wrong. I doubled checked the values of phase_diagram file including the range and the units and I don't see any problem on that side (I checked against TestPD of the LaMEM main repository).
I also checked a phase diagram that works when running the C version of LaMEM and if I use the same diagram with LaMEM.jl the densities are wrong. The values are way too high, it seems like the melt fraction is not taken into account :/

Could there be a rounding/type issue?

[boriskaus]

The statement that they are „wrong“ is not so helpful. Can you perhaps upload a simple test case and upload a 1D plot of how it should be versus how it is? That would simplify the process of identifying the issue.

is the pressure perhaps not set to zero at the top of the domain?

[boriskaus]

Update: we looked into this and noted that it is important that this parameter:

 SolutionParams(mfmax=1.0)

is set.
The default in LaMEM is 0.15, which implies that only 15% of the melt is taken into account while computing effective density. That does seem to interact with the phase diagrams; I think it would be better to set this back to 1.0 by default in the main LaMEM as well.