## phase change comparison elmer/Comsol

Numerical methods and mathematical models of Elmer
julien givernaud
Posts: 71
Joined: 18 Nov 2014, 18:18
Antispam: Yes

### phase change comparison elmer/Comsol

Hi,

I work on an industrial phase change problem of silicon in induction melting process.
I work with Comsol but I am very interested in Elmer.
I make some tests with transient simulation to compare Elmer and Comsol. I have seen many tests on this forum for phase change with many issue

First I simulate a 2D problem with solid-solid phase change (spatial2). It's works very good, Elmer and Comsol give same results.

See sif file:

Code: Select all

``````Header
CHECK KEYWORDS Warn
Mesh DB "./Mesh" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Max Iterations = 1
Output Intervals = 100
Timestep intervals = 2000
Timestep Sizes = 10
Timestepping Method = BDF
BDF Order = 2
Solver Input File = case.sif
!Post File = case.ep
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "Body 1"
Equation = 1
Material = 1
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "Body 2"
Equation = 1
Material = 1
Initial condition = 1
End

Solver 1
Equation = SaveScalars
Filename Numbering = True
Operator 2 = mean
Variable 2 = Temperature
Operator 1 = boundary sum
Variable 3 = Time
!Operator 1 = diffusive flux
!Variable 1 = Temperature
!Coefficient 1 = "Heat Conductivity"
Filename = scal_flux_diri.dat
Procedure = "SaveData" "SaveScalars"
Exec Solver = After Timestep
End

Solver 5
Calculate Flux = Logical True
Calculate Flux Abs = Logical True
Target Variable = String "Temperature"
Flux Coefficient = String "Heat Conductivity"
Procedure = "FluxSolver" "FluxSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 3
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 3 = Temperature grad
Vector Field 1 = Temperature flux
Single Precision = True
Exec Solver = After saving
End

Solver 2
Equation = SaveLine
Procedure = "SaveData" "SaveLine"
Filename = saveline_flux_diri.dat
Filename Numbering = True
Save Flux = Logical True
!Flux Variable = String Temperature
!Flux Coefficient = String "Heat Conductivity"
Variable 1 = Coordinate 1
Variable 2 = Coordinate 2
Variable 3 = String Temperature
Variable 4 = Temperature flux 1
Variable 5 = Temperature flux 2
Variable 6 = Time
Polyline Coordinates (2,2) = Real 0.0 0.05 0.1 0.06
Exec Solver = After saving
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 4
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Equation 1
Name = "Equation 1"
Phase Change Model = Spatial 2
Check Latent Heat Release = True
Active Solvers(5) = 4 1 3 2 5
End

Material 1
Name = "Material 1"
Enthalpy = Variable Temperature; Real; 298.15 0.00E+000;350 8.89E+007;401.85 1.82E+008;453.7 2.79E+008;505.55 3.79E+008;557.4 4.80E+008;609.25 5.84E+008;661.1 6.88E+008;712.95 7.95E+008;764.8 9.02E+008;816.65 1.01E+009;868.5 1.12E+009;920.35 1.23E+009;
972.2 1.34E+009;1024.05 1.46E+009;1075.9 1.57E+009;1127.75 1.69E+009;1179.6 1.80E+009;1231.45 1.92E+009;1283.3 2.04E+009;1335.15 2.16E+009;1387 2.28E+009;1438.85 2.40E+009;1490.7 2.52E+009;1542.55 2.64E+009;1594.4 2.77E+009;
1605 2.79E+009;1615.6 2.82E+009;1626.2 2.84E+009;1636.8 2.87E+009;1647.4 2.89E+009;1658 2.92E+009;1668.6 2.94E+009;1679.2 2.97E+009;1680.2 2.97E+009;1696 7.34E+009;1750 7.47E+009;1804 7.59E+009;1858 7.71E+009;1912 7.83E+009;1966 7.95E+009;2020 8.07E+009;2074 8.20E+009;2128 8.32E+009;2182 8.44E+009;2236 8.56E+009;2290 8.68E+009;
2344 8.81E+009;2398 8.93E+009;2452 9.05E+009;2506 9.17E+009;2560 9.29E+009;2614 9.41E+009;2668 9.54E+009;2722 9.66E+009;2776 9.78E+009;2830 9.90E+009;2884 1.00E+010;2938 1.01E+010;2992 1.03E+010;
End
Phase Change Intervals(2,1) = 1680.2 1696
Heat Conductivity = 20
Heat Capacity = 1000
Density = 2330
End

Initial Condition 1
Name = "InitialCondition 1"
Temperature = 300
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1750
Save Scalars = Logical True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
External Temperature = 300
Heat Transfer Coefficient = 15
Save Scalars = Logical True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
External Temperature = 300
Heat Transfer Coefficient = 15
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Heat Flux = -10000
Save Scalars = Logical True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Heat Flux = -10000
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Save Scalars = Logical True
Save Line = Logical True
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_iso"
Save Scalars = Logical True
End
``````
After I had Navier Stockes for a liquid-solid phase change with a change of viscosity between liquid and solid phase

Code: Select all

``````Header
CHECK KEYWORDS Warn
Mesh DB "./Mesh" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Max Iterations = 1
Output Intervals = 20
Timestep intervals = 400
Timestep Sizes = 2
Timestepping Method = BDF
BDF Order = 2
Solver Input File = case.sif
!Post File = case.ep
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "Body 1"
Equation = 1
Material = 2
Body Force = 1
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "Body 2"
Equation = 1
Material = 2
Body Force = 1
Initial condition = 2
End

Solver 1
Equation = SaveScalars
Filename Numbering = True
Operator 2 = mean
Variable 2 = Temperature
Operator 1 = boundary sum
Variable 3 = Time
!Operator 1 = diffusive flux
!Variable 1 = Temperature
!Coefficient 1 = "Heat Conductivity"
Filename = scal_flux_diri.dat
Procedure = "SaveData" "SaveScalars"
Exec Solver = After Timestep
End

Solver 5
Calculate Flux = Logical True
Calculate Flux Abs = Logical True
Target Variable = String "Temperature"
Flux Coefficient = String "Heat Conductivity"
Procedure = "FluxSolver" "FluxSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 3
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 3 = Temperature grad
Vector Field 1 = Temperature flux
Scalar Field 4 = Pressure
Vector Field 2 = Velocity
Single Precision = True
Exec Solver = After saving
End

Solver 2
Equation = SaveLine
Procedure = "SaveData" "SaveLine"
Filename = saveline_flux_diri.dat
Filename Numbering = True
Save Flux = Logical True
!Flux Variable = String Temperature
!Flux Coefficient = String "Heat Conductivity"
Variable 1 = Coordinate 1
Variable 2 = Coordinate 2
Variable 3 = String Temperature
Variable 4 = Temperature flux 1
Variable 5 = Temperature flux 2
Variable 6 = Velocity 1
Variable 7 = Velocity 2
Variable 8 = Pressure
Variable 9 = Time
Polyline Coordinates (2,2) = Real 0.0 0.05 0.1 0.06
Exec Solver = After saving
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 4
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 6
Equation = Navier-Stokes
Variable = Flow Solution[Velocity:2 Pressure:1]
Procedure = "FlowSolve" "FlowSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Equation 1
Name = "Equation 1"
Phase Change Model = Spatial 2
Check Latent Heat Release = True
Convection = Computed
Active Solvers(6) = 4 1 3 2 5 6
End

Material 2
Name = "Silicium_liquide"
Reference Temperature = 293
Viscosity = Variable Temperature; Real; 293.15 1.00E+007;1680 1.00E+007;1680.2 9750000;1680.6 9500000;1681 9250000;1681.4 9000000;1681.8 8750000;1682.2 8500000;
1682.6 8250000;1683 8000000;1683.4 7750000;1683.8 7500000;1684.2 7250000;1684.6 7000000;1685 6750000;1685.4 6500000;
1685.8 6250000;1686.2 6000000;1686.6 5750000;1687 5500000;1687.4 5250000;1687.8 5000000;1688.2 4750000;1688.6 4500000;
1689 4250000;1689.4 4000000;1689.8 3750000;1690.2 3500000;1690.6 3250000;1691 3000000;1691.4 2750000;1691.8 2500000;
1692.2 2250000;1692.6 2000000;1693 1750000;1693.4 1500000;1693.8 1250000;1694.2 1000000;1694.6 750000;1695 500000;
1695.4 250000;1695.8 8.80E-004;3000 8.80E-004;
End
Heat expansion Coefficient = 1.3725e-4
Enthalpy = Variable Temperature; Real; 298.15 0.00E+000;350 8.89E+007;401.85 1.82E+008;453.7 2.79E+008;505.55 3.79E+008;557.4 4.80E+008;609.25 5.84E+008;661.1 6.88E+008;712.95 7.95E+008;764.8 9.02E+008;816.65 1.01E+009;868.5 1.12E+009;920.35 1.23E+009;
972.2 1.34E+009;1024.05 1.46E+009;1075.9 1.57E+009;1127.75 1.69E+009;1179.6 1.80E+009;1231.45 1.92E+009;1283.3 2.04E+009;1335.15 2.16E+009;1387 2.28E+009;1438.85 2.40E+009;1490.7 2.52E+009;1542.55 2.64E+009;1594.4 2.77E+009;
1605 2.79E+009;1615.6 2.82E+009;1626.2 2.84E+009;1636.8 2.87E+009;1647.4 2.89E+009;1658 2.92E+009;1668.6 2.94E+009;1679.2 2.97E+009;1680.2 2.97E+009;1696 7.34E+009;1750 7.47E+009;1804 7.59E+009;1858 7.71E+009;1912 7.83E+009;1966 7.95E+009;2020 8.07E+009;2074 8.20E+009;2128 8.32E+009;2182 8.44E+009;2236 8.56E+009;2290 8.68E+009;
2344 8.81E+009;2398 8.93E+009;2452 9.05E+009;2506 9.17E+009;2560 9.29E+009;2614 9.41E+009;2668 9.54E+009;2722 9.66E+009;2776 9.78E+009;2830 9.90E+009;2884 1.00E+010;2938 1.01E+010;2992 1.03E+010;
End
Phase Change Intervals(2,1) = 1680.2 1696
Compressibility Model = Incompressible
Reference Pressure = 0
Specific Heat Ratio = 1.4
Heat Conductivity = 20
Heat Capacity = 1000
Density = 2550
End

Body Force 1
Name = "Natural convection"
Boussinesq = True
End

Initial Condition 1
Name = "InitialCondition 1"
Velocity 2 = 0
Pressure = 0
Velocity 1 = 0
Temperature = 1715
End

Initial Condition 2
Name = "InitialCondition 2"
Velocity 2 = 0
Pressure = 0
Velocity 1 = 0
Temperature = 1655
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1715
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
External Temperature = 300
Heat Transfer Coefficient = 15
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
External Temperature = 300
Heat Transfer Coefficient = 15
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Heat Flux = -10000
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Heat Flux = -10000
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Save Scalars = Logical True
Save Line = Logical True
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_iso"
Temperature = 1655
Noslip wall BC = True
Save Scalars = Logical True
End
``````
It works very well too, computation time is more than 10 times faster with Elmer!

Now I want to try phasechange solver in a transient problem because I have an isothermal phase change problem. I have some issue of dimensions I think.
In spatial2, enthalpy is given in [J/m^3] right and Heat capacity in [J/(kg.K)]?
In phasechange solver what are units? for Heat capacity and Latent Heat?
I try Heat capacity in [J/(kg.K)] and Latent Heat [J/m^3]
or Heat capacity in [J/(m^3.K)] and Latent Heat [J/m^3], there is still a problem, solid liquid interface doesn't move. If I decrease Latent Heat by some factors interface move too fast...
See my sif file adapted from github phasechange3

Code: Select all

``````! Compute transient phase change problem by updating the
! interface based on the melting speed.
\$ velo = 1.0
\$ melt = 1685
\$ trans = 5
\$ cap = 2.33e6
\$ lat = 4.19e9
\$ cond = 20

Mesh DB "." "PhaseChange"
Results Directory "./Results"
End

Simulation
Max Output Level = 3
Coordinate System = "Axi Symmetric"
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Min Iterations = 1
Steady State Max Iterations = 1
Timestepping Method = Implicit Euler
Timestep Sizes = 10
Timestep Intervals = 300
!Timestep Intervals = 100
Output Intervals = 10
Post File = "data.ep"
Output File = "data.result"
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
End

Body 1
Name = "solid"
Equation = 1
Material = 1
Initial Condition = 1
End

Body 2
Name = "melt"
Equation = 1
Material = 2
Initial Condition = 1
End

Body 3
Name = "heater"
Equation = 2
Material = 3
Initial Condition = 1
!Body Force = 1
End

Body 4
Name = "interface"
Equation = 3
Material = 1
Initial Condition = 1
End

Equation 1
Active Solvers(2) = 1 2
Convection = Constant
End

Equation 2
Active Solvers(1) = 2
Convection = Constant
End

Equation 3
Active Solvers(1) = 3
End

Solver 1
Equation = "Mesh Update"
Linear System Solver = "Direct"
Linear System Direct Method = "umfpack"
Steady State Convergence Tolerance = 1.0e-4
End

Solver 2
Equation = "Heat Equation"
Linear System Solver = "Direct"
Linear System Direct Method = "umfpack"
Nonlinear System Convergence Tolerance = 1.0e-5
Nonlinear System Max Iterations = 1
Nonlinear System Relaxation Factor = 1.0
Nonlinear System Newton After Iterations = 0
Nonlinear System Newton After Tolerance = 1.0e-2
Steady State Convergence Tolerance = 1.0e-4
Stabilize = Logical True
End

Solver 3
Variable = PhaseSurface
Equation = "Phase Surface"
! Procedure = "PhaseChangeSolve" "PhaseChangeSolve"
Procedure = "TransientPhaseChange" "TransientPhaseChange"
Nonlinear System Relaxation Factor = 1.0
Nonlinear System Newton After Iterations = 10
! Lumped Newton After Iterations = Integer 10
Steady State Convergence Tolerance = 1.0e-3
! Use Triple Point for Melting Point = Logical True
Surface Smoothing Factor = Real 0.0
Transient Speedup = Real 1.0
Velocity Smoothing Factor = Real 0.05
Velocity Relaxation Factor = Real 0.5
End

Solver 4
Exec Solver = After Timestep
Equation = String SaveLine
Procedure = File "SaveData" "SaveLine"
Filename = File "ss.dat"
File Append = Logical True
End

Solver 5
Exec Solver = After Timestep
Equation = String SaveScalars
Procedure = File "SaveData" "SaveScalars"
Filename = File "fs.dat"
Variable 1 = String PhaseSurface
Operator 1 = String max
Variable 2 = String PhaseSurface
Operator 2 = String min
Variable 3 = String Temperature
Operator 3 = String mean
Variable 4 = String time
End

Solver 6
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = case
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 3 = Temperature grad
Scalar Field 4 = PhaseSurface
Vector Field 1 = Temperature flux
Single Precision = True
Exec Solver = After saving
End

Body Force 1
! Heat Source = Real 1.0
! Smart Heater Control = Logical True
End

Material 1
Solid = Logical True
Melting Point = Real \$ melt
Density = 2330
Heat Capacity = \$ cap
Heat Conductivity = \$ cond
!Youngs Modulus = 1.0
!Poisson Ratio = 0.3
Latent Heat = Real \$ lat
!Convection Velocity 1 = 0.0
!Convection Velocity 2 = Real \$ velo
End

Material 2
Liquid = Logical True
Melting Point = Real \$ melt
Density = 2330
Heat Capacity = \$ cap
Heat Conductivity = \$ cond
!Youngs Modulus = 1.0
!Poisson Ratio = 0.3
!Convection Velocity 1 = 0.0
!Convection Velocity 2 = Real \$ velo
End

Material 3
Density = 2330
Heat Capacity = \$ cap
Heat Conductivity = 20.0
End

Initial Condition 1
Temperature = Variable Coordinate 2
Real
0.0 1700.0
2.0 1650
End

PhaseSurface = Real 0.0
Mesh Update 1 = Real 0.0
Mesh Update 2 = Real 0.0
End

Boundary Condition 1
Name = "melt_crystal"
Target Boundaries = 1
Temperature = Real \$ melt
Mesh Update 1 = 0
Mesh Update 2 = Equals PhaseSurface
Save Line = Logical True
Normal Target Body = Integer 1
Body Id = Integer 4
End

Boundary Condition 2
Name = "melt_heater"
Target Boundaries = 2
Mesh Update 1 = 0
Mesh Update 2 = 0
Temperature = 1700
End

Boundary Condition 3
Name = "solid_top"
Target Boundaries = 3
Mesh Update 1 = 0
Mesh Update 2 = 0
Temperature = 1650
End

Boundary Condition 4
Name = "right_walls"
Target Boundaries = 4
Mesh Update 1 = 0
!Heat Flux BC = Logical True
Heat Transfer Coefficient = Real \$ trans
External Temperature = Real 300
Phase Change Side = Logical True
End

Boundary Condition 5
Name = "axis"
Target Boundaries = 5
Mesh Update 1 = 0
End

Boundary Condition 6
Name = "heater_bottom"
Target Boundaries = 6
End``````
Hope someone can help

Julien
Attachments
phase_change_elmer_navier_stockes_20000s.png
elmer phasechange solid-liquid spatial2
phase_change_comsol_navier_stockes_20000s.jpg
Comsol phasechange solid-liquid apparent heat capicity model
phase_change_elmer_20000s.png
Elmer phasechange solid-solid spatial2
julien givernaud
Posts: 71
Joined: 18 Nov 2014, 18:18
Antispam: Yes

### Re: phase change comparison elmer/Comsol

Sorry for the poor presentation with sif file. I do it better:
First I simulate a 2D problem with solid-solid phase change (spatial2). It's works very good, Elmer and Comsol give same results.

Code: Select all

``````Header
CHECK KEYWORDS Warn
Mesh DB "./Mesh" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Max Iterations = 1
Output Intervals = 100
Timestep intervals = 2000
Timestep Sizes = 10
Timestepping Method = BDF
BDF Order = 2
Solver Input File = case.sif
!Post File = case.ep
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "Body 1"
Equation = 1
Material = 1
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "Body 2"
Equation = 1
Material = 1
Initial condition = 1
End

Solver 1
Equation = SaveScalars
Filename Numbering = True
Operator 2 = mean
Variable 2 = Temperature
Operator 1 = boundary sum
Variable 3 = Time
!Operator 1 = diffusive flux
!Variable 1 = Temperature
!Coefficient 1 = "Heat Conductivity"
Filename = scal_flux_diri.dat
Procedure = "SaveData" "SaveScalars"
Exec Solver = After Timestep
End

Solver 5
Calculate Flux = Logical True
Calculate Flux Abs = Logical True
Target Variable = String "Temperature"
Flux Coefficient = String "Heat Conductivity"
Procedure = "FluxSolver" "FluxSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 3
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 3 = Temperature grad
Vector Field 1 = Temperature flux
Single Precision = True
Exec Solver = After saving
End

Solver 2
Equation = SaveLine
Procedure = "SaveData" "SaveLine"
Filename = saveline_flux_diri.dat
Filename Numbering = True
Save Flux = Logical True
!Flux Variable = String Temperature
!Flux Coefficient = String "Heat Conductivity"
Variable 1 = Coordinate 1
Variable 2 = Coordinate 2
Variable 3 = String Temperature
Variable 4 = Temperature flux 1
Variable 5 = Temperature flux 2
Variable 6 = Time
Polyline Coordinates (2,2) = Real 0.0 0.05 0.1 0.06
Exec Solver = After saving
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 4
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Equation 1
Name = "Equation 1"
Phase Change Model = Spatial 2
Check Latent Heat Release = True
Active Solvers(5) = 4 1 3 2 5
End

Material 1
Name = "Material 1"
Enthalpy = Variable Temperature; Real; 298.15 0.00E+000;350 8.89E+007;401.85 1.82E+008;453.7 2.79E+008;505.55 3.79E+008;557.4 4.80E+008;609.25 5.84E+008;661.1 6.88E+008;712.95 7.95E+008;764.8 9.02E+008;816.65 1.01E+009;868.5 1.12E+009;920.35 1.23E+009;
972.2 1.34E+009;1024.05 1.46E+009;1075.9 1.57E+009;1127.75 1.69E+009;1179.6 1.80E+009;1231.45 1.92E+009;1283.3 2.04E+009;1335.15 2.16E+009;1387 2.28E+009;1438.85 2.40E+009;1490.7 2.52E+009;1542.55 2.64E+009;1594.4 2.77E+009;
1605 2.79E+009;1615.6 2.82E+009;1626.2 2.84E+009;1636.8 2.87E+009;1647.4 2.89E+009;1658 2.92E+009;1668.6 2.94E+009;1679.2 2.97E+009;1680.2 2.97E+009;1696 7.34E+009;1750 7.47E+009;1804 7.59E+009;1858 7.71E+009;1912 7.83E+009;1966 7.95E+009;2020 8.07E+009;2074 8.20E+009;2128 8.32E+009;2182 8.44E+009;2236 8.56E+009;2290 8.68E+009;
2344 8.81E+009;2398 8.93E+009;2452 9.05E+009;2506 9.17E+009;2560 9.29E+009;2614 9.41E+009;2668 9.54E+009;2722 9.66E+009;2776 9.78E+009;2830 9.90E+009;2884 1.00E+010;2938 1.01E+010;2992 1.03E+010;
End
Phase Change Intervals(2,1) = 1680.2 1696
Heat Conductivity = 20
Heat Capacity = 1000
Density = 2330
End

Initial Condition 1
Name = "InitialCondition 1"
Temperature = 300
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1750
Save Scalars = Logical True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
External Temperature = 300
Heat Transfer Coefficient = 15
Save Scalars = Logical True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
External Temperature = 300
Heat Transfer Coefficient = 15
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Heat Flux = -10000
Save Scalars = Logical True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Heat Flux = -10000
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Save Scalars = Logical True
Save Line = Logical True
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_iso"
Save Scalars = Logical True
End``````
After I had Navier Stockes for a liquid-solid phase change with a change of viscosity between liquid and solid phase and natural convection for the body force

Code: Select all

``````Header
CHECK KEYWORDS Warn
Mesh DB "./Mesh" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Max Iterations = 1
Output Intervals = 20
Timestep intervals = 400
Timestep Sizes = 2
Timestepping Method = BDF
BDF Order = 2
Solver Input File = case.sif
!Post File = case.ep
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "Body 1"
Equation = 1
Material = 2
Body Force = 1
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "Body 2"
Equation = 1
Material = 2
Body Force = 1
Initial condition = 2
End

Solver 1
Equation = SaveScalars
Filename Numbering = True
Operator 2 = mean
Variable 2 = Temperature
Operator 1 = boundary sum
Variable 3 = Time
!Operator 1 = diffusive flux
!Variable 1 = Temperature
!Coefficient 1 = "Heat Conductivity"
Filename = scal_flux_diri.dat
Procedure = "SaveData" "SaveScalars"
Exec Solver = After Timestep
End

Solver 5
Calculate Flux = Logical True
Calculate Flux Abs = Logical True
Target Variable = String "Temperature"
Flux Coefficient = String "Heat Conductivity"
Procedure = "FluxSolver" "FluxSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 3
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 3 = Temperature grad
Vector Field 1 = Temperature flux
Scalar Field 4 = Pressure
Vector Field 2 = Velocity
Single Precision = True
Exec Solver = After saving
End

Solver 2
Equation = SaveLine
Procedure = "SaveData" "SaveLine"
Filename = saveline_flux_diri.dat
Filename Numbering = True
Save Flux = Logical True
!Flux Variable = String Temperature
!Flux Coefficient = String "Heat Conductivity"
Variable 1 = Coordinate 1
Variable 2 = Coordinate 2
Variable 3 = String Temperature
Variable 4 = Temperature flux 1
Variable 5 = Temperature flux 2
Variable 6 = Velocity 1
Variable 7 = Velocity 2
Variable 8 = Pressure
Variable 9 = Time
Polyline Coordinates (2,2) = Real 0.0 0.05 0.1 0.06
Exec Solver = After saving
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 4
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 6
Equation = Navier-Stokes
Variable = Flow Solution[Velocity:2 Pressure:1]
Procedure = "FlowSolve" "FlowSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Equation 1
Name = "Equation 1"
Phase Change Model = Spatial 2
Check Latent Heat Release = True
Convection = Computed
Active Solvers(6) = 4 1 3 2 5 6
End

Material 2
Name = "Silicium_liquide"
Reference Temperature = 293
Viscosity = Variable Temperature; Real; 293.15 1.00E+007;1680 1.00E+007;1680.2 9750000;1680.6 9500000;1681 9250000;1681.4 9000000;1681.8 8750000;1682.2 8500000;
1682.6 8250000;1683 8000000;1683.4 7750000;1683.8 7500000;1684.2 7250000;1684.6 7000000;1685 6750000;1685.4 6500000;
1685.8 6250000;1686.2 6000000;1686.6 5750000;1687 5500000;1687.4 5250000;1687.8 5000000;1688.2 4750000;1688.6 4500000;
1689 4250000;1689.4 4000000;1689.8 3750000;1690.2 3500000;1690.6 3250000;1691 3000000;1691.4 2750000;1691.8 2500000;
1692.2 2250000;1692.6 2000000;1693 1750000;1693.4 1500000;1693.8 1250000;1694.2 1000000;1694.6 750000;1695 500000;
1695.4 250000;1695.8 8.80E-004;3000 8.80E-004;
End
Heat expansion Coefficient = 1.3725e-4
Enthalpy = Variable Temperature; Real; 298.15 0.00E+000;350 8.89E+007;401.85 1.82E+008;453.7 2.79E+008;505.55 3.79E+008;557.4 4.80E+008;609.25 5.84E+008;661.1 6.88E+008;712.95 7.95E+008;764.8 9.02E+008;816.65 1.01E+009;868.5 1.12E+009;920.35 1.23E+009;
972.2 1.34E+009;1024.05 1.46E+009;1075.9 1.57E+009;1127.75 1.69E+009;1179.6 1.80E+009;1231.45 1.92E+009;1283.3 2.04E+009;1335.15 2.16E+009;1387 2.28E+009;1438.85 2.40E+009;1490.7 2.52E+009;1542.55 2.64E+009;1594.4 2.77E+009;
1605 2.79E+009;1615.6 2.82E+009;1626.2 2.84E+009;1636.8 2.87E+009;1647.4 2.89E+009;1658 2.92E+009;1668.6 2.94E+009;1679.2 2.97E+009;1680.2 2.97E+009;1696 7.34E+009;1750 7.47E+009;1804 7.59E+009;1858 7.71E+009;1912 7.83E+009;1966 7.95E+009;2020 8.07E+009;2074 8.20E+009;2128 8.32E+009;2182 8.44E+009;2236 8.56E+009;2290 8.68E+009;
2344 8.81E+009;2398 8.93E+009;2452 9.05E+009;2506 9.17E+009;2560 9.29E+009;2614 9.41E+009;2668 9.54E+009;2722 9.66E+009;2776 9.78E+009;2830 9.90E+009;2884 1.00E+010;2938 1.01E+010;2992 1.03E+010;
End
Phase Change Intervals(2,1) = 1680.2 1696
Compressibility Model = Incompressible
Reference Pressure = 0
Specific Heat Ratio = 1.4
Heat Conductivity = 20
Heat Capacity = 1000
Density = 2550
End

Body Force 1
Name = "Natural convection"
Boussinesq = True
End

Initial Condition 1
Name = "InitialCondition 1"
Velocity 2 = 0
Pressure = 0
Velocity 1 = 0
Temperature = 1715
End

Initial Condition 2
Name = "InitialCondition 2"
Velocity 2 = 0
Pressure = 0
Velocity 1 = 0
Temperature = 1655
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1715
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
External Temperature = 300
Heat Transfer Coefficient = 15
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
External Temperature = 300
Heat Transfer Coefficient = 15
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Heat Flux = -10000
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Heat Flux = -10000
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Save Scalars = Logical True
Save Line = Logical True
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_iso"
Temperature = 1655
Noslip wall BC = True
Save Scalars = Logical True
End``````
It works very well too, computation time is more than 10 times faster with Elmer!

Now I want to try phasechange solver in a transient problem because I have an isothermal phase change problem. I have some issue of dimensions I think.
In spatial2, enthalpy is given in [J/m^3] right and Heat capacity in [J/(kg.K)]?
In phasechange solver what are units? for Heat capacity and Latent Heat?
I try Heat capacity in [J/(kg.K)] and Latent Heat [J/m^3]
or Heat capacity in [J/(m^3.K)] and Latent Heat [J/m^3], there is still a problem, solid liquid interface doesn't move. If I decrease Latent Heat by some factors interface move too fast...
See my sif file adapted from github phasechange3

Code: Select all

`````` Compute transient phase change problem by updating the
! interface based on the melting speed.
\$ velo = 1.0
\$ melt = 1685
\$ trans = 5
\$ cap = 2.33e6
\$ lat = 4.19e9
\$ cond = 20

Mesh DB "." "PhaseChange"
Results Directory "./Results"
End

Simulation
Max Output Level = 3
Coordinate System = "Axi Symmetric"
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Min Iterations = 1
Steady State Max Iterations = 1
Timestepping Method = Implicit Euler
Timestep Sizes = 10
Timestep Intervals = 300
!Timestep Intervals = 100
Output Intervals = 10
Post File = "data.ep"
Output File = "data.result"
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
End

Body 1
Name = "solid"
Equation = 1
Material = 1
Initial Condition = 1
End

Body 2
Name = "melt"
Equation = 1
Material = 2
Initial Condition = 1
End

Body 3
Name = "heater"
Equation = 2
Material = 3
Initial Condition = 1
!Body Force = 1
End

Body 4
Name = "interface"
Equation = 3
Material = 1
Initial Condition = 1
End

Equation 1
Active Solvers(2) = 1 2
Convection = Constant
End

Equation 2
Active Solvers(1) = 2
Convection = Constant
End

Equation 3
Active Solvers(1) = 3
End

Solver 1
Equation = "Mesh Update"
Linear System Solver = "Direct"
Linear System Direct Method = "umfpack"
Steady State Convergence Tolerance = 1.0e-4
End

Solver 2
Equation = "Heat Equation"
Linear System Solver = "Direct"
Linear System Direct Method = "umfpack"
Nonlinear System Convergence Tolerance = 1.0e-5
Nonlinear System Max Iterations = 1
Nonlinear System Relaxation Factor = 1.0
Nonlinear System Newton After Iterations = 0
Nonlinear System Newton After Tolerance = 1.0e-2
Steady State Convergence Tolerance = 1.0e-4
Stabilize = Logical True
End

Solver 3
Variable = PhaseSurface
Equation = "Phase Surface"
! Procedure = "PhaseChangeSolve" "PhaseChangeSolve"
Procedure = "TransientPhaseChange" "TransientPhaseChange"
Nonlinear System Relaxation Factor = 1.0
Nonlinear System Newton After Iterations = 10
! Lumped Newton After Iterations = Integer 10
Steady State Convergence Tolerance = 1.0e-3
! Use Triple Point for Melting Point = Logical True
Surface Smoothing Factor = Real 0.0
Transient Speedup = Real 1.0
Velocity Smoothing Factor = Real 0.05
Velocity Relaxation Factor = Real 0.5
End

Solver 4
Exec Solver = After Timestep
Equation = String SaveLine
Procedure = File "SaveData" "SaveLine"
Filename = File "ss.dat"
File Append = Logical True
End

Solver 5
Exec Solver = After Timestep
Equation = String SaveScalars
Procedure = File "SaveData" "SaveScalars"
Filename = File "fs.dat"
Variable 1 = String PhaseSurface
Operator 1 = String max
Variable 2 = String PhaseSurface
Operator 2 = String min
Variable 3 = String Temperature
Operator 3 = String mean
Variable 4 = String time
End

Solver 6
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = case
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 3 = Temperature grad
Scalar Field 4 = PhaseSurface
Vector Field 1 = Temperature flux
Single Precision = True
Exec Solver = After saving
End

Body Force 1
! Heat Source = Real 1.0
! Smart Heater Control = Logical True
End

Material 1
Solid = Logical True
Melting Point = Real \$ melt
Density = 2330
Heat Capacity = \$ cap
Heat Conductivity = \$ cond
!Youngs Modulus = 1.0
!Poisson Ratio = 0.3
Latent Heat = Real \$ lat
!Convection Velocity 1 = 0.0
!Convection Velocity 2 = Real \$ velo
End

Material 2
Liquid = Logical True
Melting Point = Real \$ melt
Density = 2330
Heat Capacity = \$ cap
Heat Conductivity = \$ cond
!Youngs Modulus = 1.0
!Poisson Ratio = 0.3
!Convection Velocity 1 = 0.0
!Convection Velocity 2 = Real \$ velo
End

Material 3
Density = 2330
Heat Capacity = \$ cap
Heat Conductivity = 20.0
End

Initial Condition 1
Temperature = Variable Coordinate 2
Real
0.0 1700.0
2.0 1650
End

PhaseSurface = Real 0.0
Mesh Update 1 = Real 0.0
Mesh Update 2 = Real 0.0
End

Boundary Condition 1
Name = "melt_crystal"
Target Boundaries = 1
Temperature = Real \$ melt
Mesh Update 1 = 0
Mesh Update 2 = Equals PhaseSurface
Save Line = Logical True
Normal Target Body = Integer 1
Body Id = Integer 4
End

Boundary Condition 2
Name = "melt_heater"
Target Boundaries = 2
Mesh Update 1 = 0
Mesh Update 2 = 0
Temperature = 1700
End

Boundary Condition 3
Name = "solid_top"
Target Boundaries = 3
Mesh Update 1 = 0
Mesh Update 2 = 0
Temperature = 1650
End

Boundary Condition 4
Name = "right_walls"
Target Boundaries = 4
Mesh Update 1 = 0
!Heat Flux BC = Logical True
Heat Transfer Coefficient = Real \$ trans
External Temperature = Real 300
Phase Change Side = Logical True
End

Boundary Condition 5
Name = "axis"
Target Boundaries = 5
Mesh Update 1 = 0
End

Boundary Condition 6
Name = "heater_bottom"
Target Boundaries = 6
End``````
Julien
julien givernaud
Posts: 71
Joined: 18 Nov 2014, 18:18
Antispam: Yes

### Re: phase change comparison elmer/Comsol

To compare transient spatial2 and phasechangesolver in a 2D case solid-solid phase change with same thermodynamical properties, I try Heat capacity in J/(kg.K) and Latent Heat in J/kg in phasechangesolver case.
Enthalpy in spatial2 is in J/m^3 (comparison ok with Comsol.
Mesh, initial conditions, boundaries conditions are the same.

For the spatial2 case I use the following sif file

Code: Select all

``````Header
CHECK KEYWORDS Warn
Mesh DB "./Mesh" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Max Iterations = 1
Output Intervals = 100
Timestep intervals = 2000
Timestep Sizes = 10
Timestepping Method = BDF
BDF Order = 2
Solver Input File = case.sif
!Post File = case.ep
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "Body 1"
Equation = 1
Material = 1
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "Body 2"
Equation = 1
Material = 1
Initial condition = 1
End

Solver 1
Equation = SaveScalars
Filename Numbering = True
Operator 2 = mean
Variable 2 = Temperature
Operator 1 = boundary sum
Variable 3 = Time
!Operator 1 = diffusive flux
!Variable 1 = Temperature
!Coefficient 1 = "Heat Conductivity"
Filename = scal_flux_diri.dat
Procedure = "SaveData" "SaveScalars"
Exec Solver = After Timestep
End

Solver 5
Calculate Flux = Logical True
Calculate Flux Abs = Logical True
Target Variable = String "Temperature"
Flux Coefficient = String "Heat Conductivity"
Procedure = "FluxSolver" "FluxSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 3
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 3 = Temperature grad
Vector Field 1 = Temperature flux
Single Precision = True
Exec Solver = After saving
End

Solver 2
Equation = SaveLine
Procedure = "SaveData" "SaveLine"
Filename = saveline_flux_diri.dat
Filename Numbering = True
Save Flux = Logical True
!Flux Variable = String Temperature
!Flux Coefficient = String "Heat Conductivity"
Variable 1 = Coordinate 1
Variable 2 = Coordinate 2
Variable 3 = String Temperature
Variable 4 = Temperature flux 1
Variable 5 = Temperature flux 2
Variable 6 = Time
Polyline Coordinates (2,2) = Real 0.0 0.05 0.1 0.06
Exec Solver = After saving
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 4
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Equation 1
Name = "Equation 1"
Phase Change Model = Spatial 2
Check Latent Heat Release = True
Active Solvers(5) = 4 1 3 2 5
End

Material 1
Name = "Material 1"
Enthalpy = Variable Temperature; Real; 298.15 0.00E+000;350 8.89E+007;401.85 1.82E+008;453.7 2.79E+008;505.55 3.79E+008;557.4 4.80E+008;609.25 5.84E+008;661.1 6.88E+008;712.95 7.95E+008;764.8 9.02E+008;816.65 1.01E+009;868.5 1.12E+009;920.35 1.23E+009;
972.2 1.34E+009;1024.05 1.46E+009;1075.9 1.57E+009;1127.75 1.69E+009;1179.6 1.80E+009;1231.45 1.92E+009;1283.3 2.04E+009;1335.15 2.16E+009;1387 2.28E+009;1438.85 2.40E+009;1490.7 2.52E+009;1542.55 2.64E+009;1594.4 2.77E+009;
1605 2.79E+009;1615.6 2.82E+009;1626.2 2.84E+009;1636.8 2.87E+009;1647.4 2.89E+009;1658 2.92E+009;1668.6 2.94E+009;1679.2 2.97E+009;1680.2 2.97E+009;1696 7.34E+009;1750 7.47E+009;1804 7.59E+009;1858 7.71E+009;1912 7.83E+009;1966 7.95E+009;2020 8.07E+009;2074 8.20E+009;2128 8.32E+009;2182 8.44E+009;2236 8.56E+009;2290 8.68E+009;
2344 8.81E+009;2398 8.93E+009;2452 9.05E+009;2506 9.17E+009;2560 9.29E+009;2614 9.41E+009;2668 9.54E+009;2722 9.66E+009;2776 9.78E+009;2830 9.90E+009;2884 1.00E+010;2938 1.01E+010;2992 1.03E+010;
End
Phase Change Intervals(2,1) = 1680.2 1696
Heat Conductivity = 20
Heat Capacity = 1000
Density = 2330
End

Initial Condition 1
Name = "InitialCondition 1"
Temperature = Variable Coordinate 1
Real
0.0 1745
0.1 1645
End
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1750
Save Scalars = Logical True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
External Temperature = 300
Heat Transfer Coefficient = 15
Save Scalars = Logical True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
External Temperature = 300
Heat Transfer Coefficient = 15
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Heat Flux = -10000
Save Scalars = Logical True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Heat Flux = -10000
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Save Scalars = Logical True
Save Line = Logical True
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_iso"
Save Scalars = Logical True
End``````
For phasechangesolver:

Code: Select all

``````Header
CHECK KEYWORDS Warn
Mesh DB "./Mesh" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Min Iterations = 1
Steady State Max Iterations = 1
Output Intervals = 100
Timestep intervals = 2000
Timestep Sizes = 10
Timestepping Method = BDF
!Timestepping Method = Implicit Euler
BDF Order = 2
Solver Input File = case.sif
!Post File = case.ep
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "liquide"
Equation = 1
Material = 2
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "solide"
Equation = 1
Material = 1
Initial condition = 1
End

Body 3
Target Bodies(1) = 4
Name = "interface"
Equation = 2
Material = 1
Initial condition = 1
End

Solver 1
Equation = SaveScalars
Filename Numbering = True
Operator 2 = mean
Variable 2 = Temperature
Operator 1 = boundary sum
Variable 3 = Time
!Operator 1 = diffusive flux
!Variable 1 = Temperature
!Coefficient 1 = "Heat Conductivity"
Filename = scal_flux_diri.dat
Procedure = "SaveData" "SaveScalars"
Exec Solver = After Timestep
End

Solver 5
Calculate Flux = Logical True
Calculate Flux Abs = Logical True
Target Variable = String "Temperature"
Flux Coefficient = String "Heat Conductivity"
Procedure = "FluxSolver" "FluxSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 3
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 3 = Temperature grad
Scalar Field 4 = PhaseSurface
Vector Field 1 = Temperature flux
Single Precision = True
Exec Solver = After saving
End

Solver 2
Equation = SaveLine
Procedure = "SaveData" "SaveLine"
Filename = saveline_flux_diri.dat
Filename Numbering = True
Save Flux = Logical True
!Flux Variable = String Temperature
!Flux Coefficient = String "Heat Conductivity"
Variable 1 = Coordinate 1
Variable 2 = Coordinate 2
Variable 3 = String Temperature
Variable 4 = Temperature flux 1
Variable 5 = Temperature flux 2
Variable 6 = Time
Polyline Coordinates (2,2) = Real 0.0 0.05 0.1 0.06
Exec Solver = After saving
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 4
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 6
Equation = "Mesh Update"
Linear System Solver = "Direct"
Linear System Direct Method = "umfpack"
Steady State Convergence Tolerance = 1.0e-4
End

Solver 7
Variable = PhaseSurface
Equation = "Phase Surface"
! Procedure = "PhaseChangeSolve" "PhaseChangeSolve"
Procedure = "TransientPhaseChange" "TransientPhaseChange"
Nonlinear System Relaxation Factor = 1.0
Nonlinear System Newton After Iterations = 10
! Lumped Newton After Iterations = Integer 10
Steady State Convergence Tolerance = 1.0e-3
! Use Triple Point for Melting Point = Logical True
Surface Smoothing Factor = Real 0.0
Transient Speedup = Real 1.0
Velocity Smoothing Factor = Real 0.05
Velocity Relaxation Factor = Real 0.5
End

Equation 1
Name = "Equation 1"
Active Solvers(6) = 4 1 3 2 5 6
!Convection = Constant
End

Equation 2
Name = "Equation 2"
Active Solvers(1) = 7
!Convection = Constant
End

Material 1
Name = "Solid"
Solid = Logical True
Melting Point = Real 1685
!Young Modulus = 150000000000
!Poisson Ratio = 0.27
Heat Conductivity = 20.0
Heat Capacity = 1000.0
Density = 2330.0
Latent Heat = Real 1800e3
!Convection Velocity 1 = 0.0
!Convection Velocity 2 = 0.0
End

Material 2
Name = "Liquid"
Liquid = Logical True
Melting Point = Real 1685
!Young Modulus = 150000000000
!Poisson Ratio = 0.27
Heat Conductivity = 20.0
Heat Capacity = 1000.0
Density = 2330.0
!Convection Velocity 1 = 0.0
!Convection Velocity 2 = 0.0
End

Initial Condition 1
Name = "InitialCondition 1"
Temperature = Variable Coordinate 1
Real
0.0 1745
0.1 1645
End

PhaseSurface = Real 0.0
Mesh Update 1 = Real 0.0
Mesh Update 2 = Real 0.0
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1750
Mesh Update 1 = 0
Mesh Update 2 = 0
Save Scalars = Logical True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
External Temperature = 300
Heat Transfer Coefficient = 15
Mesh Update 2 = 0
!Phase Change Side = Logical True
Save Scalars = Logical True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
External Temperature = 300
Heat Transfer Coefficient = 15
Mesh Update 2 = 0
!Phase Change Side = Logical True
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Heat Flux = -10000
Mesh Update 2 = 0
!Phase Change Side = Logical True
Save Scalars = Logical True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Heat Flux = -10000
Mesh Update 2 = 0
!Phase Change Side = Logical True
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Temperature = Real 1685.0
Mesh Update 1 = Equals PhaseSurface
Mesh Update 2 = 0
Normal Target Body = Integer 2
!Normal Target Body = Integer 1
Save Scalars = Logical True
Save Line = Logical True
Body Id = Integer 3
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_cold"
!Temperature = 1665
Mesh Update 1 = 0
Mesh Update 2 = 0
Save Scalars = Logical True
End
``````
Results are different but they should be indentical.

Julien
Attachments
Mesh.tar.gz
phasechangesolver_300s_solide-solide.png
spatial2_300s_solide-solide.png
julien givernaud
Posts: 71
Joined: 18 Nov 2014, 18:18
Antispam: Yes

### Re: phase change comparison elmer/Comsol

Ok, I have find the problem.
Perpendicular thermal boundaries conditions to the solidification interface can't be defined as easier in phasechange solver than phase change spatial2 model.*
After verification, same problem occurs in Comsol, results are the same for the 2 softwares.

Have a good day

Julien
raback
Posts: 3774
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

### Re: phase change comparison elmer/Comsol

Hi Julien

I had a quick look at your results and must gratulate you on a good work. These solvers were developed fairly long time ago when we had active projects in the area of crystal growth. They are not used much currently by the developer team but the consistency test should ensure that nothing gets broken. Combination of N-S and phase change has not been extensively studied and it is great that you got some results that were agreeable with Comsol.

-Peter
julien givernaud
Posts: 71
Joined: 18 Nov 2014, 18:18
Antispam: Yes

### Re: phase change comparison elmer/Comsol

Hi Peter,

Thank you for your response. I will post a little summary of my results of phase change tests with different solvers and coupling (Navier Stockes or not) associated with comsol comparisons for other users interested in phase change problems.

Julien
annier
Posts: 1141
Joined: 27 Aug 2013, 13:51
Antispam: Yes

### Re: phase change comparison elmer/Comsol

Hi Julien,
Would you please post some theoretical descriptions of the physical equations that you used for modeling in Elmer? This documentation would help the Elmer users in understanding the application fields of Elmer.

Yours
Anil Kunwar
Last edited by annier on 13 Mar 2015, 12:35, edited 1 time in total.
Anil Kunwar
Department of Materials Engineering, KU Leuven, Belgium
julien givernaud
Posts: 71
Joined: 18 Nov 2014, 18:18
Antispam: Yes

### Re: phase change comparison elmer/Comsol

Hi Anil,

Yes I will. I can make a document with a geometry case test with different phasechange models (Elmer vs Comsol), equations and sif files joined.

Julien
annier
Posts: 1141
Joined: 27 Aug 2013, 13:51
Antispam: Yes

### Re: phase change comparison elmer/Comsol

Dear Julien,
Thank you very much.

Yours
Anil Kunwar
Anil Kunwar
Department of Materials Engineering, KU Leuven, Belgium
julien givernaud
Posts: 71
Joined: 18 Nov 2014, 18:18
Antispam: Yes

### Re: phase change comparison elmer/Comsol

Hi,

I attach a summary of phase change comparison between Elmer and Comsol with different phase change strategies (enthalpy, ALE) for pure silicon

Here are the SIF files:
"Phase change solid-solid"

Code: Select all

``````! Compute transient phase change solid solid with enthalpy formulation

CHECK KEYWORDS Warn
Mesh DB "./Mesh" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Max Iterations = 1
Output Intervals = 100
Timestep intervals = 2000
Timestep Sizes = 10
Timestepping Method = BDF
BDF Order = 2
Solver Input File = case.sif
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "Body 1"
Equation = 1
Material = 1
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "Body 2"
Equation = 1
Material = 1
Initial condition = 1
End

Solver 1
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 2
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Single Precision = True
Exec Solver = After saving
End

Equation 1
Name = "Equation 1"
Phase Change Model = Spatial 2
Check Latent Heat Release = True
Active Solvers(2) = 1 2
End

Material 1
Name = "Material 1"
Enthalpy = Variable Temperature; Real; 298.15 0.00E+000;350 8.89E+007;401.85 1.82E+008;453.7 2.79E+008;505.55 3.79E+008;557.4 4.80E+008;609.25 5.84E+008;661.1 6.88E+008;712.95 7.95E+008;764.8 9.02E+008;816.65 1.01E+009;868.5 1.12E+009;920.35 1.23E+009;
972.2 1.34E+009;1024.05 1.46E+009;1075.9 1.57E+009;1127.75 1.69E+009;1179.6 1.80E+009;1231.45 1.92E+009;1283.3 2.04E+009;1335.15 2.16E+009;1387 2.28E+009;1438.85 2.40E+009;1490.7 2.52E+009;1542.55 2.64E+009;1594.4 2.77E+009;
1605 2.79E+009;1615.6 2.82E+009;1626.2 2.84E+009;1636.8 2.87E+009;1647.4 2.89E+009;1658 2.92E+009;1668.6 2.94E+009;1679.2 2.97E+009;1680.2 2.97E+009;1696 7.34E+009;1750 7.47E+009;1804 7.59E+009;1858 7.71E+009;1912 7.83E+009;1966 7.95E+009;2020 8.07E+009;2074 8.20E+009;2128 8.32E+009;2182 8.44E+009;2236 8.56E+009;2290 8.68E+009;
2344 8.81E+009;2398 8.93E+009;2452 9.05E+009;2506 9.17E+009;2560 9.29E+009;2614 9.41E+009;2668 9.54E+009;2722 9.66E+009;2776 9.78E+009;2830 9.90E+009;2884 1.00E+010;2938 1.01E+010;2992 1.03E+010;
End
Phase Change Intervals(2,1) = 1680.2 1696
Heat Conductivity = 20
Heat Capacity = 1000
Density = 2330
End

Initial Condition 1
Name = "InitialCondition 1"
Temperature = 300
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1750
Save Scalars = Logical True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
External Temperature = 300
Heat Transfer Coefficient = 15
Save Scalars = Logical True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
External Temperature = 300
Heat Transfer Coefficient = 15
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Heat Flux = -10000
Save Scalars = Logical True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Heat Flux = -10000
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Save Scalars = Logical True
Save Line = Logical True
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_iso"
Save Scalars = Logical True
End
[code]

"Phase change solid-liquid"
``````
! Compute transient solid liquid phase change with enthalpy formulation

CHECK KEYWORDS Warn
Mesh DB "./Mesh" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Max Iterations = 1
Output Intervals = 20
Timestep intervals = 400
Timestep Sizes = 2
Timestepping Method = BDF
BDF Order = 2
Solver Input File = case.sif
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "Body 1"
Equation = 1
Material = 2
Body Force = 1
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "Body 2"
Equation = 1
Material = 2
Body Force = 1
Initial condition = 2
End

Solver 1
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 2
Equation = Navier-Stokes
Variable = Flow Solution[Velocity:2 Pressure:1]
Procedure = "FlowSolve" "FlowSolver"
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 3
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 1 = Temperature loads
Scalar Field 2 = Temperature
Scalar Field 4 = Pressure
Vector Field 1 = Velocity
Single Precision = True
Exec Solver = After saving
End

Equation 1
Name = "Equation 1"
Phase Change Model = Spatial 2
Check Latent Heat Release = True
Convection = Computed
Active Solvers(3) = 1 2 3
End

Material 2
Name = "Silicium_liquide"
Reference Temperature = 293
Viscosity = Variable Temperature; Real; 293.15 1.00E+007;1680 1.00E+007;1680.2 9750000;1680.6 9500000;1681 9250000;1681.4 9000000;1681.8 8750000;1682.2 8500000;
1682.6 8250000;1683 8000000;1683.4 7750000;1683.8 7500000;1684.2 7250000;1684.6 7000000;1685 6750000;1685.4 6500000;
1685.8 6250000;1686.2 6000000;1686.6 5750000;1687 5500000;1687.4 5250000;1687.8 5000000;1688.2 4750000;1688.6 4500000;
1689 4250000;1689.4 4000000;1689.8 3750000;1690.2 3500000;1690.6 3250000;1691 3000000;1691.4 2750000;1691.8 2500000;
1692.2 2250000;1692.6 2000000;1693 1750000;1693.4 1500000;1693.8 1250000;1694.2 1000000;1694.6 750000;1695 500000;
1695.4 250000;1695.8 8.80E-004;3000 8.80E-004;
End
Heat expansion Coefficient = 1.3725e-4
Enthalpy = Variable Temperature; Real; 298.15 0.00E+000;350 8.89E+007;401.85 1.82E+008;453.7 2.79E+008;505.55 3.79E+008;557.4 4.80E+008;609.25 5.84E+008;661.1 6.88E+008;712.95 7.95E+008;764.8 9.02E+008;816.65 1.01E+009;868.5 1.12E+009;920.35 1.23E+009;
972.2 1.34E+009;1024.05 1.46E+009;1075.9 1.57E+009;1127.75 1.69E+009;1179.6 1.80E+009;1231.45 1.92E+009;1283.3 2.04E+009;1335.15 2.16E+009;1387 2.28E+009;1438.85 2.40E+009;1490.7 2.52E+009;1542.55 2.64E+009;1594.4 2.77E+009;
1605 2.79E+009;1615.6 2.82E+009;1626.2 2.84E+009;1636.8 2.87E+009;1647.4 2.89E+009;1658 2.92E+009;1668.6 2.94E+009;1679.2 2.97E+009;1680.2 2.97E+009;1696 7.34E+009;1750 7.47E+009;1804 7.59E+009;1858 7.71E+009;1912 7.83E+009;1966 7.95E+009;2020 8.07E+009;2074 8.20E+009;2128 8.32E+009;2182 8.44E+009;2236 8.56E+009;2290 8.68E+009;
2344 8.81E+009;2398 8.93E+009;2452 9.05E+009;2506 9.17E+009;2560 9.29E+009;2614 9.41E+009;2668 9.54E+009;2722 9.66E+009;2776 9.78E+009;2830 9.90E+009;2884 1.00E+010;2938 1.01E+010;2992 1.03E+010;
End
Phase Change Intervals(2,1) = 1680.2 1696
Compressibility Model = Incompressible
Reference Pressure = 0
Specific Heat Ratio = 1.4
Heat Conductivity = 20
Heat Capacity = 1000
Density = 2550
End

Body Force 1
Name = "Natural convection"
Boussinesq = True
End

Initial Condition 1
Name = "InitialCondition 1"
Velocity 2 = 0
Pressure = 0
Velocity 1 = 0
Temperature = 1715
End

Initial Condition 2
Name = "InitialCondition 2"
Velocity 2 = 0
Pressure = 0
Velocity 1 = 0
Temperature = 1655
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1715
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
External Temperature = 300
Heat Transfer Coefficient = 15
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
External Temperature = 300
Heat Transfer Coefficient = 15
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Heat Flux = -10000
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Heat Flux = -10000
Noslip wall BC = True
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Save Scalars = Logical True
Save Line = Logical True
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_iso"
Temperature = 1655
Noslip wall BC = True
Save Scalars = Logical True
End
[/code]

"Phase change solver solid-liquid"

Code: Select all

``````! Compute ALE moving mesh with phase change solid liquid

CHECK KEYWORDS Warn
Mesh DB "./Mesh2" "."
Include Path ""
Results Directory "./Results"
End

Simulation
Max Output Level = 5
Coordinate System = Cartesian
Coordinate Mapping(3) = 1 2 3
Simulation Type = Transient
Steady State Min Iterations = 1
Steady State Max Iterations = 1
Output Intervals = 100
Timestep intervals = 3000
Timestep Sizes = 0.1
Timestepping Method = BDF
BDF Order = 2
Solver Input File = case.sif
End

Constants
Gravity(4) = 0 -1 0 9.82
Stefan Boltzmann = 5.67e-08
Permittivity of Vacuum = 8.8542e-12
Boltzmann Constant = 1.3807e-23
Unit Charge = 1.602e-19
End

Body 1
Target Bodies(1) = 1
Name = "liquide"
Equation = 1
Material = 2
Body Force = 1
Initial condition = 1
End

Body 2
Target Bodies(1) = 2
Name = "solide"
Equation = 2
Material = 1
Initial condition = 2
End

Body 3
Target Bodies(1) = 4
Name = "interface"
Equation = 3
Material = 1
Initial condition = 2
End

Body Force 1
Boussinesq = Logical True
End

Solver 1
Equation = Heat Equation
Procedure = "HeatSolve" "HeatSolver"
Variable = Temperature
Exec Solver = Always
Stabilize = True
Bubbles = False
Lumped Mass Matrix = False
Optimize Bandwidth = True
Steady State Convergence Tolerance = 1.0e-5
Nonlinear System Convergence Tolerance = 1.0e-7
Nonlinear System Max Iterations = 20
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 500
Linear System Convergence Tolerance = 1.0e-10
Linear System Preconditioning = ILU0
Linear System ILUT Tolerance = 1.0e-3
Linear System Abort Not Converged = False
Linear System Residual Output = 1
Linear System Precondition Recompute = 1
End

Solver 2
Equation = "Mesh Update"
Linear System Solver = "Iterative"
Linear System Iterative Method = "BiCGStab"
Linear System Convergence Tolerance = 1.0e-8
Linear System Max Iterations = 300
Steady State Convergence Tolerance = 1.0e-4
Stabilize = Logical True
Linear System Preconditioning = ILU1
End

Solver 3
Equation = string "Phase Surface"
Procedure = File "TransientPhaseChange" "TransientPhaseChange"
Variable = string PhaseSurface
Variable DOFs = Integer 1
Use nodal loads = logical true
Linear System Solver = Direct
Linear System Direct Method = Banded
!Velocity Smoothing Factor = real 0.0000001 ! the velocity diffusion factor of the interface, numerical reason, in order not to oscilate
!Velocity Relaxation Factor = real 1.0
!Surface Smoothing Factor = Real 0.0001
Steady State Convergence Tolerance = 1.0e-3
Stabilize = logical true
solver test = logical false
End

Solver 4
Equation = "Navier-Stockes"
Variable = Flow Solution[Velocity:2 Pressure:1]
Procedure = "FlowSolve" "FlowSolver"
Exec Solver = Always
Steady State Convergence Tolerance = 1.0e-5
Stabilize = Logical True
Nonlinear System Convergence Tolerance = 1.0e-6
Nonlinear System Max Iterations = 1
Nonlinear System Newton After Iterations = 3
Nonlinear System Newton After Tolerance = 1.0e-3
Nonlinear System Relaxation Factor = 1
Linear System Solver = Iterative
Linear System Iterative Method = BiCGStab
Linear System Max Iterations = 200
Linear System Convergence Tolerance = 1.0e-8
Linear System Preconditioning = ILU5
End

Solver 5
Equation = Result Output
Save Geometry Ids = True
Output Format = Vtu
Procedure = "ResultOutputSolve" "ResultOutputSolver"
Output File Name = comp_elmer_comsol
Scalar Field 2 = Temperature
Scalar Field 4 = PhaseSurface
Scalar Field 5 = Pressure
Vector Field 2 = Velocity
Single Precision = True
Exec Solver = After saving
End

Equation 1
Name = "Equation 1"
Active Solvers(4) = 1 2 4 5
Convection = "Computed"
Heat Equation = True
Mesh Update = Logical True
End

Equation 2
Name = "Equation 2"
Active Solvers(3) = 1 2 5
!Convection = "Constant"
Heat Equation = True
Mesh Update = Logical True
End

Equation 3
Name = "Equation 2"
Active Solvers(1) = 3
!Convection = "None"
End

Material 1
Name = "Solid"
Solid = Logical True
Melting Point = Real 1685
Youngs Modulus = 150000000000
Poisson Ratio = 0.27
Heat Conductivity = 20.0
Heat Capacity = 1000.0
Density = 2330.0
Reference Temperature = 1685
Heat Expansion Coefficient = 1.3725e-6
Latent Heat = Real 1800e3
Convection Velocity 1 = 0.0
Convection Velocity 2 = 0.0
End

Material 2
Name = "Liquid"
Liquid = Logical True
Melting Point = Real 1685
Youngs Modulus = 150000000000
Poisson Ratio = 0.27
Heat Conductivity = 20.0
Heat Capacity = 1000.0
Density = 2330.0
Viscosity = 8.8e-4
Reference Temperature = 1685
Heat Expansion Coefficient = 1.3725e-4
End

Initial Condition 1
Name = "InitialCondition 1"
Temperature = Variable Coordinate 1
Real
0.0 1745
0.1 1645
End
Velocity 1 = 0.0
Velocity 2 = 0.0
Pressure = 0.0
PhaseSurface = Real 0.0
Mesh Update 1 = Real 0.0
Mesh Update 2 = Real 0.0
End

Initial Condition 2
Name = "InitialCondition 2"
Temperature = Variable Coordinate 1
Real
0.0 1745
0.1 1645
End
PhaseSurface = Real 0.0
Mesh Update 1 = Real 0.0
Mesh Update 2 = Real 0.0
End

Boundary Condition 1
Target Boundaries(1) = 1
Name = "bound_hot"
Temperature = 1745
Mesh Update 1 = 0
Mesh Update 2 = 0
Save Scalars = Logical True
Noslip wall BC = True
End

Boundary Condition 2
Target Boundaries(1) = 3
Name = "bound_up3"
Mesh Update 2 = 0
!Phase Change Side = Logical True
Save Scalars = Logical True
Noslip wall BC = True
End

Boundary Condition 3
Target Boundaries(1) = 6
Name = "bound_up6"
Mesh Update 2 = 0
!Phase Change Side = Logical True
Save Scalars = Logical True
End

Boundary Condition 4
Target Boundaries(1) = 2
Name = "bound_bottom2"
Mesh Update 2 = 0
!Phase Change Side = Logical True
Save Scalars = Logical True
Noslip wall BC = True
End

Boundary Condition 5
Target Boundaries(1) = 5
Name = "bound_bottom5"
Mesh Update 2 = 0
!Phase Change Side = Logical True
Save Scalars = Logical True
End

Boundary Condition 6
Target Boundaries(1) = 4
Name = "bound_interne"
Temperature = Real 1685.0
Mesh Update 1 = Equals PhaseSurface
Mesh Update 2 = 0
Normal Target Body = Integer 1
Save Scalars = Logical True
Save Line = Logical True
Body Id = Integer 3
Noslip wall BC = True
End

Boundary Condition 7
Target Boundaries(1) = 7
Name = "bound_cold"
Temperature = 1645
Mesh Update 1 = 0
Mesh Update 2 = 0
Save Scalars = Logical True
End

``````
Julien
Attachments
Elmer_Comsol_phasechange_summary.pdf