Hello,

I'm making my first steps with FiPy. I studied Chemical
Engineering a while back, but am mostly working as a
software developer these days, mainly with Python.

I want to formulate a 2D diffusion problem and started from
example

  
http://www.ctcms.nist.gov/fipy/examples/diffusion/generated/examples.diffusion.mesh20x20.html

after having read the FiPy introduction chapter at

  http://www.ctcms.nist.gov/fipy/documentation/introduction.html .

This code, very similar to the mesh20x20 example, works as
expected:

  def test():
      # Set up problem.
      mesh = fipy.Grid2D(dx=1.0, dy=1.0, nx=10, ny=10)
      initial_value = 2.0
      concentration = fipy.CellVariable(
        mesh, "concentration", value=initial_value)
      eq = fipy.TransientTerm() == fipy.DiffusionTerm(coeff=1.0)
      viewer = fipy.Viewer(vars=concentration, datamin=0.0, datamax=10.0)
      # Solve.
      time_step = 0.01
      steps = 10
      for step in range(steps):
          eq.solve(var=concentration, dt=time_step)
          viewer.plot()
          raw_input()

apart from printing a warning

RuntimeWarning: invalid value encountered in double_scalars if
  (numerix.sqrt(numerix.sum(errorVector**2)) / error0) <= self.tolerance

Obviously, the above problem isn't very interesting. I start
with a concentration of 2.0 everywhere and due to the
default boundary conditions, which don't allow any mass
transfer at the boundaries, the concentration stays the same
throughout the simulation.

Now I want to set different initial concentration values for
different mesh points. How do I do this?

By reading some of the FiPy code in `meshVariable.py`, I
found out that I can pass an array as the `value` argument
of the `CellVariable` constructor, like this:

  def test():
      # Set up problem.
      mesh = fipy.Grid2D(dx=1.0, dy=1.0, nx=10, ny=10)
      initial_value = fipy.numerix.zeros((10, 10))
      initial_value[:3, :3] = 5
      concentration = fipy.CellVariable(
        mesh, "concentration", value=initial_value, rank=0)
      eq = fipy.TransientTerm() == fipy.DiffusionTerm(coeff=1.0)
      viewer = fipy.Viewer(vars=concentration, datamin=0.0, datamax=10.0)
      # Solve.
      time_step = 0.01
      steps = 10
      for step in range(steps):
          eq.solve(var=concentration, dt=time_step)
          viewer.plot()
          raw_input()

However, this gives me

  Traceback (most recent call last):
    File "fipy_experiment.py", line 45, in <module>
      test()
    File "fipy_experiment.py", line 34, in test
      viewer = fipy.Viewer(vars=concentration, datamin=0.0, datamax=10.0)
    File 
"/home/schwa/.virtualenvs/diffusion/lib/python2.7/site-packages/fipy/viewers/__init__.py",
 line 136, in Viewer
      raise ImportError, "Failed to import a viewer: %s" % str(errors)
  ImportError: Failed to import a viewer: ['matplotlib: The mesh must be a 
Mesh2D instance', 'gist: __call__() takes at least 2 arguments (1 given)', 
'gnuplot: __call__() takes at least 2 arguments (1 given)', 'mayavi: No module 
named enthought.tvtk.api']

When I set the rank argument to 0 explicitly, as I
understand from the `_MeshVariable.__init__` docstring, I
still get the same traceback. (I understand the word
"element" in the mentioned docstring as referring to the
concentration value in my case, which is a scalar.) I also
tried setting `elementShape` to `()` explicitly, but this
gives me

  Traceback (most recent call last):
    File "fipy_experiment.py", line 46, in <module>
      test()
    File "fipy_experiment.py", line 33, in test
      elementshape=())
    File 
"/home/schwa/.virtualenvs/diffusion/lib/python2.7/site-packages/fipy/variables/cellVariable.py",
 line 69, in __init__
      rank=rank, elementshape=elementshape, unit=unit)
    File 
"/home/schwa/.virtualenvs/diffusion/lib/python2.7/site-packages/fipy/variables/meshVariable.py",
 line 133, in __init__
      array=array, cached=cached)
    File 
"/home/schwa/.virtualenvs/diffusion/lib/python2.7/site-packages/fipy/variables/variable.py",
 line 112, in __init__
      self._setValueInternal(value=value, unit=unit, array=array)
    File 
"/home/schwa/.virtualenvs/diffusion/lib/python2.7/site-packages/fipy/variables/variable.py",
 line 666, in _setValueInternal
      self._value = self._makeValue(value=value, unit=unit, array=array)
    File 
"/home/schwa/.virtualenvs/diffusion/lib/python2.7/site-packages/fipy/variables/variable.py",
 line 698, in _makeValue
      array[:] = value
  ValueError: output operand requires a reduction, but reduction is not enabled

I don't know what this is supposed to tell me. :-)

In summary: How can I specify different scalar initial
values for the points in the mesh?

P. S.: I just found

  http://article.gmane.org/gmane.comp.python.fipy/2803

Is this the answer to my question? If yes, could this be
made somewhat easier? For example, would it be feasible to
wrap that functionality in some convenience API of FiPy?

Stefan

_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
  [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]

Reply via email to