Which ones need to be cleaned up? I'll be working on getting tabulate_tensor working. Let me know if I should help with the cleaning up.
-- Anders
--- Begin Message --------------------------------------------------------------- revno: 1412 committer: Marie E. Rognes <m...@simula.no> branch nick: ffc-unstable timestamp: Thu 2010-01-07 14:16:49 +0100 message: Implemented interpolate_vertex_values for scalar CG_1. Time to do some clean-ups with regard to code generation utilities before continuing. modified: ffc/codegeneration.py ffc/fiatinterface.py ffc/representation.py -- lp:~ffc-core/ffc/dev https://code.launchpad.net/~ffc-core/ffc/dev You are subscribed to branch lp:~ffc-core/ffc/dev. To unsubscribe from this branch go to https://code.launchpad.net/~ffc-core/ffc/dev/+edit-subscription.=== modified file 'ffc/codegeneration.py' --- ffc/codegeneration.py 2010-01-06 18:32:59 +0000 +++ ffc/codegeneration.py 2010-01-07 13:16:49 +0000 @@ -11,7 +11,7 @@ __copyright__ = "Copyright (C) 2009 " + __author__ __license__ = "GNU GPL version 3 or any later version" -# Last changed: 2010-01-06 +# Last changed: 2010-01-07 # FFC modules from ffc.log import begin, end, debug_code @@ -20,6 +20,7 @@ # FFC code generation modules from ffc.evaluatebasis import _evaluate_basis from ffc.evaluatedof import _evaluate_dof, _evaluate_dofs, affine_weights +from ffc.codesnippets import jacobian # FFC specialized code generation modules #from ffc.quadrature import generate_quadrature_integrals @@ -98,7 +99,7 @@ code["evaluate_basis_derivatives_all"] = not_implemented code["evaluate_dof"] = _evaluate_dof(ir["evaluate_dof"]) code["evaluate_dofs"] = _evaluate_dofs(ir["evaluate_dofs"]) - code["interpolate_vertex_values"] = not_implemented + code["interpolate_vertex_values"] = _interpolate_vertex_values(ir["interpolate_vertex_values"]) code["num_sub_elements"] = ret(ir["num_sub_elements"]) code["create_sub_element"] = _create_foo(ir["create_sub_element"], prefix, "finite_element") @@ -264,31 +265,64 @@ def _tabulate_coordinates(ir): + # Raise error if tabulate_coordinates is ill-defined if ir is None: return "// Raise error here" # Extract formats: - add = format["add"] - multiply = format["multiply"] + add, multiply = format["add"], format["multiply"] precision = format["float"] # Extract coordinates and cell dimension - coordinates = ir - cell_dim = len(coordinates[0]) + cell_dim = len(ir[0]) # Aid mapping points from reference to physical element coefficients = affine_weights(cell_dim) + # Generate code for each point and each component code = ["const double * const * x = c.coordinates;"] - - for (i, c) in enumerate(coordinates): - w = coefficients(c) + for (i, coordinate) in enumerate(ir): + w = coefficients(coordinate) for j in range(cell_dim): value = add([multiply([precision(w[k]), "x[%d][%d]" % (k, j)]) for k in range(cell_dim + 1)]) code += ["coordinates[%d][%d] = %s;" % (i, j, value)] return "\n".join(code) +def _interpolate_vertex_values(ir): + + code = [] + + # Add code for Jacobian if necessary + if ir["needs_jacobian"]: + code += [jacobian(ir["cell_dim"])] + + # Extract formats + add, multiply = format["add"], format["multiply"] + precision = format["float"] + + # Iterate over the elements + for (k, data) in enumerate(ir["element_data"]): + + # Extract vertex values for all basis functions + vertex_values = data["basis_values"] + + # Create code for each vertex + for j in range(len(vertex_values)): + + # Map basis values according to mapping + vertex_value = [precision(v) for v in vertex_values[j]] + + # Contract basis values and coefficients + value = add([multiply(["dof_values[%d]" % i, v]) + for (i, v) in enumerate(vertex_value)]) + + # Construct correct vertex value label + name = "vertex_values[%d]" % j + code += [name + " = " + value + ";"] + + return "\n".join(code) + #--- Utility functioins --- def _create_foo(numbers, prefix, class_name): === modified file 'ffc/fiatinterface.py' --- ffc/fiatinterface.py 2010-01-04 12:48:48 +0000 +++ ffc/fiatinterface.py 2010-01-07 13:16:49 +0000 @@ -5,7 +5,7 @@ # Modified by Garth N. Wells, 2009. # Modified by Marie Rognes, 2009-2010. -# Last changed: 2010-01-04 +# Last changed: 2010-01-07 # UFL modules from ufl import FiniteElement as UFLFiniteElement @@ -44,6 +44,11 @@ # Mapping from dimension to number of mesh sub-entities: entities_per_dim = {1: [2, 1], 2: [3, 3, 1], 3: [4, 6, 4, 1]} +def reference_cell(dim): + if isinstance(dim, int): + return ufc_simplex(dim) + else: + return ufc_simplex(domain2dim[dim]) def create_element(ufl_element): @@ -70,7 +75,7 @@ "Create FIAT element corresponding to given UFL finite element." ElementClass = family2class[ufl_element.family()] - reference_cell = ufc_simplex(domain2dim[ufl_element.cell().domain()]) - element = ElementClass(reference_cell, ufl_element.degree()) + cell = reference_cell(ufl_element.cell().domain()) + element = ElementClass(cell, ufl_element.degree()) return element === modified file 'ffc/representation.py' --- ffc/representation.py 2010-01-06 18:32:59 +0000 +++ ffc/representation.py 2010-01-07 13:16:49 +0000 @@ -28,7 +28,7 @@ # FFC modules from ffc.utils import compute_permutations from ffc.log import info, error, begin, end, debug_ir, ffc_assert -from ffc.fiatinterface import create_element, entities_per_dim +from ffc.fiatinterface import create_element, entities_per_dim, reference_cell from ffc.mixedelement import MixedElement @@ -63,11 +63,12 @@ # Create FIAT element element = create_element(ufl_element) + cell = ufl_element.cell() # Compute data for each function ir = {} ir["signature"] = repr(ufl_element) - ir["cell_shape"] = ufl_element.cell().domain() + ir["cell_shape"] = cell.domain() ir["space_dimension"] = element.space_dimension() ir["value_rank"] = len(ufl_element.value_shape()) ir["value_dimension"] = ufl_element.value_shape() @@ -75,11 +76,12 @@ ir["evaluate_basis_all"] = not_implemented #element.get_coeffs() ir["evaluate_basis_derivatives"] = element ir["evaluate_basis_derivatives_all"] = not_implemented #element.get_coeffs() - ir["evaluate_dof"] = _evaluate_dof(element, ufl_element.cell()) - ir["evaluate_dofs"] = None - ir["interpolate_vertex_values"] = None + ir["evaluate_dof"] = _evaluate_dof(element, cell) + ir["evaluate_dofs"] = not_implemented + ir["interpolate_vertex_values"] = _interpolate_vertex_values(element, cell) ir["num_sub_elements"] = ufl_element.num_sub_elements() - ir["create_sub_element"] = [form_data.element_map[e] for e in ufl_element.sub_elements()] + ir["create_sub_element"] = [form_data.element_map[e] + for e in ufl_element.sub_elements()] debug_ir(ir, "finite_element") @@ -116,7 +118,8 @@ ir["tabulate_entity_dofs"] = not_implemented ir["tabulate_coordinates"] = _tabulate_coordinates(element) ir["num_sub_dof_maps"] = ufl_element.num_sub_elements() - ir["create_sub_dof_map"] = [form_data.element_map[e] for e in ufl_element.sub_elements()] + ir["create_sub_dof_map"] = [form_data.element_map[e] + for e in ufl_element.sub_elements()] debug_ir(ir, "dofmap") @@ -202,24 +205,18 @@ def _evaluate_dof(element, cell): "Compute intermediate representation of evaluate_dof." - # Setup ir - ir = {"mappings": element.mapping(), - "value_dim": _value_dimension(element), - "cell_dimension": cell.geometric_dimension(), - "dofs": [L.pt_dict for L in element.dual_basis()]} - # Generate offsets: i.e value offset for each basis function - if not isinstance(element, MixedElement): - ir["offsets"] = [0]*element.space_dimension() - else: - offsets = [] - offset = 0 - for e in element.elements(): - offsets += [offset]*e.space_dimension() - offset += _value_dimension(e) - ir["offsets"] = offsets + offsets = [] + offset = 0 + for e in all_elements(element): + offsets += [offset]*e.space_dimension() + offset += _value_dimension(e) - return ir + return {"mappings": element.mapping(), + "value_dim": _value_dimension(element), + "cell_dimension": cell.geometric_dimension(), + "dofs": [L.pt_dict for L in element.dual_basis()], + "offsets": offsets} def _tabulate_coordinates(element): "Compute intermediate representation of tabulate_coordinates." @@ -230,16 +227,10 @@ def _tabulate_dofs(element, cell): "Compute intermediate representation of tabulate_dofs." - if isinstance(element, MixedElement): - ir = [{"entites_per_dim": entities_per_dim[cell.geometric_dimension()], + elements = all_elements(element) + return [{"entites_per_dim": entities_per_dim[cell.geometric_dimension()], "num_dofs_per_entity": _num_dofs_per_entity(e)} - for e in element.elements()] - else: - ir = [{"entites_per_dim": entities_per_dim[cell.geometric_dimension()], - "num_dofs_per_entity": _num_dofs_per_entity(element)}] - - return ir - + for e in elements] def _tabulate_facet_dofs(element, cell): "Compute intermediate representation of tabulate_facet_dofs." @@ -271,8 +262,39 @@ return facet_dofs +def _interpolate_vertex_values(element, cell): + + ir = {} + ir["cell_dim"] = cell.geometric_dimension() + + # Check whether computing the Jacobian is necessary + mappings = element.mapping() + ir["needs_jacobian"] = any("piola" in m for m in mappings) + + # Get vertices of reference cell + vertices = reference_cell(cell.domain()).get_vertices() + + # Compute data for each constituent element + extract = lambda values: [a for a in values[(0, 0)]] + + + ir["element_data"] = [{"space_dim": e.space_dimension(), + "value_dim": _value_dimension(e), + "basis_values": extract(e.tabulate(0, vertices)), + "mapping": e.mapping()} + for e in all_elements(element)] + + return ir + + #--- Utility functions --- +def all_elements(element): + try: + return element.elements() + except: + return [element] + def _num_dofs_per_entity(element): """ Compute list of integers representing the number of dofs
--- End Message ---
signature.asc
Description: Digital signature
_______________________________________________ Mailing list: https://launchpad.net/~ffc Post to : ffc@lists.launchpad.net Unsubscribe : https://launchpad.net/~ffc More help : https://help.launchpad.net/ListHelp
