Hi TJ,

Can you use the direct UFL interpolation feature instead? It's documented in the manual.

Cheers,

David 

On Wed, 18 May 2016 at 20:48 Sun, Tianjiao <tianjiao.sun14@imperial.ac.uk> wrote:
Hi team,

I run into a "NotImplemented" exception when trying to assign values to
an function on extruded mesh of 1 layer. I'm wondering is it really not
implemented or (more likely) just my bad usage? The situation is that I
have this function on DG(0) x CG(1), and I just need to assign values to
the top and bottom surfaces.

Many thanks!
-TJ

code:
---------------------------------------------------------------------------
m = RectangleMesh(5,5,5,5,quadrilateral=True)
mesh = ExtrudedMesh(m, layers=1)
W_h = FiniteElement("DG", "quadrilateral", 0)
W_v = FiniteElement("CG", "interval", 1)
W_elt = TensorProductElement(W_h, W_v)
W = FunctionSpace(mesh, W_elt)

# This Expression reads values from an array
class Ex_intensity(Expression):
     def eval(self, value, X):
         # this is not quite completed, there will be more logic here
         x0 = int(X[0])
         x1 = int(X[1])
         value[:] = p0[x1,x0]

I = Function(W)
I.interpolate(Ex_intensity())
I.at(2.5,2.5,0)

Error message:
---------------------------------------------------------------------------
NotImplementedError                       Traceback (most recent call last)
<ipython-input-12-055ffdeb8731> in <module>()
----> 1 I.at(2.5,2.5,0)

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/firedrake/function.pyc
in at(self, arg, *args, **kwargs)
     492                 else:
     493                     p_result = np.empty(value_shape, dtype=float)
--> 494                     single_eval(points[i:i+1], p_result)
     495                     l_result.append((i, p_result))
     496             except PointNotInDomainError:

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/firedrake/function.pyc
in single_eval(x, buf)
     470         def single_eval(x, buf):
     471             """Helper function to evaluate at a single point."""
--> 472             err = self._c_evaluate(self._ctypes,
     473 x.ctypes.data_as(POINTER(c_double)),
     474 buf.ctypes.data_as(POINTER(c_double)))

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/pyop2/utils.pyc
in __get__(self, obj, cls)
      62         if obj is None:
      63             return self
---> 64         obj.__dict__[self.__name__] = result = self.fget(obj)
      65         return result
      66

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/firedrake/function.pyc
in _ctypes(self)
     406         else:
     407             c_function.n_layers = 1
--> 408         c_function.coords =
coordinates.dat.data.ctypes.data_as(POINTER(c_double))
     409         c_function.coords_map =
coordinates_space.cell_node_list.ctypes.data_as(POINTER(c_int))
     410         c_function.f =
self.dat.data.ctypes.data_as(POINTER(c_double))

<decorator-gen-136> in data(self)

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/pyop2/versioning.pyc
in modifies(method, self, *args, **kwargs)
     113     _force_copies(self)
     114
--> 115     retval = method(self, *args, **kwargs)
     116
     117     self._version_bump()

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/pyop2/base.pyc
in data(self)
    1759
    1760         """
-> 1761         _trace.evaluate(set([self]), set([self]))
    1762         if self.dataset.total_size > 0 and self._data.size == 0
and self.cdim > 0:
    1763             raise RuntimeError("Illegal access: no data
associated with this Dat!")

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/pyop2/base.pyc
in evaluate(self, reads, writes)
     169             to_run = fuse('from_trace', to_run, 0)
     170         for comp in to_run:
--> 171             comp._run()
     172
     173

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/pyop2/base.pyc
in _run(self)
    4073
    4074     def _run(self):
-> 4075         return self.compute()
    4076
    4077     def prepare_arglist(self, iterset, *args):

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/pyop2/base.pyc
in compute(self)
    4115         arglist = self.prepare_arglist(iterset, *self.args)
    4116         fun = self._jitmodule
-> 4117         self._compute(iterset.core_part, fun, *arglist)
    4118         self.halo_exchange_end()
    4119         self._compute(iterset.owned_part, fun, *arglist)

/home/ts2914/projects/firedrake/local/lib/python2.7/site-packages/pyop2/pyparloop.pyc
in _compute(self, part, *arglist)
     106     def _compute(self, part, *arglist):
     107         if part.set._extruded:
--> 108             raise NotImplementedError
     109         subset = isinstance(self._it_space._iterset, base.Subset)
     110

NotImplementedError:



_______________________________________________
firedrake mailing list
firedrake@imperial.ac.uk
https://mailman.ic.ac.uk/mailman/listinfo/firedrake