/usr/share/pyshared/dolfin/functions/function.py

"""This module handles the Function class in Python.
"""
# Copyright (C) 2009 Johan Hake
#
# This file is part of DOLFIN.
#
# DOLFIN is free software: you can redistribute it and/or modify
# it under the terms of the GNU Lesser General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# DOLFIN is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
#
# First added:  2009-10-06
# Last changed: 2011-04-18

__all__ = ["Function", "TestFunction", "TrialFunction", "Argument",
           "TestFunctions", "TrialFunctions"]

import types

# Import UFL and SWIG-generated extension module (DOLFIN C++)
import ufl
import dolfin.cpp as cpp
import numpy

from dolfin.functions.functionspace import FunctionSpaceBase

class MetaNoEvalOverloading(type):
    def __init__(mcs, name, bases, dictionary):
        if "eval" in dictionary:
            raise TypeError, "cannot overload 'eval'"

class Function(ufl.Coefficient, cpp.Function):
    """
    This class represents a function :math:`u_h` in a finite
    element function space :math:`V_h`, given by

    .. math::

        u_h = \sum_{i=1}^n U_i \phi_i,

    where :math:`\{\phi_i\}_{i=1}^n` is a basis for :math:`V_h`,
    and :math:`U` is a vector of expansion coefficients for
    :math:`u_h`.

    *Arguments*
        There is a maximum of two arguments. The first argument must be a
        Function or a :py:class:`FunctionSpace
        <dolfin.functions.functionspace.FunctionSpace>`.

        If instantiated from another Function, the (optional)
        second argument must be an integer denoting the number
        of sub functions to extract.

    *Examples*
        Create a Function:

        - from a :py:class:`FunctionSpace
          <dolfin.functions.functionspace.FunctionSpace>` ``V``

          .. code-block:: python

              f = Function(V)

        - from another Function ``f``

          .. code-block:: python

              g = Function(f)

        - from a :py:class:`FunctionSpace
          <dolfin.functions.functionspace.FunctionSpace>` ``V`` and a
          :py:class:`GenericVector <dolfin.cpp.GenericVector>` ``v``

          .. code-block:: python

              g = Function(V, v)

        - from a :py:class:`FunctionSpace
          <dolfin.functions.functionspace.FunctionSpace>` and a
          filename containg a :py:class:`GenericVector
          <dolfin.cpp.GenericVector>`

          .. code-block:: python

              g = Function(V, 'MyVectorValues.xml')

    """

    __metaclass__ = MetaNoEvalOverloading

    def __init__(self, *args):
        """Initialize Function."""
        # Check arguments
        if len(args) == 0:
            raise TypeError, "expected 1 or more arguments"

        if not isinstance(args[0], (FunctionSpaceBase, Function)):
            raise TypeError, "expected a FunctionSpace or a Function as argument 1"

        # If using the copy constuctor
        if isinstance(args[0], Function):
            other = args[0]
            # If using the copy constuctor
            if len(args) == 1:
                # Instantiate base classes
                cpp.Function.__init__(self, other)
                ufl.Coefficient.__init__(self, other._element)
                return

            # If using sub-function constructor
            elif len(args) == 2 and isinstance(args[1], int):
                i = args[1]
                num_sub_spaces = other.function_space().num_sub_spaces()

                if num_sub_spaces == 1:
                    raise RuntimeError, "No subfunctions to extract"
                if not i < num_sub_spaces:
                    raise RuntimeError, "Can only extract subfunctions with i = 0..%d"% num_sub_spaces
                cpp.Function.__init__(self, other, i)
                ufl.Coefficient.__init__(self, self.function_space().ufl_element())
                return
            else:
                raise TypeError, "expected one or two arguments when instantiating from another Function"

        V = args[0]

        # Instantiate ufl base class
        ufl.Coefficient.__init__(self, V.ufl_element())

        # Passing only the FunctionSpace
        if len(args) == 1:
            # Instantiate cpp base classes
            cpp.Function.__init__(self, V)
        elif len(args) == 2:
            # If passing FunctionSpace together with cpp.Function
            # Attached passed FunctionSpace and initialize the cpp.Function
            # using the passed Function
            if isinstance(args[1], cpp.Function):
                if args[1].function_space().dim() != V.dim():
                    raise ValueError, "non matching dimensions on passed FunctionSpaces"

                cpp.Function.__init__(self, args[1])
            else:
                cpp.Function.__init__(self, *args)
        else:
            raise TypeError, "too many arguments"

    def _sub(self, i, deepcopy = False):
        """Return a sub function.

        The sub functions are numbered from i = 0..N-1, where N is the
        total number of sub spaces.

        *Arguments*
            i
                The number of the sub function

        """
        if not isinstance(i, int):
            raise TypeError, "expects an 'int' as first argument"
        num_sub_spaces = self.function_space().num_sub_spaces()
        if num_sub_spaces == 1:
            raise RuntimeError, "No subfunctions to extract"
        if not i < num_sub_spaces:
            raise RuntimeError, "Can only extract subfunctions with i = 0..%d"% num_sub_spaces

        # Create and instantiate the Function
        if deepcopy:
            return Function(self.function_space().sub(i), cpp.Function._sub(self, i))
        else:
            return Function(self, i)

    def split(self, deepcopy=False):
        """
        Extract any sub functions.

        A sub function can be extracted from a discrete function that
        is in a :py:class:`MixedFunctionSpace
        <dolfin.functions.functionspace.MixedFunctionSpace>` or in a
        :py:class:`VectorFunctionSpace
        <dolfin.functions.functionspace.VectorFunctionSpace>`. The sub
        function resides in the subspace of the mixed space.

        *Arguments*
            deepcopy
                Copy sub function vector instead of sharing

        """

        num_sub_spaces = self.function_space().num_sub_spaces()
        if num_sub_spaces == 1:
            raise RuntimeError, "No subfunctions to extract"
        return tuple(self._sub(i, deepcopy) for i in xrange(num_sub_spaces))

    def ufl_element(self):
        """Return ufl element"""
        return self._element

    def __str__(self):
        """Return a pertty print representation of it self.

        Use ufl.__str__ is used for this.
        """
        return ufl.Coefficient.__str__(self)

    def __repr__(self):
        """Return a str repr of it self.

        Must use ufl.__repr__ for this"""
        return ufl.Coefficient.__repr__(self)

    def str(self, verose=False):
        """Return an informative str representation of itself"""
        # FIXME: We might change this using rank and dimension instead
        return "<Function in %s>" % str(self.function_space())

    def ufl_evaluate(self, x, component, derivatives):
        """Function used by ufl to evaluate the Function"""
        import numpy
        import ufl
        assert derivatives == () # TODO: Handle derivatives

        if component:
            shape = self.shape()
            assert len(shape) == len(component)
            value_size = ufl.common.product(shape)
            index = ufl.common.component_to_index(component, shape)
            values = numpy.zeros(value_size)
            self(*x, values=values)
            return values[index]
        else:
            # Scalar evaluation
            return self(*x)

    def __call__(self, *args, **kwargs):
        """
        Evaluates the Function.

        *Examples*
            1) Using an iterable as x:

              .. code-block:: python

                  fs = Expression("sin(x[0])*cos(x[1])*sin(x[3])")
                  x0 = (1.,0.5,0.5)
                  x1 = [1.,0.5,0.5]
                  x2 = numpy.array([1.,0.5,0.5])
                  v0 = fs(x0)
                  v1 = fs(x1)
                  v2 = fs(x2)

            2) Using multiple scalar args for x, interpreted as a
            point coordinate

              .. code-block:: python

                  v0 = f(1.,0.5,0.5)

            3) Using a Point

              .. code-block:: python

                  p0 = Point(1.,0.5,0.5)
                  v0 = f(p0)

            3) Passing return array

              .. code-block:: python

                  fv = Expression(("sin(x[0])*cos(x[1])*sin(x[3])",
                               "2.0","0.0"))
                  x0 = numpy.array([1.,0.5,0.5])
                  v0 = numpy.zeros(3)
                  fv(x0, values = v0)

              .. note::

                  A longer values array may be passed. In this way one can fast
                  fill up an array with different evaluations.

              .. code-block:: python

                  values = numpy.zeros(9)
                  for i in xrange(0,10,3):
                      fv(x[i:i+3], values = values[i:i+3])

        """

        if len(args)==0:
            raise TypeError, "expected at least 1 argument"

        # Test for ufl restriction
        if len(args) == 1 and args[0] in ('+','-'):
            return ufl.Coefficient.__call__(self, *args)

        # Test for ufl mapping
        if len(args) == 2 and isinstance(args[1], dict) and self in args[1]:
            return ufl.Coefficient.__call__(self, *args)

        # Some help variables
        value_size = ufl.common.product(self.ufl_element().value_shape())

        # If values (return argument) is passed, check the type and length
        values = kwargs.get("values", None)
        if values is not None:
            if not isinstance(values, numpy.ndarray):
                raise TypeError, "expected a NumPy array for 'values'"
            if len(values) != value_size or \
                   not numpy.issubdtype(values.dtype, 'd'):
                raise TypeError, "expected a double NumPy array of length"\
                      " %d for return values."%value_size
            values_provided = True
        else:
            values_provided = False
            values = numpy.zeros(value_size, dtype='d')

        # Get the dimension of the cell
        dim = self.ufl_element().cell().geometric_dimension()

        # Assume all args are x argument
        x = args

        # If only one x argument has been provided, unpack it if it's an iterable
        if len(x) == 1:
            if isinstance(x[0], cpp.Point):
                x = (x[0][i] for i in xrange(dim))
            elif hasattr(x[0], '__iter__'):
                x = x[0]

        # Convert it to an 1D numpy array
        try:
            x = numpy.fromiter(x, 'd')
        except (TypeError, ValueError, AssertionError):
            raise TypeError, "expected scalar arguments for the coordinates"

        if len(x) == 0:
            raise TypeError, "coordinate argument too short"

        if len(x) != dim:
            raise TypeError, "expected the geometry argument to be of "\
                  "length %d"%dim

        # The actual evaluation
        self.eval(values, x)

        # If scalar return statement, return scalar value.
        if value_size == 1 and not values_provided:
            return values[0]

        return values

#--- Subclassing of ufl.{Basis, Trial, Test}Function ---

class Argument(ufl.Argument):
    """An Argument represents a possibly differentiated component
    of an argument on the reference cell.
    """
    def __init__(self, V, index=None):
        if not isinstance(V, FunctionSpaceBase):
            raise TypeError, "Illegal argument for creation of Argument, not a FunctionSpace: " + str(V)
        ufl.Argument.__init__(self, V.ufl_element(), index)
        self._V = V

    def function_space(self):
        "Return the FunctionSpace"
        return self._V

def TrialFunction(V):
    """A TrialFunction is the Argument with the next lowest primary
    index. We simply pick an index lower than almost all others (-1).
    """
    return Argument(V, -1)

def TestFunction(V):
    """A TestFunction is the Argument with the lowest primary
    index. We simply pick an index lower than all others (-2).
    """
    return Argument(V, -2)

#--- TestFunctions and TrialFunctions ---

def TestFunctions(V):
    "Create test functions from mixed function space."
    return ufl.split(TestFunction(V))

def TrialFunctions(V):
    "Create trial functions from mixed function space."
    return ufl.split(TrialFunction(V))
python-dolfin 1.0.0-1 / usr / share / pyshared / dolfin / functions / function.py
