Heat source for microwave heating

Numerical methods and mathematical models of Elmer
elelel
Posts: 32
Joined: 23 May 2022, 17:37
Antispam: Yes

Re: Heat source for microwave heating

Post by elelel »

Hi Peter,
I feel honored to upload this case.
What should I do? Do I need to follow you in the github?
elelel
Posts: 32
Joined: 23 May 2022, 17:37
Antispam: Yes

Re: Heat source for microwave heating

Post by elelel »

This case has heated a potato in a straight waveguide.In the sif, the microwave heating power P can be set. P is one of the factors influencing the magnitude of the amplitude Em values. Em determines the speed of heating.
In Solver1,2, set Exec Solver = Before Simulation(Peter's advice) can greatly speed up the calculation speed.
This case has two heat sources, one heat source is the reaction heat, assumed to stay open, one is the microwave heating, on / off according to the temperature of a point.In Solver 3 ,there is a really cool thing,Apply Explicit Control.Its details can be found in this place:https://github.com/ElmerCSC/elmerfem/tr ... olExplicit It determines the success of this case
This is the sif.

Code: Select all

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

Simulation
  Max Output Level = 5
  Coordinate System = Cartesian
  Coordinate Mapping(3) = 1 2 3
  Simulation Type = Transient
  Steady State Max Iterations = 20
  Output Intervals = 1
  Timestepping Method = BDF
  BDF Order = 2
  Timestep intervals = 20
  Timestep Sizes = 1
  Solver Input File = case.sif
  Post File = case.vtu
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 = 2
  Body Force = 1
  Initial condition = 1
End


Solver 1
  Equation = Vector Helmholtz Equation
  Variable = E[E re:1 E im:1]
  Procedure = "VectorHelmholtz" "VectorHelmholtzSolver"
! seems to converge ok without
!  Linear System Preconditioning Damp Coefficient im = 1.0
  Exec Solver = Before Simulation

  Nonlinear System Max Iterations = 1
  Linear System Solver = Iterative
  Linear System Iterative Method = BiCGStabl
  Linear System Max Iterations = 1000
  Linear System Convergence Tolerance = 1.0e-7
  BiCGstabl polynomial degree = 6
  Linear System Preconditioning = vanka
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 20
  Linear System Precondition Recompute = 1
End


Solver 2
  Equation = Vector Helmholtz Post Process
  Field Variable = Electric field E[Electric field re E:3 Electric field im E:3]
  Calculate Electric Field = True
  Calculate Magnetic Field Strength = True
  Calculate Energy Functional = True
  Calculate Poynting Vector = True
  Calculate Div of Poynting Vector = True

  Calculate Nodal Fields = False
  Calculate Elemental Fields = True

  Procedure = "VectorHelmholtz" "VectorHelmholtzCalcFields"
  Exec Solver = Before Simulation

  Nonlinear System Max Iterations = 1

  Linear System Solver = Iterative
  Linear System Iterative Method = BiCGStabl
  Linear System Max Iterations = 500
  Linear System Convergence Tolerance = 1.0e-10
  BiCGstabl polynomial degree = 4
  Linear System Preconditioning = ILU0
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 20
  Linear System Precondition Recompute = 1
End

Solver 3
  Equation = Heat Equation
  Procedure = "HeatSolve" "HeatSolver"
  Variable = Temperature
  Exec Solver = Always
  Stabilize = True
  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
  Nonlinear System Consistent Norm = True

  Linear System Solver = Iterative
  Linear System Iterative Method = BiCGStab
  Linear System Max Iterations = 500
  Linear System Convergence Tolerance = 1.0e-10
  BiCGstabl polynomial degree = 6
  Linear System Preconditioning = ILU0
  Linear System ILUT Tolerance = 1.0e-3
  Linear System Abort Not Converged = False
  Linear System Residual Output = 10
  Linear System Precondition Recompute = 1

  Apply Explicit Control = Logical True   
  Number Of Controls = Integer 1    
  Control Node Coordinates(1,3) = Real 0.04 0.0225 0.102
 ! Control Variable = string "Temperature"
End

Solver 4
  Equation = SaveScalars
  Procedure = "SaveData" "SaveScalars"
  Parallel Reduce = Logical True  ! Save the control value   
  Variable 1 = cpar   
  Variable 2 = Temperature   
  Operator 2 = min    
  Variable 3 = Temperature   
  Operator 3 = max    
  Filename = f.dat   
  Show Norm Index = 1
End

Equation 1
  Name = "Equation 1"
$ a = 0.086 
$ b = 0.043
$ P = 10000
$ mu0 = 4*pi*10^-7
$ mur = 1
$ epsilon0 = 8.854e-12
$ epsilonr = 1
$ f = 2.45e9
$ vp0 = 3e8
!*************Em*********************
$ mu = mu0*mur
$ epsilon = epsilon0*epsilonr
$ fc10 = 1/(2*a*sqrt(mu*epsilon))
$ vp10 = vp0/sqrt(1-(fc10/f)^2)
$ eta = sqrt(mu/epsilon)
$ ZTE10 = eta/sqrt(1-(fc10/f)^2)
$ Em = sqrt(P*4*ZTE10/(a*b))
!************beta0*************************
$ c0 = 1/sqrt(mu0*epsilon0) 
$ omega=2*pi*f 
$ k0 = omega/c0 
$ kc = pi/a 
$ beta0 = sqrt(k0^2-kc^2)
  Angular Frequency = $2*pi*2.45e9
  Active Solvers(4) = 1 2 3 4 
End

Body Force 1
  Name = "Body Force 1"
  Heat Source = Variable "cpar, Div Poynting Vector re e"
   Real Procedure "AlternatingSource" "AlternatingSource"
End

Material 1
  Name = "Air (room temperature)"
  Relative Permittivity = 1.00059
  Sound speed = 343.0
  Viscosity = 1.983e-5
  Density = 1.205
  Heat Capacity = 0
  Heat Conductivity = 0
  Heat expansion Coefficient = 0
End

Material 2
  Name = "tudou"
  Relative Permittivity im = real 20
  Relative Permittivity = real 65
  Heat Conductivity = 0.55
  Relative Reluctivity = real 1
  Density = 1050
  Heat Capacity = 3.64e3
  electric conductivity = real 0
End


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

Boundary Condition 1
  Target Boundaries(1) = 3 
  Name = "Inport"
  Electric Robin Coefficient im = $ beta0
  Magnetic Boundary Load im 2 = Variable Coordinate 1 
    Real MATC "-2*beta0*Em*sin(kc*(tx))"
End

Boundary Condition 2
  Target Boundaries(1) = 4 
  Name = "Outport"
  Electric Robin Coefficient im = $ beta0
End

Boundary Condition 3
  Target Boundaries(4) = 1 2 5 12 
  Name = "PEC"
  E re {e} = 0
  E im {e} = 0
End
This is the udf

Code: Select all

FUNCTION AlternatingSource(Model, n, x) RESULT(ht)
  USE DefUtils
  
  IMPLICIT None
  
  TYPE(Model_t) :: Model
  INTEGER :: n
  REAL(KIND=dp) :: x(2),ht
  
  REAL(KIND=dp) :: hvar,tlimit,hcon,time,tpoint,Poy

  tpoint = x(1)
  Poy = x(2)
  
 
  hvar = abs(Poy)/1050  
  hcon = -4761.9
  tlimit = 600.0
  ht = hcon + hvar
  if (tpoint <= tlimit) then
  ht = hcon + hvar
  else
  ht = hcon 
  end if 
    
END FUNCTION AlternatingSource
This is the comparison of the results from Elmer and Comsol. In Comsol ,the default potato boundary is thermal insulation. I don't know how Elmer does thermal insulation, so set the heat transfer coefficient,Heat Conductivity and Heat Capacity of air to 0.So Elmer's potato boundary looks a little weird.The temperature is a little different, maybe the control point temperature or the reaction heat is a little different.
comsol0-20sec.png
comsol0-20sec.png (116.72 KiB) Viewed 502 times
Elmer0-20sec.png
Elmer0-20sec.png (109.4 KiB) Viewed 502 times
0.001m_Mesh.7z
(574.18 KiB) Downloaded 50 times
elelel
Posts: 32
Joined: 23 May 2022, 17:37
Antispam: Yes

Re: Heat source for microwave heating

Post by elelel »

The upper microwave power P is 10kw, probably heating too fast, so the result is a little weird.It should not use so much power in the real world.Our two-meter-long reaction was planned with 5kw,this waveguide is far too small.So I simulated P=1kw.
Conditions:P=1kw,Max mesh of potato=0.0004m,Reaction heat is expressed by the heat consumption rate,Q0=P0/V=-10W/(0.01m*0.01m*0.02m)=-5e6W/m3=-4761.9W/kg,Tlimit=400K
This is the results.
comsol0-20sec.png
comsol0-20sec.png (149.31 KiB) Viewed 484 times
Elmer0-20sec.png
Elmer0-20sec.png (127.71 KiB) Viewed 484 times
This is tpoint at different time steps
tpoint.png
tpoint.png (67.34 KiB) Viewed 470 times
Then I'll try a finer mesh, and now I'm using a virtual machine, it can't count a finer grid. If I can get better results, I will update them here.
Post Reply