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  import sfc 
 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   
 45  from sfc.representation import GeometryRepresentation 
 46  from sfc.representation import ElementRepresentation 
 47  from sfc.representation import FormRepresentation 
 48   
 49  from sfc.codegeneration.formcg import FormCG 
 50  from sfc.codegeneration.dofmapcg import DofMapCG 
 51  from sfc.codegeneration.finiteelementcg import FiniteElementCG 
 52   
 53  # FIXME: Move these somewhere else 
54 -def integral_handlers(domaintype, options):
55 # Switch family of classes based on options 56 if options.code.integral.use_uflacs: 57 import sfc.representation.integralrep_uflacs as r 58 import sfc.codegeneration.integralcg_uflacs as c 59 else: 60 import sfc.representation.integralrep_regular as r 61 import sfc.codegeneration.integralcg_regular as c 62 63 # Select integral classes based on domain type 64 if domaintype == Measure.CELL: 65 return (r.CellIntegralRepresentation, 66 c.CellIntegralCG) 67 if domaintype == Measure.EXTERIOR_FACET: 68 return (r.ExteriorFacetIntegralRepresentation, 69 c.ExteriorFacetIntegralCG) 70 if domaintype == Measure.INTERIOR_FACET: 71 return (r.InteriorFacetIntegralRepresentation, 72 c.InteriorFacetIntegralCG)
73 74 75 # version %(version)s 76 _header_template = r"""/* 77 * %(name)s.h 78 * 79 * This file was automatically generated by SFC. 80 * 81 * http://www.fenics.org/syfi/ 82 * 83 */ 84 85 #ifndef GENERATED_SFC_CODE_%(name)s 86 #define GENERATED_SFC_CODE_%(name)s 87 88 #include "ufc.h" 89 90 %(hincludes)s 91 92 //namespace sfc 93 //{ 94 %(hbody)s 95 //} 96 97 #endif 98 """ 99 100 # version %(version)s 101 _implementation_template = r"""/* 102 * %(name)s.cpp 103 * 104 * This file was automatically generated by SFC. 105 * 106 * http://www.fenics.org/syfi/ 107 * 108 */ 109 110 #include "ufc.h" 111 112 #include "%(name)s.h" 113 114 %(cincludes)s 115 116 //namespace // Anonymous namespace 117 //{ 118 %(supportcode)s 119 //} 120 121 //namespace sfc 122 //{ 123 %(cbody)s 124 //} 125 """ 126
127 -def common_system_includes():
128 return ("iostream", "cmath", "stdexcept", "cstring")
129
130 -def apply_code_dict(format_string, code_dict):
131 "Template formatting with improved checking for template argument mismatch." 132 expected_keys = set(re.findall('%\((\w+)\)', format_string)) 133 actual_keys = set(code_dict.keys()) 134 135 missing_keys = expected_keys.difference(actual_keys) 136 if missing_keys: 137 print "When formatting code using template, missing keys:" 138 print "\n".join(sorted(missing_keys)) 139 print 140 141 #superfluous_keys = actual_keys.difference(expected_keys) 142 #if superfluous_keys: 143 # print "When formatting code using template, got superfluous keys:" 144 # print "\n".join(sorted(superfluous_keys)) 145 # print 146 147 return format_string % code_dict
148
149 -def generate_class_files(cg):
150 sfc_debug("Entering generate_class_files") 151 152 code_dict = cg.generate_code_dict() 153 classname = code_dict["classname"] 154 155 hbody = apply_code_dict(cg.header_template(), code_dict) 156 supportcode = cg.generate_support_code() 157 cbody = apply_code_dict(cg.impl_template(), code_dict) 158 159 # TODO: Probably don't need all these always, could speed up compilation a bit 160 system_includes = common_system_includes() 161 162 hincludes = "\n".join('#include "%s"' % inc for inc in cg.hincludes()) 163 164 cincludes = "\n".join('#include <%s>' % f for f in system_includes) 165 cincludes += "\n" 166 cincludes += "\n".join('#include "%s"' % inc 167 for inc in chain(cg.cincludes(), cg.local_includes())) 168 169 hdict = { "name": classname, 170 "hbody": hbody, 171 "hincludes": hincludes } 172 cdict = { "name": classname, 173 "supportcode": supportcode, 174 "cbody": cbody, 175 "cincludes": cincludes } 176 177 hcode = _header_template % hdict 178 ccode = _implementation_template % cdict 179 180 # External includes, i.e. ones not generated by sfc 181 includes = cg.hincludes() + cg.cincludes() 182 183 sfc_debug("Leaving generate_class_files") 184 return classname, (hcode, ccode), includes
185
186 -def write_file(filename, contents):
187 "Write contents to a file." 188 #print "WRITING FILE", filename 189 # TODO: Setup some kind of hook here for testing without writing files? 190 f = open(filename, "w") 191 f.write(contents) 192 f.close()
193
194 -def write_code(classname, code):
195 "Write code for a class to a header file or header + impl file." 196 sfc_debug("Entering write_code") 197 if isinstance(code, tuple): 198 # Code is split in a header and implementation file 199 hcode, ccode = code 200 hname = classname + ".h" 201 cname = classname + ".cpp" 202 write_file(hname, hcode) 203 write_file(cname, ccode) 204 ret = (hname, cname) 205 else: 206 # All code is combined in a header file 207 name = classname + ".h" 208 write_file(name, code) 209 ret = (name,) 210 sfc_debug("Leaving write_code") 211 return ret
212
213 -def compiler_input(input, objects=None, common_cell=None):
214 """Map different kinds of input to a list of 215 UFL elements and a list of FormData instances. 216 217 The following input formats are allowed: 218 - ufl.Form 219 - ufl.algorithms.FormData 220 - ufl.FiniteElementBase 221 - list of the above 222 223 Returns: 224 elements, formdatas 225 """ 226 sfc_debug("Entering compiler_input") 227 sfc_assert(input, "Got no input!") 228 fd = [] 229 fe = [] 230 if not isinstance(input, list): 231 input = [input] 232 for d in input: 233 if isinstance(d, ufl.form.Form): 234 d = d.compute_form_data(object_names=objects, common_cell=common_cell) 235 if isinstance(d, ufl.algorithms.formdata.FormData): 236 fd.append(d) 237 fe.extend(d.sub_elements) 238 elif isinstance(d, ufl.FiniteElementBase): 239 if d.cell().is_undefined(): 240 sfc_assert(common_cell and not common_cell.is_undefined(), 241 "Element no defined cell, cannot compile this.") 242 d = d.reconstruct(cell=common_cell) 243 sfc_assert(not d.cell().is_undefined(), "Still undefined?") 244 fe.append(d) 245 else: 246 sfc_error("Not a FormData or FiniteElementBase object: %s" % str(type(d))) 247 fe = sorted(set(fe)) 248 sfc_debug("Leaving compiler_input") 249 for x in fe: 250 sfc_assert(not x.cell().is_undefined(), 251 "Can never have undefined cells at this point!") 252 for x in fd: 253 sfc_assert(not x.cell.is_undefined(), 254 "Can never have undefined cells at this point!") 255 return fe, fd
256 257 dolfin_header_template = """/* 258 * DOLFIN wrapper code generated by the SyFi Form Compiler. 259 */ 260 261 %s 262 """
263 -def generate_dolfin_header_file(prefix, formnames, hfilenames):
264 if formnames: 265 sfc_info("Generating DOLFIN wrapper code, formnames are %s." %\ 266 "\n".join(map(str,formnames))) 267 268 header = dolfin_header_template % "\n".join('#include "%s"' % h 269 for h in hfilenames) 270 dolfin_code = generate_dolfin_code(prefix, header, formnames) 271 272 dolfin_filename = "%s.h" % prefix 273 write_file(dolfin_filename, dolfin_code) 274 return dolfin_filename 275 else: 276 # TODO: Generate FunctionSpaces for elements? 277 sfc_info("NOT generating dolfin wrappers, no forms to wrap!") 278 return None
279
280 -def update_global_syfi_nsd(nsd):
281 if nsd == update_global_syfi_nsd.last_nsd: 282 return # TODO: Check actual dimension? 283 else: 284 update_global_syfi_nsd.last_nsd = nsd 285 SyFi.initSyFi(nsd)
286 update_global_syfi_nsd.last_nsd = None 287
288 -def update_element_representations(e, geomrep, element_reps, generated_elements, options):
289 # Already visited this element? Skip it if possible. 290 if e in element_reps: 291 return [] 292 293 # Construct ElementRepresentation objects for all elements 294 quad_rule = None # FIXME: What to do with this one? 295 assert not "Quadrature" in e.family() # No quadrature rule defined! 296 297 # Create representation and cache for reuse 298 erep = ElementRepresentation(e, geomrep, quad_rule=quad_rule, 299 options=options, cache=element_reps) 300 element_reps[e] = erep 301 302 # Build flat list of all subelements 303 todo = [erep] 304 ereps = [] 305 while todo: 306 erep = todo.pop() 307 # Skip already generated subelements! 308 if not erep.ufl_element in generated_elements: 309 generated_elements.add(erep.ufl_element) 310 ereps.append(erep) 311 # Queue subelements for inspection 312 for subrep in erep.sub_elements: 313 todo.append(subrep) 314 315 return ereps
316
317 -def generate_element_ufc_code(element_reps_list, options):
318 # Generate code for all element representations 319 generated_classes = [] 320 for rep in element_reps_list: 321 # Generate code for finite_element 322 cg = FiniteElementCG(rep) 323 generated_classes.append(generate_class_files(cg)) 324 325 # Generate code for dof_map 326 cg = DofMapCG(rep) 327 generated_classes.append(generate_class_files(cg)) 328 return generated_classes
329
330 -def generate_form_ufc_code(formrep, options):
331 332 pf = formdata = formrep.formdata.preprocessed_form 333 ig = pf.integral_groups() 334 335 generated_classes = [] 336 337 # Generate code for all types of integral classes 338 for domaintype in (Measure.CELL, Measure.EXTERIOR_FACET, Measure.INTERIOR_FACET): 339 # Generate code for this particular integral type 340 IntegralRepresentation, IntegralCG = integral_handlers(domaintype, options) 341 342 for subdomain in pf.domains(domaintype): 343 integrals = ig[subdomain] 344 if integrals: 345 # Create symbolic representation for this integral 346 itgrep = IntegralRepresentation(integrals, formrep) 347 348 # Create code generation object for this integral class 349 cg = IntegralCG(itgrep) 350 351 # Generate code for this class 352 generated_classes.append(generate_class_files(cg)) 353 354 # Generate code for form class, which basically collects everything 355 cg = FormCG(formrep) 356 generated_classes.append(generate_class_files(cg)) 357 358 return generated_classes
359
360 -def create_geometry_rep(domain, geometry_reps):
361 nsd = ufl.geometry.domain2dim[domain] 362 363 # Initialize global variables in SyFi with the 364 # right space dimension (not very nice! need a 365 # significant refactoring of SyFi to fix this...) 366 update_global_syfi_nsd(nsd) 367 368 if not domain in geometry_reps: 369 geomrep = GeometryRepresentation(domain) 370 geometry_reps[domain] = geomrep 371 372 return geometry_reps[domain]
373
374 -def generate_code(input, objects, options=None, common_cell=None):
375 """Generate code from input and options. 376 377 @param input: 378 TODO 379 @param options: 380 TODO 381 """ 382 sfc_debug("Entering generate_code") 383 sfc_assert(input, "Got no input!") 384 385 if options is None: 386 options = default_parameters() 387 388 DEBUG = 0 389 if DEBUG: print "Calling compiler_input" 390 ufl_elements, formdatas = compiler_input(input, objects, common_cell=common_cell) 391 392 generated_classes = [] 393 geometry_reps = {} 394 element_reps = {} 395 generated_elements = set() 396 397 # Generate UFC code for elements 398 if DEBUG: print "Making element representations" 399 for e in ufl_elements: 400 # Create a shared geometry representation 401 # (sets the global space dimension in SyFi as a side effect) 402 geomrep = create_geometry_rep(e.cell().domain(), geometry_reps) 403 404 # Create symbolic representations for all new elements 405 element_reps_list = update_element_representations(e, 406 geomrep, 407 element_reps, 408 generated_elements, 409 options) 410 411 # Generate code for all new elements (finite_element and dof_map) 412 gc = generate_element_ufc_code(element_reps_list, options) 413 generated_classes.extend(gc) 414 415 # Generate UFC code for forms 416 formnames = [] 417 for formdata in formdatas: 418 if DEBUG: print "Starting processing another form." 419 420 # Create a shared geometry representation 421 # (sets the global space dimension in SyFi as a side effect) 422 geomrep = create_geometry_rep(formdata.cell.domain(), geometry_reps) 423 424 # Build symbolic representation 425 if DEBUG: print "Making form representation" 426 formrep = FormRepresentation(formdata, element_reps, geomrep, options) 427 428 # Generate code for form class 429 if DEBUG: print "Generating form code" 430 gc = generate_form_ufc_code(formrep, options) 431 generated_classes.extend(gc) 432 433 # Collect classnames for use with dolfin wrappers 434 namespace = "" # "sfc::" 435 fn = UFCFormNames(formdata.name, 436 formdata.coefficient_names, 437 namespace + formrep.classname, 438 [namespace + name for name in formrep.fe_names], 439 [namespace + name for name in formrep.dm_names]) 440 formnames.append(fn) 441 442 if DEBUG: print "Done with this form." 443 444 # Write code to files 445 filenames = [] 446 needed_files = [] 447 for classname, code, includes in generated_classes: 448 fn = write_code(classname, code) 449 filenames.extend(fn) 450 needed_files.extend(includes) 451 452 # Fetch other needed files: 453 if needed_files: 454 raise NotImplementedError("FIXME: Implement fetching non-ufc-class "+\ 455 "files like DofPtv and quadrature rule files.") 456 457 filenames = list(unique(chain(filenames, needed_files))) 458 hfilenames = [f for f in filenames if f.endswith(".h")] 459 cfilenames = [f for f in filenames if f.endswith(".cpp")] 460 461 # Generate DOLFIN wrapper code 462 if options.code.dolfin_wrappers: 463 hfn = generate_dolfin_header_file(options.code.prefix, formnames, hfilenames) 464 if hfn: hfilenames.append(hfn) 465 466 sfc_debug("Leaving generate_code") 467 return hfilenames, cfilenames
468