CFD Online Discussion Forums - Melting of Steel via 3D Goldak Heat Source

Hello together,

I'm trying to model a 3D volumetric Goldak heat source that is moving over a simple steel plate.
The problem I’m facing is that the simulation stops after a short amount of iterations (about 0.2s simulated time). I used an expression for thermal expansion and phase density difference buoyancy forces, the Kozeny-carman permeability model with min/max clipping. The setup is shown in the attached pictures. The used CEL expressions are :

LIBRARY:
CEL:
EXPRESSIONS:
BetaLiq = 4.2333e-05[K^-1]
BuoySrcRho = ((LiqFrac*rhoL+(1-LiqFrac)*rhoS)-rhoRef)*g
BuoySrcT = -BetaLiq*rhoL*(T-Tref)*LiqFrac*g
Kperm = Kstar*(1-LiqFrac)^2/max(LiqFrac^3, 1e-10)
Kperm2 = max(min(KpermMax,Kperm), 1e-15 [m^-2])
KpermMax = 1.0E+9[kg s^-1 m^-3] / Visc
Kstar = 1e10[m^-2]
L = 246800[J/kg]
LiqFrac = Liquid Steel.Mass Fraction
PL = 2400[W]
Q = if(qf1*qf+qr1*qr<=0[kg m^-1 s^-3],1[kg m^-1 s^-3],qf1*qf+qr1*qr)
Tcold = 1700[K]
Thot = 1800[K]
Tmelt = 1800[K]
Tref = Tcold
Visc = Liquid Steel.Dynamic Viscosity
a = 0.004[m]
b = 0.005[m]
cf = 0.0025[m]
cr = 0.0075[m]
ff = 0.5
fr = 1.5
kf = (6*sqrt(3)*ff*PL)/(a*b*cf*pi*sqrt(pi))
kr = (6*sqrt(3)*fr*PL)/(a*b*cr*pi*sqrt(pi))
qf = kf*exp(-3*((x-vl*t)^2/cf^2+y^2/a^2+z^2/b^2))
qf1 = if(x>vl*t,1,0)
qr = kr*exp(-3*((x-vl*t)^2/cr^2+y^2/a^2+z^2/b^2))
qr1 = if(x<vl*t,1,0)
rhoL = Liquid Steel.Density
rhoRef = rhoS
rhoS = 7854[kg/m^3]
vl = 0.005[m s^-1]
END
END

Something is wrong with this model, but I don’t know what it is.
Maybe someone has a clue what is wrong with the simulation.

Best regards
Joe T.