Dear all,
we (my student and I) have a probably simple problem:
For a 2D problem, we would like to define a velocity in the convection
coefficient, which is a function of the spatial position (as known for simple
laminar flow problems).
Something like u=[f(x),f(y)].
For testing i just thought about: u=[x,y] for each cell
However we fail.
It is no problem if a fixed value for u is given (by the way, if a constant is
used, in which direction does the flow go in a 2D system?).
Also a list like [0.,-.1] is no problem.
However, if I try to use the positions or the transposed positions from
pos= mesh.getCellCenters() or
posT=mesh.getCellCenters()
in eq1 = 0 == ExponentialConvectionTerm(coeff = pos)\
+ ImplicitDiffusionTerm(coeff=a)
We get errors not inthis line, but later, while solving. In the first case it
is:
-----------------------(1 start)
Traceback (most recent call last):
File "D:/Atakan/Cantera/PyWorkBA/FiPyWork/Meshtest0.py", line 31, in <module>
eq1.solve(var=temp, boundaryConditions=bc)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\term.py", line 156, in solve
solver, matrix, RHSvector = self._prepareLinearSystem(var, solver,
boundaryConditions, dt)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\term.py", line 138, in
_prepareLinearSystem
matrix, RHSvector = self.__buildMatrix(var, solver._getMatrixClass(),
boundaryConditions, dt)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\term.py", line 116, in
__buildMatrix
matrix, RHSvector = self._buildMatrix(var, SparseMatrix,
boundaryConditions, dt)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\equation.py", line 94, in
_buildMatrix
dt, self)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\faceTerm.py", line 165, in
_buildMatrix
weight = self._getWeight(mesh, equation=equation)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\convectionTerm.py", line 143, in
_getWeight
alpha = self._Alpha(-self._getGeomCoeff(mesh) / diffCoeff)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\term.py", line 483, in
_getGeomCoeff
self.geomCoeff = self._calcGeomCoeff(mesh)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\convectionTerm.py", line 123, in
_calcGeomCoeff
projectedCoefficients = self.coeff * mesh._getOrientedAreaProjections()
ValueError: shape mismatch: objects cannot be broadcast to a single shape
-----------------------(1 end)
With the other one I get:
-----------------------(2 start)
Traceback (most recent call last):
File "D:/Atakan/Cantera/PyWorkBA/FiPyWork/Meshtest0.py", line 31, in <module>
eq1.solve(var=temp, boundaryConditions=bc)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\term.py", line 156, in solve
solver, matrix, RHSvector = self._prepareLinearSystem(var, solver,
boundaryConditions, dt)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\term.py", line 138, in
_prepareLinearSystem
matrix, RHSvector = self.__buildMatrix(var, solver._getMatrixClass(),
boundaryConditions, dt)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\term.py", line 116, in
__buildMatrix
matrix, RHSvector = self._buildMatrix(var, SparseMatrix,
boundaryConditions, dt)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\equation.py", line 94, in
_buildMatrix
dt, self)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\faceTerm.py", line 168, in
_buildMatrix
self._implicitBuildMatrix(SparseMatrix, L, id1, id2, b, weight['implicit'],
mesh, boundaryConditions, interiorFaces, dt)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\terms\faceTerm.py", line 65, in
_implicitBuildMatrix
L.addAt(numerix.take(coeffMatrix['cell 1 diag'], interiorFaces), id1,
id1)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\tools\numerix.py", line 1082, in take
taken = a.take(indices, axis=axis)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\variables\variable.py", line 1341, in
take
return numerix.take(self.getValue(), ids, axis)
File "D:\Atakan\Cantera\FiPy-2.0\fipy\tools\numerix.py", line 1106, in take
taken = NUMERIX.take(a, indices, axis=axis)
File "C:\Python25\lib\site-packages\numpy\core\fromnumeric.py", line 103, in
take
return take(indices, axis, out, mode)
IndexError: index out of range for array
-----------------------(2 end)
The code is the following:
-----------------------(3 start)
from fipy import *
steps=5
nx = steps
ny = steps
Laenge=1.0
dx = Laenge/steps
dy = Laenge/steps
mesh = Grid2D(dx=dx, dy=dy, nx=nx, ny=ny)
print mesh
positionen=mesh.getCellCenters()
posT=positionen.transpose()
T0=300.
a=3e-3
temp = CellVariable(name = "Temperature",
mesh = mesh,
value = T0)
bc=(FixedValue(faces=mesh.getFacesBottom(), value=T0),\
#FixedValue(faces=mesh.getFacesTop(), value=T0),\
FixedValue(faces=mesh.getFacesLeft(), value=T0),\
FixedValue(faces=mesh.getFacesRight(), value=T0*2.))
print positionen
eq1 = 0 == ExponentialConvectionTerm(coeff = posT)\
+ ImplicitDiffusionTerm(coeff=a)
eq1.solve(var=temp, boundaryConditions=bc)
if __name__ == '__main__':
viewer=Viewer(title='Temperaturverteilung', vars=(temp,), datamax=2000,
datamin=0)
viewer.plot()
-------- (3 end)
So the question is: what is wrong with our approach? How should such a position
dependent velocity profile be implemented?
Thanks in advance!
Burak
