Re: Phase Field Modeling Using moose

287 views
Skip to first unread message

Roma Gurung

unread,
Jul 18, 2015, 4:33:29 PM7/18/15
to moose...@googlegroups.com
Hi Michael,
With your strong instructions, i have learnt a lot about he temperature/voltage driven shape evolution in void. 

Role of the Diffusivity
After you have mentioned  the role of Diffusivity, I tried to study about it and as D = D_o *exp(-Em/kbT), i have observed from the simulation run that , the smaller value of enthalpy of vacancy migration (Em), would (1) make the solution converge in one hand , and (2) make the evolution of surface more pronounced. Hence, greater diffusivity ensures the prominent shape evolution.

Length Scale and Time Scale
Now, i am trying to understand the role of length scales and time scales in the output. For a geometry in microscale, i first tried to define _length_scale = 1.0e-6 )(um), the solution did not converge. Then i again put _length_scale = 1.0e-9 (nm) whereas defined mesh geometry 10^3 times bigger which eventually shifted the computational time from few minutes to several hours; however this made the solution converged.
Also, for slower evolutions, i increased the time scale as _time_scale = 1.0e3. So, that the evolution after a duration of 2 days can be tracked with smaller timesteps (number of time steps). What happens to the interface if the time scale is too large? Because i noted the spots having same color of the interface (betweeen void and metal) scattered in move in the domain when i made the time scale large.

Yours Sincerely
Anil Kunwar

On Thu, Jul 9, 2015 at 4:30 AM, Tonks, Michael R <michae...@inl.gov> wrote:
If the solve is converging but the concentration appears to not evolve, then your time step is too small or your diffusivity is too small (increasing either one should have the same effect).

If the concentration is evolving but it isn't impact by the voltage, then your driving force in SplitCHVoltage may be too small.

On Tue, Jul 7, 2015 at 12:46 PM, Roma Gurung <romagu...@gmail.com> wrote:
Hi,
I would like to see the shape evolution of surface voids due to voltage gradient and i gave a try using the following input file:
[Mesh]
type = GeneratedMesh
dim = 2
nx = 60
ny = 60
xmax = 500
ymax = 500
elem_type = QUAD4
[]
[GlobalParams]
polynomial_order = 8
[]
[Variables]
[./c]
[../]
[./w]
scaling = 1.0e2
[../]
[./volt]
initial_condition = 1000
scaling = 1.0e5
[../]
[]
[ICs]
[./c_IC]
type = SmoothCircleIC
x1 = 125.0
y1 = 250.0
radius = 60.0
invalue = 1.0
outvalue = 0.1
int_width = 100.0
variable = c
[../]
[]
[Kernels]
[./electrical_potential]
type = MatDiffusion
variable = volt
D_name = electrical_conductivity
[../]
#[./electric_potential]
# type = ElectricPotential
# variable = volt
# conductivity = 73 # (W/m K) From NIST leadfree solder database
# [../]
[./c_res]
type = SplitCHParsed
variable = c
kappa_name = kappa
w = w
f_name = F
[../]
[./w_res]
type = SplitCHWRes
variable = w
mob_name = M
[../]
[./w_res_voltage]
type = SplitCHVoltage
variable = w
c = c
volt = volt
diff_name = D
z_name = zeff
T = T
[../]
[./time]
type = CoupledTimeDerivative
variable = w
v = c
[../]
[]
[BCs]
[./bottom_volt]
type = DirichletBC
variable = volt
boundary = bottom
value = 1024
[../]
[./top_volt]
type = DirichletBC
variable = volt
boundary = top
value = 1000
[../]
[]
[Materials]
[./Copper]
type = VoltPFParamsPolyFreeEnergy
block = 0
c = c
volt = volt
T = T # K
int_width = 60.0
length_scale = 1.0e-9
time_scale = 1.0e-9
D0 = 3.1e-5 # m^2/s, from Brown1980
Em = 0.71 # in eV, from Balluffi1978 Table 2
Ef = 1.28 # in eV, from Balluffi1978 Table 2
surface_energy = 0.708 # Total guess
[../]
[./elcond]
type = ParsedMaterial
block = 0
args = 'c'
function = 'if(c>0.7,1.0e-8,5.5e-8)'
f_name = electrical_conductivity
outputs = exodus
[../]
[./free_energy]
type = PolynomialFreeEnergy
block = 0
c = c
derivative_order = 2
[../]
[]
[Preconditioning]
[./SMP]
type = SMP
full = true
[../]
[]
[Executioner]
type = Transient
scheme = 'bdf2'
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -ksp_grmres_restart -sub_ksp_type -sub_pc_type -pc_asm_overlap'
petsc_options_value = 'asm 31 preonly lu 1'
l_max_its = 30
l_tol = 1.0e-4
nl_max_its = 25
nl_rel_tol = 1.0e-9
num_steps = 60
dt = 20.0
[]
[Outputs]
output_initial = true
interval = 1
exodus = true
print_perf_log = true
[]
Where, VoltPFParamsPolyFreeEnergy is the modified Class of  PFParamsPolyFreeEnergy and
SplitCHVoltage is the modified class for SoretDiffusion.
However, the shape evolution of void is not seen pronounced. What are the numerical parameters that need to be checked to observe it:
1) surface diffusivity (D_s)
2) kappa ?


Yours Sincerely
Anil Kunwar

On Mon, Jul 6, 2015 at 1:50 AM, Roma Gurung <romagu...@gmail.com> wrote:
Hi Michael,
Looking at your input file, i noted that i had not put "Initial Condition" on temperature variable. Putting it , the solver seem to be more converged than before.
Moreover, i need to try with the length scales, time scales and other suggestions given by you to ensure the convergence of the solve.I am happy to learn many concepts about MOOSE framework.

Yours Sincerely
Anil Kunwar

On Sun, Jul 5, 2015 at 6:45 AM, Tonks, Michael R <michae...@inl.gov> wrote:
Roma,

I think that your best bet would be to start with phase_field/tests/SoretDiffusion/split_temp_test.i and modify it. I took that test and made it 2D in the attached input file.

Other things you should think about:

- The units are set in the material block, and the units in PFParamsPolyFreeEnergy are Kelvin for temperature, "length_scale" for length and "time_scale" for time. So, you need to make sure that your thermal conductivity is in the same units (in your input file you used W/mK for the thermal conductivity even though length_scale = 1e-9 (nm) and time_scale = 1e-6 (microsec).
- If you decrease the temperature, the diffusivity goes down exponentially, so you will need to increase the time scale or evolution will occur VERY slowly
- Your temperature boundary conditions need to be in the same units as in your material (K in this case)



On Sat, Jul 4, 2015 at 12:09 PM, Roma Gurung <romagu...@gmail.com> wrote:
Hi,
As an exercise, i have tried to learn SoretDiffusion tutorial by coupling Heat Conduction Equation with the Cahn Hilliard Equation (Phase FIeld Equation).
The input file was modified from the SoretDiffusion test input file and has the following structure:
[Mesh]
 type = GeneratedMesh
 dim = 2
 xmax = 1000
 ymax = 1000
 nx = 50
 ny = 50
[]
[GlobalParams]
 polynomial_order = 8
[]
[Variables]
 [./T]
 [../]
 [./c]
 [../]
 [./w]
     scaling = 1.0e3
 [../]
[]
[ICs]
 [./c_IC]
    type = SmoothCircleIC
    x1 = 175.0
    y1 = 0.0
    radius = 100
    invalue = 1.0
    outvalue = 0.01
    int_width = 100.0
    variable = c
 [../]
[]
#[AuxVariables]
#[./T]
#[../]
#[]
[Kernels]
 [./heat_conduction]
    type = HeatConduction
    variable = T
 [../]
 [./c_res]
    type = SplitCHParsed
    variable = c
    kappa_name = kappa
    w = w
    f_name = F
 [../]
 [./w_res]
    type = SplitCHWRes
    variable = w
    mob_name = M
 [../]
 [./w_res_soret]
    type = SoretDiffusion
    variable = w
    c = c
    T = T
    diff_name = D
    Q_name = Qstar
 [../]
 [./time]
    type = CoupledTimeDerivative
    variable = w
    v = c
 [../]
[]
#[AuxKernels]
# [./Temp]
#    type = FunctionAux
#    variable = T
#    function = 1000.0+0.025*x
# [../]
#[]
[BCs]
  [./bottom_temperature]
    type = DirichletBC
    variable = T
    boundary = bottom
    value = 350 # (C)
  [../]
  [./top_temperature]
    type = DirichletBC
    variable = T
    boundary = top
    value = 348 # (C)
  [../]
[Materials]
 [./Copper]
    type = PFParamsPolyFreeEnergy
    block = 0
    c = c
    T = T # K
    int_width = 80.0
    length_scale = 1.0e-9
    time_scale = 1.0e-6
    D0 = 3.1e-5 # m^2/s, from Brown1980
    Em = 0.71 # in eV, from Balluffi1978 Table 2
    Ef = 1.28 # in eV, from Balluffi1978 Table 2
    surface_energy = 0.708 # Total guess
 [../]
 [./Cu_cond]
    type = GenericConstantMaterial
    block = 0
    prop_names = thermal_conductivity
    prop_values = 73 # K: (W/m*K) from wikipedia @296K
  [../]
 [./free_energy]
    type = PolynomialFreeEnergy
    block = 0
    c = c
    outputs = exodus
    derivative_order = 2
 [../]
[]
[Preconditioning]
 [./SMP]
    type = SMP
    full = true
 [../]
[]
[Executioner]
    type = Transient
    scheme = 'bdf2'
    solve_type = 'NEWTON'
    petsc_options_iname = '-pc_type -ksp_grmres_restart -sub_ksp_type -sub_pc_type -pc_asm_overlap'
    petsc_options_value = 'asm 31 preonly lu 1'
    l_max_its = 10
    l_tol = 1.0e-4
    nl_max_its = 25
    nl_rel_tol = 1.0e-9
    start_time = 0.0
    num_steps = 20
    dt = 3
[]
[Outputs]
    output_initial = true
    exodus = true
    print_perf_log = true
[]
When i tried to run the input file, the solve did not converge with the following message.
The log_file is attached herewith. What is the cause of the "convergence" issue?

Yours Sincerely
Anil Kunwar 

On Wednesday, July 1, 2015 at 11:40:57 PM UTC+8, michael.tonks wrote:
As a quick update, SplitCHSoret has changed names. We adjusted it to work for both the split and direct solve of the Cahn-HIlliard equation, so it is now call SoretDiffusion. See phase_field/tests/SoretDiffusion.

Mike

On Tue, Jun 30, 2015 at 3:49 PM, Tonks, Michael R <michae...@inl.gov> wrote:
Anil,

When you have a temperature dependent mobility and free energy, temperature gradient terms come automatically out of the Cahn-Hilliard equation. In that case, you would only need to SplitCHParsed. However, when you are using a polynomial free energy without temperature dependence, you have to artificially add the temperature gradient driving force, as we did here.

On Tue, Jun 30, 2015 at 2:03 PM, Daniel Schwen <dan...@schwen.de> wrote:
5. A newer question: Earlier Daniel had commented somewhere that as fields like Temperature (T) and Potential (V) have some free-energy contribution associated with them. Is there a possibility to derive a material class from the DerivativeFunctionMaterial
 class for the effect of temperature field and couple it with the heat conduction equation rather than use the function = A + B.T, where T is temperature in [AuxKernels] block?

You can certainly derive a new material class from DerivativeFunctionMaterialBase [1], but I strongly suggest that you try implementing the model with a DerivativeParsedMaterial first.
Daniel

--
You received this message because you are subscribed to the Google Groups "moose-users" group.
 To unsubscribe from this group and stop receiving emails from it, send an email to moose-users...@googlegroups.com.
Visit this group at http://groups.google.com/group/moose-users.
 To view this discussion on the web visit https://groups.google.com/d/msgid/moose-users/CAPxoKqcz_5God2fE0CaxL_0bN3sR1ziwnrp9_BdopzFBbTKD9g%40mail.gmail.com.

For more options, visit https://groups.google.com/d/optout.



--

*************************

Michael R. Tonks

Fuel Modeling and Simulation Department

Idaho National Laboratory

208-526-6319

michae...@inl.gov




--

*************************

Michael R. Tonks

Fuel Modeling and Simulation Department

Idaho National Laboratory

208-526-6319

michae...@inl.gov

--
You received this message because you are subscribed to the Google Groups "moose-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to moose-users...@googlegroups.com.
Visit this group at http://groups.google.com/group/moose-users.
 To view this discussion on the web visit https://groups.google.com/d/msgid/moose-users/931d66ef-6d03-440b-ae51-d298f460cdf7%40googlegroups.com.

For more options, visit https://groups.google.com/d/optout.



--

*************************

Michael R. Tonks

Fuel Modeling and Simulation Department

Idaho National Laboratory

208-526-6319

michae...@inl.gov

--
You received this message because you are subscribed to a topic in the Google Groups "moose-users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/moose-users/4qzSmJJlI_g/unsubscribe.
To unsubscribe from this group and all its topics, send an email to moose-users...@googlegroups.com.
Visit this group at http://groups.google.com/group/moose-users.
 To view this discussion on the web visit https://groups.google.com/d/msgid/moose-users/CA%2BEq24_U2WPwrTWtxvdv5OKVh8GYB%2BkDQgQ5KLDLhHOawyPE%2BQ%40mail.gmail.com.

For more options, visit https://groups.google.com/d/optout.


--
You received this message because you are subscribed to the Google Groups "moose-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to moose-users...@googlegroups.com.
Visit this group at http://groups.google.com/group/moose-users.
To view this discussion on the web visit https://groups.google.com/d/msgid/moose-users/CAFuDzTxyZ%2Bqc1xv1_tAUtRie7e34fSXs4qdySdrOF_ULFKCtHg%40mail.gmail.com.
For more options, visit https://groups.google.com/d/optout.



--

*************************

Michael R. Tonks

Fuel Modeling and Simulation Department

Idaho National Laboratory

208-526-6319

michae...@inl.gov

--
You received this message because you are subscribed to a topic in the Google Groups "moose-users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/moose-users/4qzSmJJlI_g/unsubscribe.
To unsubscribe from this group and all its topics, send an email to moose-users...@googlegroups.com.
Visit this group at http://groups.google.com/group/moose-users.
To view this discussion on the web visit https://groups.google.com/d/msgid/moose-users/CA%2BEq24_iiAHgpExoaDpOM2uxSmW0benMmEoax1Lia4t8txSirQ%40mail.gmail.com.
For more options, visit https://groups.google.com/d/optout.

john.m...@uconn.edu

unread,
Jul 19, 2015, 12:23:44 PM7/19/15
to moose...@googlegroups.com
Interesting,

What is the difference between CoupledTimeDerivative and TimeDerivative kernels?

John
...

Roma Gurung

unread,
Jul 19, 2015, 4:07:48 PM7/19/15
to moose...@googlegroups.com
Hi John,
In principle they are same.
TimeDerivative kernel represents dc/dt where c is the variable of the FEM simulation.
CoupledTimeDerivative is also the time derivative kernel utilized mostly in coupled form of two residual equations. Previously, it was given the Kernel name CoupledImplicitEuler

For example,
(Please refer the corresponding tables in the above website.)

1. The  residual for the direct solution of Cahn-Hilliard Equation and the Residual for the Allen-Cahn equation uses "TimeDerivative" Kernel.

2. The residual for the split form of Cahn-Hilliard Equation uses the "CoupledTImeDerivative" Kernel.


Yours Sincerely
Anil Kunwar





...

john.m...@uconn.edu

unread,
Jul 19, 2015, 4:25:54 PM7/19/15
to moose...@googlegroups.com
Oh ok.

I have a coupled Cahn-Hilliard-like system that we're solving and I have gotten some moderate success just using TimeDerivative so was hoping that wasn't something obvious that I'm missing.

Also perhaps related to your post: sometimes if you are having troubles with the nonlinear convergence using the length_scale method to scale your kernels to the mesh size, using a certain unit system encourages precision loss and thus makes it hard for the nonlinear solve. I had a problem where we were using meters for coefficients in the input file (LGD type expansion of the order parameter), length_scale = 1e-9, and it was not converging; but when we went to an entirely different unit system where coefficients fed to the kernels were around order 1, the problems with the nonlinear convergence disappeared.

Cheers,
John
...

Derek Gaston

unread,
Jul 20, 2015, 10:00:21 AM7/20/15
to moose...@googlegroups.com
BTW: If you go look at the code it's pretty easy to see the difference between the two TimeDerivative classes (pay no attention to the Jacobian stuff... just look at the residuals and the coupling):



TimeDerivative is just the normal du/dt for _this_ variable.  CoupledTimeDerivative pulls the time derivative of another variable into this variable's equation.

So if you have two equations:

du/dt + dv/dt = f
dv/dt = g

In the first equation you would use TimeDerivative for du/dt and CoupledTimeDerivative for dv/dt.

Derek



--
You received this message because you are subscribed to the Google Groups "moose-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to moose-users...@googlegroups.com.
Visit this group at http://groups.google.com/group/moose-users.
To view this discussion on the web visit https://groups.google.com/d/msgid/moose-users/6b2f139d-9dd6-440b-b13d-d3585da706be%40googlegroups.com.

Tonks, Michael R

unread,
Jul 20, 2015, 10:55:17 AM7/20/15
to moose...@googlegroups.com
Roma,

When you change the length and time scales, you are really just changing the magnitude of the Mobility terms in the equations (the mobilities have units involving time and length). I recommend that you look at the mobility material properties as you change the unit scales and see what is really going on. The easiest way of doing this is to add the line

outputs = exodus

in the materials block that is calculating the mobility. Then, you should be able to see the mobility value in the exodus file that is outputted each time step. Note, it will be a element value, not a nodal value. 


For more options, visit https://groups.google.com/d/optout.

Tonks, Michael R

unread,
Jul 20, 2015, 11:03:40 AM7/20/15
to moose...@googlegroups.com
I would like to clarify a little what Roma said about the two time derivative kernels.

Both add a time derivative of a variable. However, remember that each kernel makes up a piece of the residual equation for one variable, denoted in the input file with the line

variable = X

The TimeDerivative kernel adds a time derivative of that variable to the residual equation.

However, often we will be solving for multiple variables simultaneously that are all coupled. So, if the residual for variable X has a term using the time derivative of variable Y, we use the CoupledTimeDerivative. It takes two inputs, 

variable = X denotes which residual equation this term is added to (but is not used for the time derivative) and

v = Y denotes the coupled variable (not the one the residual is solving for) for which you need the time derivative.
Reply all
Reply to author
Forward
0 new messages