Учет теплоты фазового перехода в задачах теплопроводности.

Здравствуйте, 

Очень вероятно, что придётся писать UDF. Обычный Expression не сможет задать вариацию источника для каждой отдельной ячейки. Там только для всего объёма можно.

Математика тепловых эффектов при промерзании есть? Решаете в стационаре или нестационаре?

Макрос DEFINE_SOURCE для уравнения энергии

Пример использования макроса DEFINE_SOURCE:

Данный код используется для задания источника тепла при образовании тумана, а точнее при диффузионном росте капли тумана или центра конденсации в зоне пересыщения.

C_UDMI(c, t, UDM_RH) - это записанная в память влажность воздуха в рамках макроса DEFINE_ADJUST.

real molefrac(real mass_frac, real M1, real M2)
{
    real x;  /* x is mole fraction */
    if(mass_frac <= 0.0)
        x = 0.0;
    else
        x = M2/( M2+M1*( (1./mass_frac) - 1.) );
    return x;
}

real H2O_Saturation_Pressure(real tsat)
{

real psat;

if (tsat >= 273.16)
    psat =100000*pow(10, (28.59051 - 8.2*log10(tsat) + 0.0024804*(tsat) - 3142.31/(tsat)));
else
    psat =100*pow(10, (10.5380997 - 2663.91/(tsat)));
 
return psat;
}

DEFINE_SOURCE(my_energy_source,c,t,dS,eqn)
{
    real source;    /* vapor source */
    real sigma = 1.0;    /* evaporation efficiency */
    real constant = sigma *
        sqrt( mw_vapor / (2.0 * M_PI* UNIVERSAL_GAS_CONSTANT * C_T(c,t)));    /* condensation relaxation time, 1/s */
    real diameter_drop = 2.0e-5; /* 10-50 mkm */
    real latent_heat = 2258000.0; /* latent_heat, J/kg */

    
    real X_vap = C_YI(c, t, 0);  /* mass fraction of vapor */
    real Y_vap = molefrac(X_vap, mw_vapor, mw_air);  /* mole fraction, alpha vapor */
    real X_drop = C_YI(c, t, 1);  /* mass fraction of droplets */

    real area_density = 6.0/diameter_drop*X_drop*C_R(c,t)/1000.0; 

    real P_sat = H2O_Saturation_Pressure(C_T(c,t));

    if (C_UDMI(c, t, UDM_RH) > 101.0)
    {
    
    if (X_drop > 0.0)
    {    source = -latent_heat*area_density*Y_vap*constant*P_sat*((100.0-C_UDMI(c, t, UDM_RH))/100.0);
    dS[eqn] = 0.0;    }
    
    
    else
    {source = 1.0e-7*latent_heat;
    dS[eqn] = 0.0;}

    }

    else 
  {
    source = 0.0;
    dS[eqn] = 0.0;
  }
      
      C_UDMI(c, t, UDM_RH_new) = source;
    return source;
}

По-моему, это безумная идея)))

Я бы на вашем месте задал бы землю как жидкость с очень высокой вязкостью и постоянной плотностью, чтобы она не двигалась при смене её температуры. 
Далее нужно задать ей теплопроводность как у грунта и уменьшенную скрытую теплоту образования с учётом того, что доля воды там не 100%. 

В итоге ко всему этому можно применить модель Solidification/melting, чтобы смоделировать то, как всё это будет послойно замерзать.

Способ очень прост и не требует программирования UDF.