Transient Phase Change Solver

Numerical methods and mathematical models of Elmer
Post Reply
Roland
Posts: 226
Joined: 12 Apr 2018, 11:29
Antispam: Yes

Transient Phase Change Solver

Post by Roland »

Hi,
I would like to use the Transient Phase Change Solver with the associated ALE front tracking method for settig up a cristal growth process in an induction furnace. In a former discussion on the forum I saw the very interesting work done by Julien Givernau (I don't know if you are still on the forum?), especially a Comsol/Elmer comparison work of the phase change in a closed box (see the attached plot).
For this the Transient Phase Change Solver is used (in addition with the convective Heat_Transfer, Navier_Stokes and Mesh_Update solvers).
I would like to better understand the use of this "front tracking" solver.
In Julien's corresponding case.sif file, the BC on the phase change moving front is following:
*************************************************************
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
*************************************************************
My questions are:
1/What is exactly the PhaseSurface variable? Is it the x coordinate of the phase surface, or both x and y coordinates of each interface node?
2/ Why is the Mesh Update 2 variable set to zero?
3/ what means the command line "Normal Target Body = Integer 1" ?
4/ How does the Transient Phase Change Solver know what side of the interface is the liquid phase and what side is the solid phase ?
5/ Must the Latent Heat be given in both the solid and the liquid materials or only in one material?
6/ Are there some such examples (or tutorials) of phase change ( in static and also in cristal pulling situation) anywhere in Elmer?

Thanks a lot in advance for your help about this since, after the different kinds of models using magnetic fields which I have setup these last months (like cold crucible, EM_levitation, with your precious help, Peter, Kevin, Rich...!) we would like to investigate the phase change physics for develloping cristal growth processes by magnetic fields with Elmer in the next future.
Roland
Attachments
Front_Tracking_Phase_Change.PNG
(375.39 KiB) Not downloaded yet
kevinarden
Posts: 2237
Joined: 25 Jan 2019, 01:28
Antispam: Yes

Re: Transient Phase Change Solver

Post by kevinarden »

More discussion here
viewtopic.php?t=4138
Roland
Posts: 226
Joined: 12 Apr 2018, 11:29
Antispam: Yes

Re: Transient Phase Change Solver

Post by Roland »

Hi Kevin,
Thank you very much, I will take a look at this.
Roland
Posts: 226
Joined: 12 Apr 2018, 11:29
Antispam: Yes

Re: Transient Phase Change Solver

Post by Roland »

Hi Kevin again,
I just took a quick look at your mentioned discussion, but I found nothing about the Lagrangian ALE front tracking method(done by the Transient Phase Change Solver which I would understand better...). There are only discussions about the enthalpy or equivalent heat capacity methods in which I am less interested.
Thanks to tell me what you think about that.
Roland
raback
Site Admin
Posts: 4812
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: Transient Phase Change Solver

Post by raback »

Hi Roland,

It is assumed that the phase surface is associated to some coordinate direction. Like in your case it is the x-coordinate and y-coordinate is then frozen. The PhaseSurface is the displacement of this active coordinate at the free surface only. Hence the MeshSolver is needed to extend the displacement to the fluid in some continuous way. It transient cases this also gives rise to Mesh Velocity which is automatically taken into account in the FlowSolver using the ALE formulation.

The initial use of this phase change solver was for crystal growth type of phenomena where you have a nicely behaving free surface. This approach fails when the free surface cannot be projected to plane.

-Peter
Roland
Posts: 226
Joined: 12 Apr 2018, 11:29
Antispam: Yes

Re: Transient Phase Change Solver

Post by Roland »

Hi Peter,
Thanks for your reply and your explanations,even if I don't yet understand all the subtilities of this very interesting phase change solver. I have some additional questions for which your help would be very precious:
1/ could you please explain the comand line "Normal target Body = integer 1" in the phase change interface of Julien Givernau's Sif file(my 3rd former question)?
2/what do you mean by "This approach fails when the free surface cannot be projected to plane"?
3/how is it possible to add a BC on the phase change interface which simulates the latent heat release on the heat transfer physics?
4/are there some examples of crystal growth with ALE front tracking method in Elmer? Is it possible to include the crystal pulling in the Transient Phase Change Solver?
Thanks in advance for your help and we keep us informed
Roland
raback
Site Admin
Posts: 4812
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: Transient Phase Change Solver

Post by raback »

Hi

1) I think you need to have consistent normals. At internal interface this is ambiguous. I think this is the reason, I'm not sure though...
2) As the free surface in 3D is h(x,y) there needs to be a unique value for each (x,y) point.
3) Add Latent Heat + velocity
4) elmeruser@elmeruser-VirtualBox:~/elmerfem/fem/tests$ grep TransientPhaseChange -r
PhaseChange3/PhaseChangeTrans.sif: Procedure = "TransientPhaseChange" "TransientPhaseChange"

-Peter
Roland
Posts: 226
Joined: 12 Apr 2018, 11:29
Antispam: Yes

Re: Transient Phase Change Solver

Post by Roland »

Hi Peter,
Thank you very much for these details.
I will try to setup a first simplified crystal growth model. Just one thing concerning my last 3rd question: in order to add the latent heat on the phase change interface, which depends on the Stefan velocity, how is it possible to get this Stefan velocity?
Thanks for you help and we keep us informed.
Roland
raback
Site Admin
Posts: 4812
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: Transient Phase Change Solver

Post by raback »

Hi Roland,

In transient case you set the interface to be at melt point, compute the heat fluxes, and from the timestep can deduce the speed of the melt interface. For crystal growth with pulling you need to in some way follow this velocity (feedback) or use some predefined pull velocity. The feedback requires some special routines. Generally it may be a difficult task because there may be different timescales in the system.

-Peter
Roland
Posts: 226
Joined: 12 Apr 2018, 11:29
Antispam: Yes

Re: Transient Phase Change Solver

Post by Roland »

Hi,
After having made with the GUI the phase change model in a box (strongly inspired from the Julien's Givernau model / still many thanks Julien!) I try now to make a phase change model in a vertical crystal growth configuration(the zipped GUI Elmer project is here attached) . Here attached also is the simplified geometry: the lower body is the liquid (body 1/liquid Silicon), the upper left body (body 2) is the crystal pulled out with a "v_cast" upward velocity, and the upper right body (body 3) is the surrounding air. The crystal (body 2) is considered as a "solid equivalent" domain since it is still treated, like the lower liquid, with the Navier Stokes physics, but with an infinitely big viscosity (7e10 Pa*s!) which makes it behave like a rigid, meaning a solid, body.

The coupled physics are:
- Heat transfer equation(say HT/priority 5)
- Navier Stokes equation (say NS/priority 4)
- Transient Phase Change Solver (say PS/priority 3) on the phase change interface ("PhaseChange front") boundary between lower liquid body 1 and upper left "solid" body 2.
- Compute Normal Solver (say CN/priority 2)) on the free surfaces of both the lower liquid body 1 and the right side boundary of the upper left solid body 2 (which is the ALE method for displacing the boundary nodes of the MeshUpdate equation )
- MeshUpdate Solver (say MU/priority 1) for updating the mesh.

The "hot" temperature is fixed on the right vertical wall of the liquid body, the melting temperature is fixed on the interface boundary between the liquid and the solid, and the "cold" temperature is fixed on the upper horizontal boundary of the solid body (which is the pulled up crystal output). The thermoconvection in the liquid is setup with the Boussinesq assumption.

In order to set the priorities of the solvers in the GUI I used a trick, meaning that for setting the PS priority at 3, I activated the ResultOutput equation and replaced its procedure by "TransientPhaseChange" "TransientPhaseChange" and entered the priority 3, and same for the CN for which I activated the Electrostatics equation and replaced its procedure by "NormalSolver" "NormalSolver" and entered the priority 2.

Unfortunately the solve stops at the 1rst step and the Solver Log (also here attached) shows something like a "NaN" for the Transient Phase Change Solver. When I disable PS, then the solve works well. Same when I disable CN. It seems that PS and CN don't want to work together. Does that perhaps come from the triple point where both PS and CN apply which is perhaps conflicting for MU?

As I am now out of ideas, it would be very kind if somebody who has more experience than me could take a look at it and help me to fix this problem. Thanks in advance for any help!
Roland
Attachments
Solver Log.txt
(9.05 KiB) Downloaded 32 times
Crystal_Growth_geom.PNG
Crystal_Growth_geom.PNG (35.16 KiB) Viewed 869 times
CG_02_Elmer.zip
(110.95 KiB) Downloaded 31 times
Post Reply