Dear FiPy Authors, Could anyone of you answer my questions in the attached QFiPy060407.pdf?
Thanks and regards, Zhiheng Huang
print "Constants definition:-\n"
numberOfCells=500
Length=numberOfCells*2.5/100.
nx=numberOfCells
ny=numberOfCells
dx=Length/nx
dy=Length/ny
radius=dx*5.
seedCenter=(0.,0.)
initialTemperature=-0.5
print "Mesh generation:-\n"
from fipy.meshes.grid2D import Grid2D
mesh=Grid2D(dx,dy,nx,ny)
print "Another set of constants:-\n"
timeStepDuration=5e-5
tau=3e-4
alpha=0.015
c=0.02
N=4.
kappa1=0.9
kappa2=20.
tempDiffusionCoeff=2.25
theta=0.
print "Cell variables assembly:-\n"
from fipy.variables.cellVariable import CellVariable
phase=CellVariable(
name='phaseField',
mesh=mesh,
value=0.,
hasOld=1)
interiorCells=mesh.getCells(filter=lambda cell: \
(cell.getCenter()[0]-seedCenter[0])**2+ \
(cell.getCenter()[1]-seedCenter[1])**2<radius**2)
phase.setValue=(1.,interiorCells)
temperature=CellVariable(
name='temperature',
mesh=mesh,
value=initialTemperature,
hasOld=1)
from fipy.tools import numerix
mVar=phase-0.5-kappa1/numerix.pi*\
numerix.arctan(kappa2*temperature)
dPhiy=phase.getFaceGrad().dot((0,1))
dPhix=phase.getFaceGrad().dot((1,0))
arc=N*numerix.arctan2(dPhiy,dPhix)+theta
Phi=numerix.tan(arc/2.)
PhiSq=Phi**2
beta=(1.-PhiSq)/(1+PhiSq)
dbdpsi=-N*2*Phi/(1.+PhiSq)
A=alpha**2*c*(1.+c*beta)*dbdpsi
D=alpha**2*(1.+c*beta)**2
dxi=phase.getFaceGrad()._take((1,0),axis=1)*(-1,1)
anisotropySource=(A*dxi).getDivergence()
print "Construction order parameter PDE:-\n"
from fipy.terms.transientTerm import TransientTerm
from fipy.terms.explicitDiffusionTerm import ExplicitDiffusionTerm
from fipy.terms.implicitSourceTerm import ImplicitSourceTerm
phaseEq=TransientTerm(tau)==ExplicitDiffusionTerm(D)+\
ImplicitSourceTerm(mVar*((mVar<0.)-phase))+\
((mVar>0.)*mVar*phase+anisotropySource)
print "Construction temperature PDE:-\n"
from fipy.terms.implicitDiffusionTerm import ImplicitDiffusionTerm
temperatureEq=TransientTerm()==\
ImplicitDiffusionTerm(tempDiffusionCoeff)+\
(phase-phase.getOld())/timeStepDuration
print "Creat views:-\n"
if __name__ == '__main__':
import fipy.viewers
phaseViewer=fipy.viewers.make(vars=phase)
temperatureViewer=fipy.viewers.make(vars=temperature,
limits={'datamin':-0.5,'datamax':0.5})
phaseViewer.plot()
temperatureViewer.plot()
steps=10000
for i in range(steps):
print i
phase.updateOld()
temperature.updateOld()
phaseEq.solve(phase,dt=timeStepDuration)
temperatureEq.solve(temperature,dt=timeStepDuration)
if i%10 == 0 and __name__ == '__main__':
phaseViewer.plot()
temperatureViewer.plot()
raw_input('Press any key to continue:-\n')
QFiPy060407.pdf
Description: Adobe PDF document
