Hi Gurus,
This is a continuation of my attempts to model a bit of silicon being
heated by a pulse of laser.
I am also including the top metal layer as a change in diffusion coef.
Modelling the heating by a FixedFlux along the right hand side, and would
like it to be non-zero before time T, and zero after time T.
Eventually I would include some convection/radiation terms along
the other 3 sides, but that is not included yet (doing this bit by bit).
I am having trouble making the flux go to zero. The signum function is
working and printing out OK, but it seems that when a BC has a time
dependent term, something goes wrong.
Here is the code:
#!/usr/bin/env python
from fipy import *
def signum ():
""" This returns 1 or 0 depending on time """
if 3.0 > t:
return 1.0
else:
return 0.0
## my mesh and subdivisions and time variable
Lx = 10.0
Ly = 10.0
Nx = 20.0
Ny = 20.0
Dx = Lx/Nx
Dy = Ly/Ny
t = Variable(value=0.0)
## create mesh
m2D = Grid2D(nx=Nx, ny=Ny, dx=Dx, dy=Dy)
x, y = m2D.getFaceCenters()
## define diffusion coefficients and different material along top edge of
box
nominalDiffusionCoef = 0.1
D = FaceVariable(mesh=m2D, value=nominalDiffusionCoef)
D.setValue(3 * nominalDiffusionCoef, where = (y > 9))
## define sources - not using at present
## boundarySource = FaceVariable(mesh=m2D, value = 0.0, rank=1)
## boundarySource.setValue(-1, where = m2D.getExteriorFaces() )
## other things I tried
## eq = TransientTerm() == ImplicitDiffusionTerm(coeff=D) + cv * t
## eq = TransientTerm() == ImplicitDiffusionTerm(coeff=D) +
boundarySource.getDivergence()
# Use the simple one
eq = TransientTerm() == ImplicitDiffusionTerm(coeff=D)
## set up boundary conditions
facesHeaters = (m2D.getFacesRight())
# want flux to go to zero after a set time
BCs = ( FixedFlux(faces=facesHeaters, value = -0.5 * signum() ))
phi = CellVariable(name = "solution variable",
mesh = m2D,
value = 1.0)
dt = 0.9 * Dx*Dy / (2 * 1.0)
## setup viewer and end-time
totalTime=10
viewer=Viewer(vars=phi)
viewer.plot()
while t() < totalTime:
eq.solve(var=phi, boundaryConditions=BCs, dt=dt)
viewer.plot()
t.setValue(t()+dt)
print signum(), t ## check
raw_input("Press <return> to proceed...")
Could you help?
Thanks a lot
--
David Wende
home +972-8-9353488
work +972-2-5886116
mobile +972-54-234-6479