~njansson/dolfin/hpc

« back to all changes in this revision

Viewing changes to site-packages/dolfin/assemble.py

  • Committer: ilmarw@gogmagog.simula.no
  • Date: 2008-05-29 10:46:14 UTC
  • mfrom: (2122.3.2 trunk)
  • mto: (2668.1.45 trunk)
  • mto: This revision was merged to the branch mainline in revision 2670.
  • Revision ID: ilmarw@gogmagog.simula.no-20080529104614-xjn0a188wtj2wlmw
Repo fix, roll back unwanted Parameter changes.

Show diffs side-by-side

added added

removed removed

Lines of Context:
9
9
rely on the dolfin::Form class which is not used on the Python side."""
10
10
 
11
11
__author__ = "Anders Logg (logg@simula.no)"
12
 
__date__ = "2007-08-15 -- 2008-05-31"
 
12
__date__ = "2007-08-15 -- 2008-04-15"
13
13
__copyright__ = "Copyright (C) 2007-2008 Anders Logg"
14
14
__license__  = "GNU LGPL Version 2.1"
15
15
 
20
20
           "AvgMeshSize",
21
21
           "LinearPDE",
22
22
           "DirichletBC",
23
 
           "PeriodicBC",
24
 
           "compile_functions"]
 
23
           "PeriodicBC"]
25
24
 
26
25
from ffc import *
27
26
from dolfin import *
28
 
from compile_functions import compile_functions
29
 
 
30
 
# Cache for dof maps
31
 
_dof_map_cache = {}
32
 
 
33
 
# Cache for tensors
34
 
_tensor_cache = {}
 
27
 
 
28
 
35
29
 
36
30
# JIT assembler
37
31
def assemble(form, mesh, coefficients=None, dof_maps=None,
39
33
    tensor=None, backend=None, return_dofmaps=False):
40
34
    "Assemble form over mesh and return tensor"
41
35
 
 
36
    
42
37
    # Create empty list of coefficients, filled below
43
 
    _coefficients = ArrayFunctionPtr()
44
 
 
 
38
    coefficients_ = ArrayFunctionPtr()
 
39
    
45
40
    # Check if we need to compile the form (JIT)
46
41
    if hasattr(form, "create_cell_integral"):
47
42
        # UFC form, no need to compile
49
44
        
50
45
        # Extract coefficients
51
46
        if not coefficients is None: 
52
 
            # Compile all strings as dolfin::Function
53
 
            string_expressions = []
54
 
            for c in coefficients:
55
 
                # Note: To allow tuples of floats or ints below, this logic becomes more involved...
56
 
                if isinstance(c, (tuple, str)):
57
 
                    string_expressions.append(c)
58
 
            if string_expressions:
59
 
                compiled_functions = compile_functions(string_expressions, mesh)
60
 
                compiled_functions.reverse()
61
 
            
62
 
            # Build list of coefficients
63
 
            for c in coefficients:
64
 
                # Note: We could generalize this to support more objects 
65
 
                # like sympy expressions, tuples for constant vectors, etc...
 
47
            for c in coefficients:
 
48
                # Note: We could generalize this to support more objects like sympy expressions, tuples for constant vectors, etc...
66
49
                if isinstance(c, (float, int)):
67
50
                    c = cpp_Function(mesh, float(c))
68
 
                elif isinstance(c, (tuple, str)):
69
 
                    c = compiled_functions.pop()
70
 
                _coefficients.push_back(c)
 
51
                coefficients_.push_back(c)
71
52
    else:
72
53
        # FFC form, call JIT compile
73
 
        optimize = dolfin_get("optimize form") or dolfin_get("optimize")
74
 
        (compiled_form, module, form_data) = jit(form, optimize=optimize)
 
54
        (compiled_form, module, form_data) = jit(form)
75
55
        
76
56
        # Extract coefficients
77
57
        for c in form_data.coefficients:
78
 
            _coefficients.push_back(c.f)
79
 
 
80
 
    # FIXME: do we need these lines now? None works fine with Assembler
81
 
    # Create dummy arguments for domains (not yet supported in Python)
 
58
            coefficients_.push_back(c.f)
 
59
    
 
60
    # Create dummy arguments for domains (not yet supported in Python) FIXME: do we need these lines now? None works fine with Assembler.
82
61
    if cell_domains is None:
83
62
        cell_domains = MeshFunction("uint")
84
63
    if exterior_facet_domains is None:
85
64
        exterior_facet_domains = MeshFunction("uint")
86
65
    if interior_facet_domains is None:
87
66
        interior_facet_domains = MeshFunction("uint")
88
 
 
 
67
    
89
68
    # Create dof map set
90
 
    dof_maps = _create_dof_map_set(form, compiled_form, mesh, dof_maps)
 
69
    if dof_maps is None:
 
70
        dof_maps = DofMapSet(compiled_form, mesh)
91
71
 
 
72
    
92
73
    # Create tensor
93
74
    rank = compiled_form.rank()
94
 
    (tensor, reset_tensor) = _create_tensor(form, rank, backend, tensor, reset_tensor)
95
 
 
96
 
    # Set default value for reset_tensor if not specified
97
 
    if reset_tensor is None:
98
 
        reset_tensor = True
99
 
 
 
75
    if tensor is None:
 
76
        if rank == 0:
 
77
            tensor = Scalar()
 
78
        elif rank == 1:
 
79
            if backend: tensor = backend.createVector()
 
80
            else:       tensor = Vector()
 
81
        elif rank == 2:
 
82
            if backend: tensor = backend.createMatrix()
 
83
            else:       tensor = Matrix()
 
84
        else:
 
85
            raise RuntimeError, "Unable to create tensors of rank %d." % rank
 
86
    
100
87
    # Assemble tensor from compiled form
101
 
    cpp_assemble(tensor, compiled_form, mesh, _coefficients, dof_maps,
102
 
                 cell_domains, exterior_facet_domains, interior_facet_domains, reset_tensor)
 
88
    cpp_assemble(tensor, compiled_form, mesh, coefficients_, dof_maps,
 
89
                 cell_domains, exterior_facet_domains, interior_facet_domains, True)
103
90
    
104
91
    # Convert to float for scalars
105
92
    if rank == 0:
111
98
    else:
112
99
        return tensor
113
100
 
114
 
def _create_dof_map_set(form, compiled_form, mesh, dof_maps):
115
 
    "Create dof map set for form"
116
 
 
117
 
    # Check if dof map set is supplied by user
118
 
    if dof_maps:
119
 
        return dof_maps
120
 
 
121
 
    # Check if we should use the cache
122
 
    use_cache = dolfin_get("optimize use dof map cache") or dolfin_get("optimize")
123
 
    if use_cache and form in _dof_map_cache:
124
 
        return _dof_map_cache[form]
125
 
 
126
 
    # Create dof map set
127
 
    dof_maps = DofMapSet(compiled_form, mesh)
128
 
 
129
 
    # Store in cache
130
 
    if use_cache:
131
 
        _dof_map_cache[form] = dof_maps
132
 
 
133
 
    return dof_maps
134
 
 
135
 
def _create_tensor(form, rank, backend, tensor, reset_tensor):
136
 
    "Create tensor for form"
137
 
 
138
 
    # Check if tensor is supplied by user
139
 
    if tensor:
140
 
        return (tensor, reset_tensor)
141
 
 
142
 
    # Decide if we should reset the tensor
143
 
    use_cache = dolfin_get("optimize use tensor cache") or dolfin_get("optimize")
144
 
    if use_cache and reset_tensor is None:
145
 
        reset_tensor = not form in _tensor_cache
146
 
 
147
 
    # Check if we should use the cache
148
 
    if use_cache and form in _tensor_cache:
149
 
        return (_tensor_cache[form], reset_tensor)
150
 
 
151
 
    # Create tensor
152
 
    if rank == 0:
153
 
        tensor = Scalar()
154
 
    elif rank == 1:
155
 
        if backend: tensor = backend.createVector()
156
 
        else:       tensor = Vector()
157
 
    elif rank == 2:
158
 
        if backend: tensor = backend.createMatrix()
159
 
        else:       tensor = Matrix()
160
 
    else:
161
 
        raise RuntimeError, "Unable to create tensors of rank %d." % rank
162
 
 
163
 
    # Store in cache
164
 
    if use_cache:
165
 
        _tensor_cache[form] = tensor
166
 
 
167
 
    return (tensor, reset_tensor)
168
 
 
169
101
# Rename FFC Function
170
102
ffc_Function = Function
171
103