Elmer Discussion Forum

Posted: **06 Mar 2015, 13:12**

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 1 = Temperature Loads
  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
  Equation = Flux and Gradient
  Calculate Flux = Logical True
  Calculate Flux Abs = Logical True
  Calculate Grad = 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"
  Calculate Loads = True
  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 1 = Temperature Loads
  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
  Equation = Flux and Gradient
  Calculate Flux = Logical True
  Calculate Flux Abs = Logical True
  Calculate Grad = 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"
  Calculate Loads = True
  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"
  Calculate Loads = True
  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

Header
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

Posted: **06 Mar 2015, 13:23**

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 1 = Temperature Loads
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
Equation = Flux and Gradient
Calculate Flux = Logical True
Calculate Flux Abs = Logical True
Calculate Grad = 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"
Calculate Loads = True
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 1 = Temperature Loads
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
Equation = Flux and Gradient
Calculate Flux = Logical True
Calculate Flux Abs = Logical True
Calculate Grad = 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"
Calculate Loads = True
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"
Calculate Loads = True
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

Header
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

Posted: **06 Mar 2015, 18:22**

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 1 = Temperature Loads
  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
  Equation = Flux and Gradient
  Calculate Flux = Logical True
  Calculate Flux Abs = Logical True
  Calculate Grad = 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"
  Calculate Loads = True
  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 1 = Temperature Loads
  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
  Equation = Flux and Gradient
  Calculate Flux = Logical True
  Calculate Flux Abs = Logical True
  Calculate Grad = 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"
  Calculate Loads = True
  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

Posted: **12 Mar 2015, 11:33**

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

Posted: **12 Mar 2015, 23:40**

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

Posted: **13 Mar 2015, 10:39**

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

Posted: **13 Mar 2015, 11:39**

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

Posted: **13 Mar 2015, 11:51**

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

Posted: **13 Mar 2015, 12:26**

Dear Julien,
Thank you very much.

Yours
Anil Kunwar

Posted: **19 Mar 2015, 13:11**

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

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
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"
  Calculate Loads = True
  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

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
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"
Calculate Loads = True
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"
Calculate Loads = True
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

Header
  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"
  Calculate Loads = True
  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"
  Calculate Loads = True
  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

Elmer Discussion Forum

phase change comparison elmer/Comsol

phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol

Re: phase change comparison elmer/Comsol