_RL const1,const2
const1=HeatCapacity_Cp*rhoNil*drF(1)
& (deconcbar(i,j,k,bi,bj)-SEAICE_freeze)*const1
& ( (deconcbar(i,j,k,bi,bj)-SEAICE_freeze)*const1
& const2 ) - (epsilonTemp*const1)