Skip to content

solidificationMeltingSource is buggy

Summary

In the constructor the alpha1 field calculating the liquid phase fraction (1 = full liquid) is initialized as zero.

    alpha1_
    (
        IOobject
        (
            IOobject::scopedName(name_, "alpha1"),
            mesh.time().timeName(),
            mesh,
            IOobject::READ_IF_PRESENT,
            IOobject::AUTO_WRITE,
            IOobject::REGISTER
        ),
        mesh,
        dimensionedScalar(dimless, Zero),
        fvPatchFieldBase::zeroGradientType()
    ),

The liquid phase fraction is calculated as follow from the temperature field :

alpha1New = alpha1_[celli] + relax_*Cpc*(Tc - Tmelt_)/L_;

For the first time step alpha1_[celli] is 0. Assuming we would like to model a pure solidification process (alpha1(t=0) = 1.0), the code will return alpha1(t=0)=alpha1New=relax_Cpc(Tc - Tmelt_)/L_ which can be < 1 whereas Tc > Tmelt. Ex:

| relax | 1 |

| Cpc | 900 |

| Tc | 1200 |

| Tmelt | 933 |

| L | 400000 |

| alpha1New | 0.60075 |

We end up by having an initial liquid phase fraction < 1.0 wich is obviously incorrect.

In case of melting the current implementation seems ok since alpha1(t=0) = 0.

Environment information

OpenFOAM version : v2306 Operating system : ubuntu Compiler : gcc

Possible fixes

I am not a specialist of melting/solidification so I am not sure of a correct fix. Maybe we could feed the constructor by an user input defining the initial fraction of liquid.

Edited by Paulin FERRO