diff --git a/python/dune/perftool/blockstructured/geometry.py b/python/dune/perftool/blockstructured/geometry.py index ce77c3d095316a71550bb46c6b3a7da243426869..63bbaa4153f81c7f2cd3dd22a7d61ad2d96b96f4 100644 --- a/python/dune/perftool/blockstructured/geometry.py +++ b/python/dune/perftool/blockstructured/geometry.py @@ -320,22 +320,30 @@ def apply_default_to_global_transformation(name, local): if isinstance(local, str): local = prim.Variable(local) - d_var = prim.Variable(component_iname('to_global')) + dim_pym = prim.Variable(component_iname('to_global')) - coeffs = bilinear_transformation_coefficients() + coeffs_pym = [prim.Subscript(prim.Variable(coeff), (dim_pym,)) for coeff in bilinear_transformation_coefficients()] + + local_pym = [prim.Subscript(local, i) for i in range(dim)] + + corner_0_pym = prim.Subscript(prim.Variable(corners), (0, dim_pym)) # global[d] = T(local)[d] if dim == 2: - a_var, b_var, c_var = [prim.Variable(coeff) for coeff in coeffs] - expr = (prim.Subscript(a_var, (d_var)) * prim.Subscript(local, (0,)) * prim.Subscript(local, (1,)) + - prim.Subscript(b_var, (d_var,)) * prim.Subscript(local, (0,)) + - prim.Subscript(c_var, (d_var,)) * prim.Subscript(local, (1,)) + - prim.Subscript(prim.Variable(corners), (0, d_var))) + a_pym, b_pym, c_pym = coeffs_pym + expr = a_pym * local_pym[0] * local_pym[1] + b_pym * local_pym[0] + c_pym * local_pym[1] + corner_0_pym + elif dim == 3: + a_pym, b_pym, c_pym, d_pym, e_pym, f_pym, g_pym = coeffs_pym + expr = a_pym * local_pym[0] * local_pym[1] * local_pym[2] + b_pym * local_pym[0] * local_pym[1] + \ + c_pym * local_pym[0] * local_pym[2] + d_pym * local_pym[1] * local_pym[2] + \ + e_pym * local_pym[0] + f_pym * local_pym[1] + g_pym * local_pym[2] + corner_0_pym + else: + raise NotImplementedError - assignee = prim.Subscript(prim.Variable(name), (d_var,)) + assignee = prim.Subscript(prim.Variable(name), (dim_pym,)) instruction(assignee=assignee, expression=expr, - within_inames=frozenset(sub_element_inames() + get_backend(interface="quad_inames")() + (d_var.name,)), + within_inames=frozenset(sub_element_inames() + get_backend(interface="quad_inames")() + (dim_pym.name,)), within_inames_is_final=True, depends_on=frozenset({Writes(local.name), Writes(corners)}) ) @@ -349,19 +357,16 @@ def apply_constant_to_global_transformation(name, local): if isinstance(local, str): local = prim.Variable(local) - d = component_iname('to_global') + dim_pym = prim.Variable(component_iname('to_global')) # global[d] = lower_left[d] + local[d] * (upper_right[d] - lower_left[d]) - expr = prim.Sum((prim.Subscript(prim.Variable(corners), (0,prim.Variable(d))), - prim.Product((prim.Subscript(local, (prim.Variable(d),)), - prim.Sum((prim.Subscript(prim.Variable(corners), (2**dim - 1, prim.Variable(d))), - -1 * prim.Subscript(prim.Variable(corners), (0, prim.Variable(d))))) - )) - )) - assignee = prim.Subscript(prim.Variable(name), (prim.Variable(d),)) + expr = (prim.Subscript(prim.Variable(corners), (0, dim_pym)) + + (prim.Subscript(local, (dim_pym,)) * (prim.Subscript(prim.Variable(corners), (2**dim - 1, dim_pym)) - + prim.Subscript(prim.Variable(corners), (0, dim_pym))))) + assignee = prim.Subscript(prim.Variable(name), (dim_pym,)) instruction(assignee=assignee, expression=expr, - within_inames=frozenset(sub_element_inames() + get_backend(interface="quad_inames")() + (d,)), + within_inames=frozenset(sub_element_inames() + get_backend(interface="quad_inames")() + (dim_pym.name,)), within_inames_is_final=True, depends_on=frozenset({Writes(local.name), Writes(corners)}) )