← Back to branch summary

~ubuntu-branches/ubuntu/wily/dolfin/wily-proposed

« back to all changes in this revision

Viewing changes to site-packages/dolfin/functions/expression.py

  • Committer: Package Import Robot
  • Author(s): Johannes Ring
  • Date: 2014-09-22 14:35:34 UTC
  • mfrom: (1.1.17) (19.1.23 sid)
  • Revision ID: package-import@ubuntu.com-20140922143534-0yi89jyuqbgdxwm9
Tags: 1.4.0+dfsg-4
http://bugs.debian.org/760658
* debian/control: Disable libcgal-dev on i386, mipsel and sparc.
* debian/rules: Remove bad directives in pkg-config file dolfin.pc
  (closes: #760658).
* Remove debian/libdolfin-dev.lintian-overrides.

Show diffs side-by-side

added added

removed removed

Lines of Context:
site-packages/dolfin/functions/expression.py
35
35
 
36
36
"""
37
37
 
38
 
# Copyright (C) 2008-2011 Johan Hake
 
38
# Copyright (C) 2008-2014 Johan Hake
39
39
#
40
40
# This file is part of DOLFIN.
41
41
#
53
53
# along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
54
54
#
55
55
# Modified by Anders Logg, 2008-2009.
56
 
#
57
 
# First added:  2008-11-03
58
 
# Last changed: 2009-12-16
 
56
# Modified by Martin S. Alnaes 2013-2014
59
57
 
60
58
__all__ = ["Expression"]
61
59
 
86
84
        assert(isinstance(expr, cpp.Expression))
87
85
        self._expr = expr
88
86
        dict.__init__(self, **kwargs)
89
 
        
 
87
 
90
88
    def __setitem__(self, name, value):
91
89
        __doc__ = dict.__setitem__.__doc__
92
90
        if name not in self:
102
100
        else:
103
101
            for name, value in other:
104
102
                self[name] = value
105
 
            
 
103
 
106
104
def create_compiled_expression_class(cpp_base):
107
105
    # Check the cpp_base
108
106
    assert(isinstance(cpp_base, (types.ClassType, type)))
109
107
 
110
 
    def __init__(self, cppcode, element=None, cell=None, \
 
108
    def __init__(self, cppcode, element=None, cell=None, domain=None, \
111
109
                 degree=None, name=None, label=None, **kwargs):
112
110
        """Initialize the Expression """
113
111
        # Initialize the cpp base class first and extract value_shape
114
112
        cpp_base.__init__(self)
115
113
        value_shape = tuple(self.value_dimension(i) \
116
 
                                  for i in range(self.value_rank()))
 
114
                            for i in range(self.value_rank()))
117
115
 
118
116
        # Store the value_shape
119
117
        self._value_shape = value_shape
120
118
 
121
119
        # Select an appropriate element if not specified
122
120
        if element is None:
123
 
            element = _auto_select_element_from_shape(value_shape, degree, cell)
 
121
            element = _auto_select_element_from_shape(value_shape, degree, cell, domain)
124
122
        else:
125
123
            # Check that we have an element
126
124
            if not isinstance(element, ufl.FiniteElementBase):
135
133
 
136
134
        # Initialize UFL base class
137
135
        self._ufl_element = element
138
 
        ufl.Coefficient.__init__(self, self._ufl_element)
 
136
        ufl.Coefficient.__init__(self, self._ufl_element, count=self.id())
139
137
 
140
138
        # Set default variables
141
139
        for member, value in kwargs.items():
188
186
 
189
187
    def __init__(self, *args, **kwargs):
190
188
        """Initialize the Expression"""
191
 
        # Get element, cell and degree
 
189
        # Get element, domain and degree
192
190
        element = kwargs.get("element", None)
193
191
        degree = kwargs.get("degree", None)
194
192
        cell = kwargs.get("cell", None)
 
193
        domain = kwargs.get("domain", None)
195
194
        name = kwargs.get("name", None)
196
195
        label = kwargs.get("label", None)
197
196
 
201
200
            from operator import add
202
201
            # First count how many valid kwargs is passed
203
202
            num_valid_kwargs = reduce(add, [kwarg is not None \
204
 
                                      for kwarg in [element, degree, cell, name, label]])
 
203
                                      for kwarg in [element, degree, cell, domain, name, label]])
205
204
            if len(kwargs) != num_valid_kwargs:
206
205
                raise TypeError, "expected only 'kwargs' from one of "\
207
 
                      "'element', 'degree', 'cell', 'name' or 'label'"
 
206
                      "'element', 'degree', 'cell', 'domain', 'name' or 'label'"
208
207
            if len(args) != 0:
209
208
                raise TypeError, "expected no 'args'"
210
209
 
211
210
        # Select an appropriate element if not specified
212
211
        if element is None:
213
212
            element = _auto_select_element_from_shape(self.value_shape(),
214
 
                                                      degree, cell)
 
213
                                                      degree, cell, domain)
215
214
        elif isinstance(element, ufl.FiniteElementBase):
216
215
            pass
217
216
        else:
218
217
            raise TypeError, "The 'element' argument must be a UFL finite"\
219
218
                  " element."
220
219
 
221
 
        # Initialize UFL base class
222
 
        self._ufl_element = element
223
 
        ufl.Coefficient.__init__(self, element)
224
 
 
225
220
        # Initialize cpp_base class
226
221
        cpp.Expression.__init__(self, list(element.value_shape()))
227
222
 
 
223
        # Initialize UFL base class
 
224
        ufl.Coefficient.__init__(self, element, count=self.id())
 
225
 
 
226
        # Why do we store the ufl element here in addition to in ufl.Coefficient? It's available as self.element()...
 
227
        self._ufl_element = element
 
228
 
228
229
        # Calling the user defined_init
229
230
        if user_init is not None:
230
231
            user_init(self, *args, **kwargs)
406
407
            {
407
408
            public:
408
409
 
409
 
              boost::shared_ptr<MeshFunction<unsigned int> > cell_data;
 
410
              std::shared_ptr<MeshFunction<std::size_t> > cell_data;
410
411
 
411
412
              MyFunc() : Expression()
412
413
              {
416
417
                      const ufc::cell& c) const
417
418
              {
418
419
                assert(cell_data);
419
 
                const Cell cell(cell_data->mesh(), c.index);
 
420
                const Cell cell(*cell_data->mesh(), c.index);
420
421
                switch ((*cell_data)[cell.index()])
421
422
                {
422
423
                case 0:
450
451
              class Delta : public Expression
451
452
              {
452
453
              public:
453
 
              
 
454
 
454
455
                Delta() : Expression() {}
455
 
              
 
456
 
456
457
                void eval(Array<double>& values, const Array<double>& data,
457
458
                          const ufc::cell& cell) const
458
459
                { }
459
 
              
460
 
                void update(const boost::shared_ptr<const Function> u, 
 
460
 
 
461
                void update(const std::shared_ptr<const Function> u,
461
462
                            double nu, double dt, double C1,
462
463
                            double U_infty, double chord)
463
464
                {
464
 
                  const boost::shared_ptr<const Mesh> mesh = u->function_space()->mesh();
465
 
                  const boost::shared_ptr<const GenericDofMap> dofmap = u->function_space()->dofmap();
 
465
                  const std::shared_ptr<const Mesh> mesh = u->function_space()->mesh();
 
466
                  const std::shared_ptr<const GenericDofMap> dofmap = u->function_space()->dofmap();
466
467
                  const uint ncells = mesh->num_cells();
467
468
                  uint ndofs_per_cell;
468
 
                  if (ncells > 0) 
 
469
                  if (ncells > 0)
469
470
                  {
470
471
                    CellIterator cell(*mesh);
471
472
                    ndofs_per_cell = dofmap->cell_dimension(cell->index());
472
 
                  } 
473
 
                  else 
 
473
                  }
 
474
                  else
474
475
                  {
475
476
                     return;
476
477
                  }
477
478
                }
478
479
              };
479
480
            }'''
480
 
            
 
481
 
481
482
            e = Expression(code)
482
483
 
483
484
    *3. User-defined expressions by subclassing*
561
562
    # Set the meta class
562
563
    __metaclass__ = ExpressionMetaClass
563
564
 
564
 
    def __new__(cls, cppcode=None, element=None, cell=None, degree=None, \
 
565
    def __new__(cls, cppcode=None, element=None, cell=None, domain=None, degree=None, \
565
566
                name=None, label=None, **kwargs):
566
567
 
567
568
        # If the __new__ function is called because we are instantiating
581
582
        cpp_base, members = compile_expressions([cppcode],
582
583
                                                [generic_function_members])
583
584
        cpp_base, members = cpp_base[0], members[0]
584
 
        
585
 
        # Check passed default arguments 
 
585
 
 
586
        # Check passed default arguments
586
587
        not_allowed = [n for n in dir(cls) if n[0] !="_"]
587
588
        not_allowed += ["cppcode", "user_parameters"]
588
589
        _check_kwargs(members, kwargs, not_allowed)
595
596
 
596
597
    # This method is only included so a user can check what arguments
597
598
    # one should use in IPython using tab completion
598
 
    def __init__(self, cppcode=None, element=None, cell=None, degree=None, \
 
599
    def __init__(self, cppcode=None, element=None, cell=None, domain=None, degree=None, \
599
600
                 name=None, label=None, **kwargs): pass
600
601
 
601
602
    # Reuse the docstring from __new__
725
726
            values_provided = False
726
727
            values = numpy.zeros(value_size, dtype='d')
727
728
 
728
 
        # Check if a cell is defined
729
 
        cell = self.ufl_element().cell()
730
 
        cell_defined = cell is not None
731
 
 
732
 
        if cell_defined:
733
 
            dim = cell.geometric_dimension()
 
729
        # Get dim if element has any domains
 
730
        domains = self.ufl_element().domains()
 
731
        dim = domains[0].geometric_dimension() if domains else None
734
732
 
735
733
        # Assume all args are x argument
736
734
        x = args
738
736
        # If only one x argument has been provided, unpack it if it's an iterable
739
737
        if len(x) == 1:
740
738
            if isinstance(x[0], cpp.Point):
741
 
                if cell_defined:
 
739
                if dim is not None:
742
740
                    x = [x[0][i] for i in xrange(dim)]
743
741
                else:
744
742
                    x = [x[0][i] for i in xrange(3)]
755
753
        if len(x) == 0:
756
754
            raise TypeError, "coordinate argument too short"
757
755
 
758
 
        if cell_defined:
 
756
        if dim is None:
 
757
            warning("Evaluating an Expression without knowing the right geometric dimension, assuming %d is correct." % len(x))
 
758
        else:
759
759
            if len(x) != dim:
760
 
                raise TypeError, "expected the geometry argument to be of "\
761
 
                      "length %d"%dim
 
760
                raise TypeError("expected the geometry argument to be of "\
 
761
                                "length %d" % dim)
762
762
 
763
763
        # The actual evaluation
764
764
        self.eval(values, x)
779
779
        raise TypeError, "Please provide a 'str', 'tuple' or 'list' for "\
780
780
              "the 'cppcode' argument."
781
781
 
782
 
def _auto_select_element_from_shape(shape, degree=None, cell=None):
 
782
def _auto_select_element_from_shape(shape, degree=None, cell=None, domain=None):
783
783
    "Automatically select an appropriate element from cppcode."
784
784
 
 
785
    # Cell to domain conversion TODO: deprecate cell argument
 
786
    if domain is None:
 
787
        if cell is None:
 
788
            pass # No domain will result in undefined element
 
789
        else:
 
790
            domain = ufl.as_domain(cell)
 
791
    else:
 
792
        if cell is None:
 
793
            domain = ufl.as_domain(domain) # Got only domain, good!
 
794
        else:
 
795
            raise ValueError("Got both domain and cell, use only domain.")
 
796
 
785
797
    # Default element, change to quadrature when working
786
798
    Family = "Lagrange"
787
799
 
788
800
    # Check if scalar, vector or tensor valued
789
801
    if len(shape) == 0:
790
 
        element = ufl.FiniteElement(Family, cell, degree)
 
802
        element = ufl.FiniteElement(Family, domain, degree)
791
803
    elif len(shape) == 1:
792
 
        element = ufl.VectorElement(Family, cell, degree, dim=shape[0])
 
804
        element = ufl.VectorElement(Family, domain, degree, dim=shape[0])
793
805
    else:
794
 
        element = ufl.TensorElement(Family, cell, degree, shape=shape)
 
806
        element = ufl.TensorElement(Family, domain, degree, shape=shape)
795
807
 
796
808
    cpp.debug("Automatic selection of expression element: " + str(element))
797
809
 
802
814
    Check that all kwargs passed is either scalars or scalar Constants,
803
815
    and that the name is allowed
804
816
    """
805
 
    
 
817
 
806
818
    # Check passed default values
807
819
    if not all(member in kwargs for member in members):
808
820
        missing = []
817
829
        if name in not_allowed:
818
830
            raise RuntimeError("Parameter name: '%s' is not allowed. It is "\
819
831
                               "part of the interface of Expression" % name)
820
 
            
 
832
 
821
833
    if not (all(isinstance(value, (int, float)) or \
822
834
                (isinstance(value, cpp.GenericFunction) and value.value_size()==1) \
823
835
                for value in kwargs.values())):
824
836
        raise TypeError, "expected default arguments for member variables "\
825
837
              "to be scalars or a scalar GenericFunctions."
826