← Back to branch summary

~fenics-core/ffc/spatder2grad

« back to all changes in this revision

Viewing changes to ffc/cpp.py

  • Committer: Marie E. Rognes
  • Date: 2013-01-08 12:10:22 UTC
  • mfrom: (1760.1.82 ffc-manifolds-trunk)
  • Revision ID: meg@simula.no-20130108121022-svq8jno62hvwp8eu
Merge work on code generation for manifolds. You can now define
elements and forms over cells with different geometric and topological
dimensions.

Show diffs side-by-side

added added

removed removed

Lines of Context:
ffc/cpp.py
179
179
    "local dof":                "dof",
180
180
    "basisvalues":              "basisvalues",
181
181
    "coefficients":             lambda i: "coefficients%d" %(i),
182
 
    "num derivatives":          "num_derivatives",
183
 
    "derivative combinations":  "combinations",
 
182
    "num derivatives":          lambda t_or_g :"num_derivatives" + t_or_g,
 
183
    "derivative combinations":  lambda t_or_g :"combinations" + t_or_g,
184
184
    "transform matrix":         "transform",
185
185
    "transform Jinv":           "Jinv",
186
186
    "dmats":                    lambda i: "dmats%s" %(i),
228
228
 
229
229
# Code snippets
230
230
from codesnippets import *
 
231
 
231
232
format.update({
232
233
    "cell coordinates":     cell_coordinates,
233
 
    "jacobian":             lambda n, r="": jacobian[n] % {"restriction": r},
234
 
    "inverse jacobian":     lambda n, r="": inverse_jacobian[n] % {"restriction": r},
235
 
    "jacobian and inverse": lambda n, r=None: format["jacobian"](n, choose_map[r]) +\
236
 
                            "\n" + format["inverse jacobian"](n, choose_map[r]),
237
 
    "facet determinant":    lambda n, r=None: facet_determinant[n] % {"restriction": choose_map[r]},
238
 
    "fiat coordinate map":  lambda n: fiat_coordinate_map[n],
239
 
    "generate normal":      lambda d, i: _generate_normal(d, i),
240
 
    "generate cell volume": lambda d, i: _generate_cell_volume(d, i),
241
 
    "generate circumradius": lambda d, i: _generate_circumradius(d, i),
242
 
    "generate facet area":  lambda d: facet_area[d],
 
234
    "jacobian":             lambda gdim, tdim, r="": jacobian[gdim][tdim] % {"restriction": r},
 
235
    "inverse jacobian":     lambda gdim, tdim, r="", oriented=False: inverse_jacobian[gdim][tdim] % {"restriction": r}  + format["orientation"](oriented, gdim, tdim, r),
 
236
    "jacobian and inverse": lambda gdim, tdim=None, r=None, oriented=False: (format["jacobian"](gdim, tdim, choose_map[r]) + \
 
237
                                                                                 "\n" + format["inverse jacobian"](gdim, tdim, choose_map[r], oriented)),
 
238
    "orientation":          lambda oriented, gdim, tdim, r="": orientation_snippet % {"restriction": r} if (oriented and gdim != tdim) else "",
 
239
    "facet determinant":    lambda gdim, tdim, r=None: facet_determinant[gdim][tdim] % {"restriction": choose_map[r]},
 
240
    "fiat coordinate map":  lambda cell, gdim: fiat_coordinate_map[cell][gdim],
 
241
    "generate normal":      lambda gdim, tdim, i: _generate_normal(gdim, tdim, i),
 
242
    "generate cell volume": lambda gdim, tdim, i: _generate_cell_volume(gdim, tdim, i),
 
243
    "generate circumradius": lambda gdim, tdim, i: _generate_circumradius(gdim, tdim, i),
 
244
    "generate facet area":  lambda gdim, tdim: facet_area[gdim][tdim],
243
245
    "generate ip coordinates":  lambda g, num_ip, name, ip, r=None: (ip_coordinates[g][0], ip_coordinates[g][1] % \
244
246
                                {"restriction": choose_map[r], "ip": ip, "name": name, "num_ip": num_ip}),
245
247
    "scale factor snippet": scale_factor,
529
531
 
530
532
    return name
531
533
 
532
 
def _generate_jacobian(cell_dimension, integral_type):
533
 
    "Generate code for computing jacobian"
534
 
 
535
 
    # Choose space dimension
536
 
    if cell_dimension == 1:
537
 
        jacobian = jacobian_1D
538
 
        facet_determinant = facet_determinant_1D
539
 
    elif cell_dimension == 2:
540
 
        jacobian = jacobian_2D
541
 
        facet_determinant = facet_determinant_2D
542
 
    else:
543
 
        jacobian = jacobian_3D
544
 
        facet_determinant = facet_determinant_3D
545
 
 
546
 
    # Check if we need to compute more than one Jacobian
547
 
    if integral_type == "cell":
548
 
        code  = jacobian % {"restriction":  ""}
549
 
        code += "\n\n"
550
 
        code += scale_factor
551
 
    elif integral_type == "exterior facet":
552
 
        code  = jacobian % {"restriction":  ""}
553
 
        code += "\n\n"
554
 
        code += facet_determinant % {"restriction": "", "facet" : "facet"}
555
 
    elif integral_type == "interior facet":
556
 
        code  = jacobian % {"restriction": choose_map["+"]}
557
 
        code += "\n\n"
558
 
        code += jacobian % {"restriction": choose_map["-"]}
559
 
        code += "\n\n"
560
 
        code += facet_determinant % {"restriction": choose_map["+"], "facet": "facet0"}
561
 
 
562
 
    return code
563
 
 
564
 
def _generate_normal(geometric_dimension, domain_type, reference_normal=False):
 
534
def _generate_normal(geometric_dimension, topological_dimension, domain_type,
 
535
                     reference_normal=False):
565
536
    "Generate code for computing normal"
566
537
 
567
538
    # Choose snippets
568
 
    direction = normal_direction[geometric_dimension]
569
 
    normal = facet_normal[geometric_dimension]
 
539
    direction = normal_direction[geometric_dimension][topological_dimension]
 
540
 
 
541
    assert (facet_normal[geometric_dimension].has_key(topological_dimension)),\
 
542
        "Facet normal not yet implemented for this gdim/tdim combo"
 
543
    normal = facet_normal[geometric_dimension][topological_dimension]
570
544
 
571
545
    # Choose restrictions
572
546
    if domain_type == "exterior_facet":
580
554
        error("Unsupported domain_type: %s" % str(domain_type))
581
555
    return code
582
556
 
583
 
def _generate_cell_volume(geometric_dimension, domain_type):
 
557
def _generate_cell_volume(gdim, tdim, domain_type):
584
558
    "Generate code for computing cell volume."
585
559
 
586
560
    # Choose snippets
587
 
    volume = cell_volume[geometric_dimension]
 
561
    volume = cell_volume[gdim][tdim]
588
562
 
589
563
    # Choose restrictions
590
564
    if domain_type in ("cell", "exterior_facet"):
596
570
        error("Unsupported domain_type: %s" % str(domain_type))
597
571
    return code
598
572
 
599
 
def _generate_circumradius(geometric_dimension, domain_type):
 
573
def _generate_circumradius(gdim, tdim, domain_type):
600
574
    "Generate code for computing a cell's circumradius."
601
575
 
602
576
    # Choose snippets
603
 
    radius = circumradius[geometric_dimension]
 
577
    radius = circumradius[gdim][tdim]
604
578
 
605
579
    # Choose restrictions
606
580
    if domain_type in ("cell", "exterior_facet"):