On May 30, 2012, at 3:56 PM, Kendall Boniface wrote:
> I'm looking for a bit of assistance in how to organize a Fipy script I am
> trying to write. I am very new to programming so please forgive me if this
> becomes overly evident as I go on.
No need to apologize. Your script was pretty close.
> I am trying to model the temperature profile in a thermohydraulic loop.
> Specifically, right now I would just like to look at the 1D transient
> conduction occurring through Inconel 600 (what the pipe is made of).
> Temperature depends on the material's thermal conductivity (obviously), but
> the thermal conductivity of Inconel 600 in turn changes slightly with
> temperature. (The specific equation is: k = (k0) * e^(omega*temp) ), where
> k0 and omega are given.
>
> I am fairly certain that interations, sweeps and time steps are required but
> I am having trouble organizing it properly to that effect.
Actually, your iteration/sweep/step loops were nearly right. I've modified your
script below, with remarks from me starting with "JEG:". Please ask if the
things I've done aren't clear.
from fipy import *
L = 1.
nx = 50.
dx = L / nx
mesh = Grid1D(nx=nx, dx=dx)
# Some definitions
omega = 0.00164
k0 = 11.83
# Create cell variable objects
temp = CellVariable(mesh=mesh, name='Temperature', hasOld=1)
# JEG: Thermal conductivity is a function of temperature, so we write it that
way
# FiPy will take care of updating it automatically as temperature changes
therm_cond = k0 * numerix.exp(omega * temp)
therm_cond.name = 'Thermal Conductivity'
# JEG: *never* use the math package with FiPy.
# Always use fipy.numerix, instead.
# The functions of numpy are needed to deal with arrays of values,
# and FiPy augments numpy in many ways (such as being able to write
# lazy expressions like above).
# FiPy's augmented version of numpy is called numerix.
# Boundary conditions
# JEG: It doesn't matter that much in 1D, but fluxes are vectors
outerFlux = [100]
outerTemp = 100
BCs = (FixedValue(faces=mesh.getFacesLeft(), value=outerTemp),\
FixedFlux(faces=mesh.getFacesRight(), value=outerFlux),)
eq = TransientTerm() == DiffusionTerm(coeff=therm_cond)
# JEG: no need to declare a solver. FiPy's defaults are pretty good.
elapsedTime = 0
totalElapsedTime = 300
timeStep = 1
desiredResidual = 1e-4
# JEG: Declare the viewers before we start solving so we can watch the evolution
# I made separate viewers because the data range is so different
if __name__ == '__main__':
Tviewer = Viewer(vars=temp)
kviewer = Viewer(vars=therm_cond)
while elapsedTime < totalElapsedTime:
# JEG: therm_cond is not a solution variable,
# so we don't need to update its old value
temp.updateOld()
# JEG: always reset the residual at each timestep or no sweeping will happen
residual = 1e100
while residual > desiredResidual:
# JEG: temp is the only solution variable (and eq can only be solved
for one variable)
# JEG: be sure to apply the boundary conditions!
residual = eq.sweep(var=temp, boundaryConditions=BCs, dt=timeStep)
elapsedTime += timeStep
# JEG: plot the results at each timestep to watch the evolution
if __name__ == '__main__':
Tviewer.plot()
kviewer.plot()
_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
