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

Source Code for Module sfc.codegeneration.codegeneration

  1  #!/usr/bin/env python 
  2  # -*- coding: utf-8 -*- 
  3  """ 
  4  Codegeneration utilities. 
  5  """ 
  6   
  7  # Copyright (C) 2008-2009 Martin Sandve Alnes and Simula Resarch Laboratory 
  8  # 
  9  # This file is part of SyFi. 
 10  # 
 11  # SyFi is free software: you can redistribute it and/or modify 
 12  # it under the terms of the GNU General Public License as published by 
 13  # the Free Software Foundation, either version 2 of the License, or 
 14  # (at your option) any later version. 
 15  # 
 16  # SyFi is distributed in the hope that it will be useful, 
 17  # but WITHOUT ANY WARRANTY; without even the implied warranty of 
 18  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 
 19  # GNU General Public License for more details. 
 20  # 
 21  # You should have received a copy of the GNU General Public License 
 22  # along with SyFi. If not, see <http://www.gnu.org/licenses/>. 
 23  # 
 24  # Modified by Kent-Andre Mardal, 2010. 
 25  # 
 26  # First added:  2008-08-13 
 27  # Last changed: 2009-05-15 
 28   
 29  from itertools import chain 
 30  import re 
 31  import ufl 
 32  import ufc_utils 
 33  import SyFi 
 34   
 35  from ufl.classes import Measure 
 36  from ufl.algorithms import FormData, preprocess 
 37   
 38  from dolfin_utils.wrappers import generate_dolfin_code, UFCFormNames 
 39   
 40  from sfc.common.utilities import unique 
 41  from sfc.common.output import sfc_assert, sfc_error, sfc_info, sfc_debug 
 42  #from sfc.common import ParameterDict 
 43  from sfc.common.options import default_parameters 
 44  from sfc.representation import ElementRepresentation, FormRepresentation, \ 
 45          CellIntegralRepresentation, ExteriorFacetIntegralRepresentation, \ 
 46          InteriorFacetIntegralRepresentation 
 47  from sfc.codegeneration.formcg import FormCG 
 48  from sfc.codegeneration.dofmapcg import DofMapCG 
 49  from sfc.codegeneration.finiteelementcg import FiniteElementCG 
 50  from sfc.codegeneration.cellintegralcg import CellIntegralCG 
 51  from sfc.codegeneration.exteriorfacetintegralcg import ExteriorFacetIntegralCG 
 52  from sfc.codegeneration.interiorfacetintegralcg import InteriorFacetIntegralCG 
 53   
 54  # version %(version)s 
 55  _header_template = r"""/* 
 56   * %(name)s.h 
 57   * 
 58   * This file was automatically generated by SFC. 
 59   * 
 60   * http://www.fenics.org/syfi/ 
 61   * 
 62   */ 
 63   
 64  #ifndef GENERATED_SFC_CODE_%(name)s 
 65  #define GENERATED_SFC_CODE_%(name)s 
 66   
 67  #include "ufc.h" 
 68   
 69  %(includes)s 
 70   
 71  //namespace sfc 
 72  //{ 
 73  %(body)s 
 74  //} 
 75   
 76  #endif 
 77  """ 
 78   
 79  # version %(version)s 
 80  _implementation_template = r"""/* 
 81   * %(name)s.cpp 
 82   * 
 83   * This file was automatically generated by SFC. 
 84   * 
 85   * http://www.fenics.org/syfi/ 
 86   * 
 87   */ 
 88   
 89  #include "ufc.h" 
 90   
 91  #include "%(name)s.h" 
 92   
 93  %(includes)s 
 94   
 95  //namespace sfc 
 96  //{ 
 97  %(body)s 
 98  //} 
 99  """ 
100   
101 -def common_system_headers():
102 return ("iostream", "cmath", "stdexcept", "cstring")
103
104 -def apply_code_dict(format_string, code_dict):
105 "Template formatting with improved checking for template argument mismatch." 106 expected_keys = set(re.findall('%\((\w+)\)', format_string)) 107 actual_keys = set(code_dict.keys()) 108 109 missing_keys = expected_keys.difference(actual_keys) 110 if missing_keys: 111 print "When formatting code using template, missing keys:" 112 print "\n".join(sorted(missing_keys)) 113 print 114 115 #superfluous_keys = actual_keys.difference(expected_keys) 116 #if superfluous_keys: 117 # print "When formatting code using template, got superfluous keys:" 118 # print "\n".join(sorted(superfluous_keys)) 119 # print 120 121 return format_string % code_dict
122
123 -def generate_finite_element_code(ferep):
124 sfc_debug("Entering generate_finite_element_code") 125 classname = ferep.finite_element_classname 126 cg = FiniteElementCG(ferep) 127 code_dict = cg.generate_code_dict() 128 supportcode = cg.generate_support_code() 129 130 hcode = apply_code_dict(ufc_utils.finite_element_header, code_dict) 131 ccode = supportcode + "\n"*3 132 ccode += apply_code_dict(ufc_utils.finite_element_implementation, code_dict) 133 134 includes = cg.hincludes() + cg.cincludes() 135 136 system_headers = common_system_headers() 137 local_headers = unique("%s.h" % e.finite_element_classname for e in ferep.sub_elements) 138 139 hincludes = "\n".join('#include "%s"' % inc for inc in cg.hincludes()) 140 cincludes = "\n".join('#include <%s>' % f for f in system_headers) 141 cincludes += "\n" 142 cincludes += "\n".join('#include "%s"' % inc for inc in chain(cg.cincludes(), local_headers)) 143 144 hcode = _header_template % { "body": hcode, "name": classname, "includes": hincludes } 145 ccode = _implementation_template % { "body": ccode, "name": classname, "includes": cincludes } 146 sfc_debug("Leaving generate_finite_element_code") 147 return classname, (hcode, ccode), includes
148
149 -def generate_dof_map_code(ferep):
150 sfc_debug("Entering generate_dof_map_code") 151 classname = ferep.dof_map_classname 152 cg = DofMapCG(ferep) 153 code_dict = cg.generate_code_dict() 154 supportcode = cg.generate_support_code() 155 156 hcode = apply_code_dict(ufc_utils.dofmap_header, code_dict) 157 ccode = supportcode + "\n"*3 158 ccode += apply_code_dict(ufc_utils.dofmap_implementation, code_dict) 159 160 includes = cg.hincludes() + cg.cincludes() 161 162 system_headers = common_system_headers() 163 local_headers = unique("%s.h" % e.dof_map_classname for e in ferep.sub_elements) 164 165 hincludes = "\n".join('#include "%s"' % inc for inc in cg.hincludes()) 166 cincludes = "\n".join('#include <%s>' % f for f in system_headers) 167 cincludes += "\n" 168 cincludes += "\n".join('#include "%s"' % inc for inc in chain(cg.cincludes(), local_headers)) 169 170 hcode = _header_template % { "body": hcode, "name": classname, "includes": hincludes } 171 ccode = _implementation_template % { "body": ccode, "name": classname, "includes": cincludes } 172 sfc_debug("Leaving generate_dof_map_code") 173 return classname, (hcode, ccode), includes
174
175 -def generate_cell_integral_code(integrals, formrep):
176 sfc_debug("Entering generate_cell_integral_code") 177 itgrep = CellIntegralRepresentation(integrals, formrep) 178 179 cg = CellIntegralCG(itgrep) 180 code_dict = cg.generate_code_dict() 181 supportcode = cg.generate_support_code() 182 183 # TODO: Support code might be placed in unnamed namespace? 184 hcode = apply_code_dict(ufc_utils.cell_integral_header, code_dict) 185 ccode = supportcode + "\n"*3 + apply_code_dict(ufc_utils.cell_integral_implementation, code_dict) 186 187 includes = cg.hincludes() + cg.cincludes() 188 189 system_headers = common_system_headers() 190 191 hincludes = "\n".join('#include "%s"' % inc for inc in cg.hincludes()) 192 cincludes = "\n".join('#include <%s>' % f for f in system_headers) 193 cincludes += "\n" 194 cincludes += "\n".join('#include "%s"' % inc for inc in cg.cincludes()) 195 196 hcode = _header_template % { "body": hcode, "name": itgrep.classname, "includes": hincludes } 197 ccode = _implementation_template % { "body": ccode, "name": itgrep.classname, "includes": cincludes } 198 sfc_debug("Leaving generate_cell_integral_code") 199 return itgrep.classname, (hcode, ccode), includes
200
201 -def generate_exterior_facet_integral_code(integrals, formrep):
202 sfc_debug("Entering generate_exterior_facet_integral_code") 203 itgrep = ExteriorFacetIntegralRepresentation(integrals, formrep) 204 205 cg = ExteriorFacetIntegralCG(itgrep) 206 code_dict = cg.generate_code_dict() 207 supportcode = cg.generate_support_code() 208 209 # TODO: Support code might be placed in unnamed namespace? 210 hcode = apply_code_dict(ufc_utils.exterior_facet_integral_header, code_dict) 211 ccode = supportcode + "\n"*3 + apply_code_dict(ufc_utils.exterior_facet_integral_implementation, code_dict) 212 213 includes = cg.hincludes() + cg.cincludes() 214 215 system_headers = common_system_headers() 216 217 hincludes = "\n".join('#include "%s"' % inc for inc in cg.hincludes()) 218 cincludes = "\n".join('#include <%s>' % f for f in system_headers) 219 cincludes += "\n" 220 cincludes += "\n".join('#include "%s"' % inc for inc in cg.cincludes()) 221 222 hcode = _header_template % { "body": hcode, "name": itgrep.classname, "includes": hincludes } 223 ccode = _implementation_template % { "body": ccode, "name": itgrep.classname, "includes": cincludes } 224 sfc_debug("Leaving generate_exterior_facet_integral_code") 225 return itgrep.classname, (hcode, ccode), includes
226
227 -def generate_interior_facet_integral_code(integrals, formrep):
228 sfc_debug("Entering generate_interior_facet_integral_code") 229 itgrep = InteriorFacetIntegralRepresentation(integrals, formrep) 230 231 cg = InteriorFacetIntegralCG(itgrep) 232 code_dict = cg.generate_code_dict() 233 supportcode = cg.generate_support_code() 234 235 hcode = apply_code_dict(ufc_utils.interior_facet_integral_header, code_dict) 236 ccode = supportcode + "\n"*3 + apply_code_dict(ufc_utils.interior_facet_integral_implementation, code_dict) 237 238 includes = cg.hincludes() + cg.cincludes() 239 240 system_headers = common_system_headers() 241 242 hincludes = "\n".join('#include "%s"' % inc for inc in cg.hincludes()) 243 cincludes = "\n".join('#include <%s>' % f for f in system_headers) 244 cincludes += "\n" 245 cincludes += "\n".join('#include "%s"' % inc for inc in cg.cincludes()) 246 247 hcode = _header_template % { "body": hcode, "name": itgrep.classname, "includes": hincludes } 248 ccode = _implementation_template % { "body": ccode, "name": itgrep.classname, "includes": cincludes } 249 sfc_debug("Leaving generate_interior_facet_integral_code") 250 return itgrep.classname, (hcode, ccode), includes
251
252 -def generate_form_code(formrep):
253 sfc_debug("Entering generate_form_code") 254 cg = FormCG(formrep) 255 code_dict = cg.generate_code_dict() 256 supportcode = cg.generate_support_code() 257 258 hcode = apply_code_dict(ufc_utils.form_header, code_dict) 259 ccode = supportcode + "\n"*3 + apply_code_dict(ufc_utils.form_implementation, code_dict) 260 261 includes = cg.hincludes() + cg.cincludes() 262 263 system_headers = common_system_headers() 264 local_headers = unique(chain(formrep.fe_names, formrep.dm_names, \ 265 sorted(formrep.itg_names.values()), cg.cincludes())) 266 267 hincludes = "\n".join('#include "%s"' % inc for inc in cg.hincludes()) 268 cincludes = "\n".join('#include <%s>' % f for f in system_headers) 269 cincludes += "\n" 270 cincludes += "\n".join('#include "%s.h"' % f for f in local_headers) 271 272 hcode = _header_template % { "body": hcode, "name": formrep.classname, "includes": hincludes } 273 ccode = _implementation_template % { "body": ccode, "name": formrep.classname, "includes": cincludes } 274 sfc_debug("Leaving generate_form_code") 275 return (hcode, ccode), includes
276
277 -def write_file(filename, text):
278 f = open(filename, "w") 279 f.write(text) 280 f.close()
281
282 -def write_code(classname, code):
283 sfc_debug("Entering write_code") 284 if isinstance(code, tuple): 285 # Code is split in a header and implementation file 286 hcode, ccode = code 287 hname = classname + ".h" 288 cname = classname + ".cpp" 289 open(hname, "w").write(hcode) 290 open(cname, "w").write(ccode) 291 ret = (hname, cname) 292 else: 293 # All code is combined in a header file 294 name = classname + ".h" 295 open(name, "w").write(code) 296 ret = (name,) 297 sfc_debug("Leaving write_code") 298 return ret
299
300 -def compiler_input(input, objects=None, common_cell=None):
301 """Map different kinds of input to a list of 302 UFL elements and a list of FormData instances. 303 304 The following input formats are allowed: 305 - ufl.Form 306 - ufl.algorithms.FormData 307 - ufl.FiniteElementBase 308 - list of the above 309 310 Returns: 311 elements, formdatas 312 """ 313 sfc_debug("Entering compiler_input") 314 sfc_assert(input, "Got no input!") 315 fd = [] 316 fe = [] 317 if not isinstance(input, list): 318 input = [input] 319 for d in input: 320 if isinstance(d, ufl.form.Form): 321 d = d.compute_form_data(object_names=objects, common_cell=common_cell) 322 if isinstance(d, ufl.algorithms.formdata.FormData): 323 fd.append(d) 324 fe.extend(d.sub_elements) 325 elif isinstance(d, ufl.FiniteElementBase): 326 if d.cell().is_undefined(): 327 sfc_assert(common_cell and not common_cell.is_undefined(), 328 "Element no defined cell, cannot compile this.") 329 d = d.reconstruct(cell=common_cell) 330 sfc_assert(not d.cell().is_undefined(), "Still undefined?") 331 fe.append(d) 332 else: 333 sfc_error("Not a FormData or FiniteElementBase object: %s" % str(type(d))) 334 fe = sorted(set(fe)) 335 sfc_debug("Leaving compiler_input") 336 for x in fe: 337 sfc_assert(not x.cell().is_undefined(), "Can never have undefined cells at this point!") 338 for x in fd: 339 sfc_assert(not x.cell.is_undefined(), "Can never have undefined cells at this point!") 340 return fe, fd
341 342 dolfin_header_template = """/* 343 * DOLFIN wrapper code generated by the SyFi Form Compiler. 344 */ 345 346 %s 347 """ 348
349 -def generate_code(input, objects, options = None, common_cell=None):
350 """Generate code from input and options. 351 352 @param input: 353 TODO 354 @param options: 355 TODO 356 """ 357 sfc_debug("Entering generate_code") 358 359 if options is None: 360 options = default_parameters() 361 362 ufl_elements, formdatas = compiler_input(input, objects, common_cell=common_cell) 363 364 filenames = [] 365 needed_files = [] 366 formnames = [] 367 368 # Generate UFC code for elements 369 element_reps = {} 370 last_nsd = None 371 generated_elements = set() 372 for e in ufl_elements: 373 # Initialize global variables in SyFi with the right space dimension (not very nice!) 374 nsd = e.cell().geometric_dimension() 375 if not nsd == last_nsd and isinstance(nsd, int): 376 SyFi.initSyFi(nsd) 377 last_nsd = nsd 378 379 # Construct ElementRepresentation objects for all elements 380 quad_rule = None # TODO: What to do with this one? 381 assert not "Quadrature" in e.family() # No quadrature rule defined! 382 if e in element_reps: 383 continue 384 385 erep = ElementRepresentation(e, quad_rule, options, element_reps) 386 element_reps[e] = erep 387 388 # Build flat list of all subelements 389 todo = [erep] 390 ereps = [] 391 while todo: 392 erep = todo.pop() 393 # Skip already generated! 394 if not erep.ufl_element in generated_elements: 395 generated_elements.add(erep.ufl_element) 396 ereps.append(erep) 397 # Queue subelements for inspection 398 for subrep in erep.sub_elements: 399 todo.append(subrep) 400 401 for erep in ereps: 402 # Generate code for finite_element 403 classname, code, includes = generate_finite_element_code(erep) 404 filenames.extend( write_code(classname, code) ) 405 needed_files.extend(includes) 406 407 # Generate code for dof_map 408 classname, code, includes = generate_dof_map_code(erep) 409 filenames.extend( write_code(classname, code) ) 410 needed_files.extend(includes) 411 412 # Generate UFC code for forms 413 for formdata in formdatas: 414 # Initialize global variables in SyFi with the right space dimension (not very nice!) 415 nsd = formdata.geometric_dimension 416 if not nsd == last_nsd: 417 SyFi.initSyFi(nsd) 418 last_nsd = nsd 419 420 # This object extracts and collects various information about the ufl form 421 formrep = FormRepresentation(formdata, element_reps, options) 422 423 pf = formdata.preprocessed_form 424 425 # TODO: Get integrals from formrep, in case of modifications there? Or rather make sure that FormRep doesn't do any such things. 426 ig = pf.integral_groups() 427 428 # Generate code for cell integrals 429 for domain in pf.domains(Measure.CELL): 430 integrals = ig[domain] 431 classname, code, includes = generate_cell_integral_code(integrals, formrep) 432 filenames.extend( write_code(classname, code) ) 433 needed_files.extend(includes) 434 435 # Generate code for exterior facet integrals 436 for domain in pf.domains(Measure.EXTERIOR_FACET): 437 integrals = ig[domain] 438 classname, code, includes = generate_exterior_facet_integral_code(integrals, formrep) 439 filenames.extend( write_code(classname, code) ) 440 needed_files.extend(includes) 441 442 # Generate code for interior facet integrals 443 for domain in pf.domains(Measure.INTERIOR_FACET): 444 integrals = ig[domain] 445 classname, code, includes = generate_interior_facet_integral_code(integrals, formrep) 446 filenames.extend( write_code(classname, code) ) 447 needed_files.extend(includes) 448 449 # Generate code for form! 450 code, includes = generate_form_code(formrep) 451 filenames.extend( write_code(formrep.classname, code) ) 452 needed_files.extend(includes) 453 454 # Collect classnames for use with dolfin wrappers 455 namespace = "" # "sfc::" 456 ufc_form_classname = namespace + formrep.classname 457 ufc_finite_element_classnames = [namespace + name for name in formrep.fe_names] 458 ufc_dof_map_classnames = [namespace + name for name in formrep.dm_names] 459 460 fn = UFCFormNames(formdata.name, formdata.coefficient_names, 461 ufc_form_classname, ufc_finite_element_classnames, ufc_dof_map_classnames) 462 formnames.append(fn) 463 464 # Get other needed files: 465 if needed_files: 466 raise NotImplementedError("FIXME: Implement fetching non-ufc-class files like DofPtv and quadrature rule files.") 467 468 filenames = list(unique(chain(filenames, needed_files))) 469 hfilenames = [f for f in filenames if f.endswith(".h")] 470 cfilenames = [f for f in filenames if f.endswith(".cpp")] 471 472 # Generate DOLFIN wrapper code 473 if options.code.dolfin_wrappers: 474 if not formnames: 475 print "NOT generating dolfin wrappers, missing forms!" # TODO: Generate FunctionSpaces for elements? 476 else: 477 prefix = options.code.prefix 478 header = dolfin_header_template % "\n".join('#include "%s"' % h for h in hfilenames) 479 sfc_info("Generating DOLFIN wrapper code, formnames are %s." % "\n".join(map(str,formnames))) 480 dolfin_code = generate_dolfin_code(prefix, header, formnames) 481 dolfin_filename = "%s.h" % prefix 482 write_file(dolfin_filename, dolfin_code) 483 hfilenames.append(dolfin_filename) 484 485 sfc_debug("Leaving generate_code") 486 return hfilenames, cfilenames
487