Package sfc :: Package codegeneration :: Module integralcg_uflacs
[hide private]
[frames] | no frames]

Source Code for Module sfc.codegeneration.integralcg_uflacs

  1  # -*- coding: utf-8 -*- 
  2  """ 
  3  This module contains uflacs based code generation tools for the ufc::*_integral classes. 
  4  """ 
  5   
  6  # Copyright (C) 2012-2012 Martin Sandve Alnes and Simula Resarch Laboratory 
  7  # 
  8  # This file is part of SyFi. 
  9  # 
 10  # SyFi is free software: you can redistribute it and/or modify 
 11  # it under the terms of the GNU General Public License as published by 
 12  # the Free Software Foundation, either version 2 of the License, or 
 13  # (at your option) any later version. 
 14  # 
 15  # SyFi is distributed in the hope that it will be useful, 
 16  # but WITHOUT ANY WARRANTY; without even the implied warranty of 
 17  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 
 18  # GNU General Public License for more details. 
 19  # 
 20  # You should have received a copy of the GNU General Public License 
 21  # along with SyFi. If not, see <http://www.gnu.org/licenses/>. 
 22  # 
 23  # First added:  2008-08-13 
 24  # Last changed: 2012-05-25 
 25   
 26  from sfc.common.output import sfc_debug, sfc_error, sfc_assert, sfc_warning 
 27  import ufc_utils 
 28  from sfc.common.utilities import indices_subset, unique 
 29  from sfc.codegeneration.support_code import (generate_sign_correction_support_code, 
 30                                               gen_A_reset_block) 
 31  from sfc.quadrature import (gen_quadrature_rule_definition, 
 32                              gen_quadrature_rule_definitions) 
 33   
 34  from sfc.codegeneration.integralcgbase import IntegralCGBase 
 35   
 36  try: 
 37      import uflacs 
 38  except: 
 39      uflacs = None 
 40  if uflacs is not None: 
 41      import uflacs.backends.cpp2 
 42      import uflacs.backends.cpp2.compiler 
 43      import uflacs.codeutils.format_code_structure 
 44      from uflacs.backends.cpp2.compiler import compile_expression_body 
 45      from uflacs.codeutils.format_code_structure import format_code_structure 
 46      from sfc.codegeneration.sfc_uflacs_formatter import SfcUflacsFormatter 
 47   
48 -class UflacsIntegralCG(IntegralCGBase):
49 - def __init__(self, itgrep):
50 IntegralCGBase.__init__(self, itgrep) 51 sfc_assert(uflacs is not None, "Couldn't import uflacs.") 52 53 # Construct an SFC C++ expression formatter 54 self.target_formatter = SfcUflacsFormatter(itgrep)
55
56 - def hincludes(self):
57 return []
58
59 - def cincludes(self):
60 return []
61
62 - def local_includes(self):
63 return []
64
65 - def header_template(self):
66 raise NotImplementedError("Missing implementation of header_template().")
67
68 - def impl_template(self):
69 raise NotImplementedError("Missing implementation of impl_template().")
70
71 - def generate_support_code(self):
72 code = [] 73 74 # Write quadrature rule if we have one: 75 fr = self.itgrep.formrep 76 # TODO: Place quadrature rules in a separate shared file? 77 if self.itgrep._on_facet: 78 if fr.facet_quad_rules is not None: 79 code += [gen_quadrature_rule_definitions(fr.facet_quad_rules)] 80 else: 81 if fr.quad_rule is not None: 82 code += [gen_quadrature_rule_definition(fr.quad_rule)] 83 84 # Write sign correction code if needed 85 code += [generate_sign_correction_support_code(self.itgrep.formrep)] 86 87 # Declare struct for precomputed variables 88 decls = self.target_formatter.facg.precomp_declarations() 89 code += ["", "struct %s::precomputed" % (self.itgrep.classname,), "{"] 90 code += [" " + decl for decl in decls] 91 code += ["};", ""] 92 93 return '\n'.join(code)
94
95 - def gen_constructor(self):
96 code = [] 97 code += ["pc = new precomputed[%d];" % self.itgrep.num_quadrature_points] 98 code += [""] 99 code += self.target_formatter.define_quadrature_for_loop() 100 code += ["{"] 101 code += self.target_formatter.define_quadrature_loop_vars() 102 lines = [''] 103 for assign in self.target_formatter.facg.precomp_assignments(): 104 lines.append(" " + assign) 105 code.extend(lines) 106 code += ["}"] 107 return '\n'.join(code)
108
110 return ""
111
112 - def gen_initializer_list(self):
113 return ""
114
115 - def gen_destructor(self):
116 code = [] 117 code += ["delete [] pc;"] 118 return '\n'.join(code)
119
120 - def gen_members(self):
121 code = ["", 122 "struct precomputed;", 123 "precomputed * pc;", 124 ""] 125 return '\n'.join(code)
126
127 - def xgen_tabulate_tensor(self):
128 code = [] 129 code += ["int num_quadrature_points = 1;"] 130 code += ["double quadrature_points[1][3] = { { 0.0, 0.0, 0.0 } };"] 131 code += ["double quadrature_weights[1] = { 1.0 };"] 132 code += ["tabulate_tensor(A, w, c,", 133 " num_quadrature_points,", 134 " quadrature_points," 135 " quadrature_weights);"] 136 return code
137
138 - def gen_tabulate_tensor(self):
139 # Get integrand expression 140 integral, = self.itgrep.integrals # TODO: Handle multiple integrals with metadata 141 expr = integral.integrand() 142 #measure = integral.measure() 143 #metadata = measure.metadata() 144 145 # Initial UFC/SFC specific setup code 146 header = gen_A_reset_block(len(self.itgrep.indices)) 147 148 # Compile body with generic compiler routine 149 body = compile_expression_body(expr, self.target_formatter, integrate=True) 150 151 code = format_code_structure([header,body]) 152 return code
153
155 return 'throw std::runtime_error("Not implemented in uflacs mode.");'
156
157 -class UflacsCellIntegralCG(UflacsIntegralCG):
158 """ 159 /// Tabulate the tensor for the contribution from a local cell 160 virtual void tabulate_tensor(double* A, 161 const double * const * w, 162 const cell& c) const = 0; 163 164 /// Tabulate the tensor for the contribution from a local cell 165 /// using the specified reference cell quadrature points/weights 166 virtual void tabulate_tensor(double* A, 167 const double * const * w, 168 const cell& c, 169 unsigned int num_quadrature_points, 170 const double * const * quadrature_points, 171 const double* quadrature_weights) const = 0; 172 """
173 - def __init__(self, itgrep):
174 UflacsIntegralCG.__init__(self, itgrep)
175
176 - def header_template(self):
177 return ufc_utils.cell_integral_header
178
179 - def impl_template(self):
180 return ufc_utils.cell_integral_implementation
181
182 -class UflacsExteriorFacetIntegralCG(UflacsIntegralCG):
183 """ 184 /// Tabulate the tensor for the contribution from a local exterior facet 185 virtual void tabulate_tensor(double* A, 186 const double * const * w, 187 const cell& c, 188 unsigned int facet) const = 0; 189 190 /// Tabulate the tensor for the contribution from a local exterior facet 191 /// using the specified reference cell quadrature points/weights 192 virtual void tabulate_tensor(double* A, 193 const double * const * w, 194 const cell& c, 195 unsigned int num_quadrature_points, 196 const double * const * quadrature_points, 197 const double* quadrature_weights) const = 0; 198 """
199 - def __init__(self, itgrep):
200 UflacsIntegralCG.__init__(self, itgrep)
201
202 - def header_template(self):
203 return ufc_utils.exterior_facet_integral_header
204
205 - def impl_template(self):
206 return ufc_utils.exterior_integral_implementation
207
209 return 'throw std::runtime_error("This function is a bug in UFC. Not implemented in uflacs mode.");'
210
211 -class UflacsInteriorFacetIntegralCG(UflacsIntegralCG):
212 """ 213 /// Tabulate the tensor for the contribution from a local interior facet 214 virtual void tabulate_tensor(double* A, 215 const double * const * w, 216 const cell& c0, 217 const cell& c1, 218 unsigned int facet0, 219 unsigned int facet1) const = 0; 220 221 /// Tabulate the tensor for the contribution from a local interior facet 222 /// using the specified reference cell quadrature points/weights 223 virtual void tabulate_tensor(double* A, 224 const double * const * w, 225 const cell& c, 226 unsigned int num_quadrature_points, 227 const double * const * quadrature_points, 228 const double* quadrature_weights) const = 0; 229 """
230 - def __init__(self, itgrep):
231 UflacsIntegralCG.__init__(self, itgrep)
232
233 - def header_template(self):
234 return ufc_utils.interior_facet_integral_header
235
236 - def impl_template(self):
237 return ufc_utils.interior_integral_implementation
238
240 return 'throw std::runtime_error("This function is a bug in UFC. Not implemented in uflacs mode.");'
241 242 # For dynamic importing of code generation strategy: 243 CellIntegralCG = UflacsCellIntegralCG 244 ExteriorFacetIntegralCG = UflacsExteriorFacetIntegralCG 245 InteriorFacetIntegralCG = UflacsInteriorFacetIntegralCG 246