Skip to content
Snippets Groups Projects
Commit 3ac1629f authored by Marcel Koch's avatar Marcel Koch
Browse files

Reverts introduction of nested loops over basis functions

Computes tensor index from single index.
parent 8385c3b1
No related branches found
No related tags found
No related merge requests found
......@@ -11,8 +11,7 @@ from dune.perftool.blockstructured.geometry import (pymbolic_jacobian_inverse_tr
pymbolic_jacobian_determinant,
pymbolic_facet_jacobian_determinant,
to_global)
from dune.perftool.blockstructured.basis import (pymbolic_reference_gradient,
pymbolic_basis)
from dune.perftool.blockstructured.basis import pymbolic_reference_gradient
from dune.perftool.blockstructured.tools import sub_element_inames
from dune.perftool.pdelab import PDELabInterface
......@@ -30,17 +29,6 @@ class BlockStructuredInterface(PDELabInterface):
def lfs_inames(self, element, restriction, number=None, context=''):
return lfs_inames(element, restriction, number, context) + sub_element_inames()
#
# Test and trial function related generator functions
#
def pymbolic_basis(self, element, restriction, number):
return pymbolic_basis(element, restriction, number)
def pymbolic_reference_gradient(self, element, restriction, number):
return pymbolic_reference_gradient(element, restriction, number)
#
# Geometry related generator functions
#
......
......@@ -80,17 +80,6 @@ def evaluate_basis(leaf_element, name, restriction):
)
@backend(interface="pymbolic_basis", name="blockstructured")
def pymbolic_basis(leaf_element, restriction, number, context=''):
assert leaf_element.num_sub_elements() == 0
name = "phi_{}".format(FEM_name_mangling(leaf_element))
name = restricted_name(name, restriction)
evaluate_basis(leaf_element, name, restriction)
iname = lfs_inames(leaf_element, restriction, number, context=context)
return prim.Subscript(prim.Variable(name), (tensor_index_to_sequential_index(iname, 2), ))
@backend(interface="evaluate_grad", name="blockstructured")
@kernel_cached
def evaluate_reference_gradient(leaf_element, name, restriction):
......
......@@ -2,7 +2,7 @@ from dune.perftool.generation import (backend,
domain)
from dune.perftool.pdelab.geometry import world_dimension
#TODO subelem_inames hinzufugen?
@backend(interface="lfs_inames", name="blockstructured")
def lfs_inames(element, restriction, count=None, context=''):
assert not ((context == '') and (count is None))
......@@ -12,11 +12,7 @@ def lfs_inames(element, restriction, count=None, context=''):
else:
context = str(count)
inames = tuple()
dim = world_dimension()
for i in range(dim):
name = "micro_{}_index_{}".format(context, i)
inames = inames + (name,)
domain(name, 2)
name = "micro_{}_index".format(context)
domain(name, pow(2,world_dimension()))
return inames
return (name, )
......@@ -32,9 +32,13 @@ def tensor_index_to_sequential_index(inames, k):
return index
# TODO better name
def to_tensor_index(iname):
return tuple(prim.RightShift(prim.BitwiseAnd((prim.Variable(iname),2**n)), n) for n in range(world_dimension()))
# TODO 3d
def micro_index_to_macro_index(inames):
inames = tensor_index_to_sequential_index(inames, get_option("number_of_blocks")+1)
def micro_index_to_macro_index(iname):
subelem_inames_orig = sub_element_inames()
......@@ -83,7 +87,8 @@ def micro_index_to_macro_index(inames):
# if local_dimension() > 2:
# modified_index = prim.Sum((modified_index, prim.Variable(z)))
modified_index = prim.Sum((inames, tensor_index_to_sequential_index(sublem_inames, get_option("number_of_blocks")+1)))
modified_index = prim.Sum(tuple(to_tensor_index(iname)[n] * (k + 1) ** n for n in range(world_dimension()))
+ (tensor_index_to_sequential_index(sublem_inames, get_option("number_of_blocks") + 1), ))
return modified_index
......@@ -139,9 +139,6 @@ def pymbolic_coefficient(container, lfs, index):
if isinstance(lfs, str):
valuearg(lfs, dtype=NumpyType("str"))
if isinstance(index, tuple):
index = index[0]
# If the LFS is not yet a pymbolic expression, make it one
from pymbolic.primitives import Expression
if not isinstance(lfs, Expression):
......
......@@ -14,7 +14,6 @@ from dune.perftool.pdelab.spaces import (lfs_child,
name_lfs,
name_lfs_bound,
type_gfs,
lfs_inames
)
from dune.perftool.pdelab.geometry import (component_iname,
world_dimension,
......@@ -80,6 +79,7 @@ def declare_cache_temporary(element, restriction, which):
return decl
@backend(interface="evaluate_basis")
@kernel_cached
def evaluate_basis(leaf_element, name, restriction):
lfs = name_leaf_lfs(leaf_element, restriction)
......@@ -97,17 +97,17 @@ def evaluate_basis(leaf_element, name, restriction):
)
@backend(interface="pymbolic_basis")
def pymbolic_basis(leaf_element, restriction, number, context=''):
assert leaf_element.num_sub_elements() == 0
name = "phi_{}".format(FEM_name_mangling(leaf_element))
name = restricted_name(name, restriction)
evaluate_basis(leaf_element, name, restriction)
iname = lfs_inames(leaf_element, restriction, number, context=context)[0]
get_backend("evaluate_basis", selector=option_switch("blockstructured"))(leaf_element, name, restriction)
iname = get_backend("lfs_inames", selector=option_switch("blockstructured"))(leaf_element, restriction, number, context=context)[0]
return Subscript(Variable(name), (Variable(iname), ))
@backend(interface="evaluate_grad")
@kernel_cached
def evaluate_reference_gradient(leaf_element, name, restriction):
lfs = name_leaf_lfs(leaf_element, restriction)
......@@ -125,13 +125,12 @@ def evaluate_reference_gradient(leaf_element, name, restriction):
)
@backend(interface="pymbolic_reference_gradient")
def pymbolic_reference_gradient(leaf_element, restriction, number, context=''):
assert leaf_element.num_sub_elements() == 0
name = "js_{}".format(FEM_name_mangling(leaf_element))
name = restricted_name(name, restriction)
evaluate_reference_gradient(leaf_element, name, restriction)
iname = lfs_inames(leaf_element, restriction, number, context=context)[0]
get_backend("evaluate_grad", selector=option_switch("blockstructured"))(leaf_element, name, restriction)
iname = get_backend("lfs_inames", selector=option_switch("blockstructured"))(leaf_element, restriction, number, context=context)[0]
return Subscript(Variable(name), (Variable(iname), 0))
......@@ -200,7 +199,7 @@ def evaluate_coefficient_gradient(element, name, container, restriction, tree_pa
# and proceed to call the necessary generator functions
temporary_variable(name, shape=shape, shape_impl=shape_impl)
lfs = name_lfs(element, restriction, tree_path)
basis = get_backend("pymbolic_reference_gradient", selector=option_switch("blockstructured"))(leaf_element, restriction, 0, context='trialgrad')
basis = pymbolic_reference_gradient(leaf_element, restriction, 0, context='trialgrad')
index, _ = get_pymbolic_indices(basis)
if isinstance(sub_element, (VectorElement, TensorElement)):
......
......@@ -200,6 +200,7 @@ def name_in_cell_geometry(restriction):
return name
# TODO is it always necessary to add quadrature inames?
def apply_in_cell_transformation(name, local, restriction):
geo = name_in_cell_geometry(restriction)
return quadrature_preamble("{} = {}.global({});".format(name,
......@@ -420,7 +421,7 @@ def apply_to_global_transformation(name, local):
depends_on=frozenset({Writes(get_pymbolic_basename(local))})
)
#TODO check backend for consistency
#TODO check backend & function name for consistency
@backend(interface="to_global", name="default")
def to_global(local):
assert isinstance(local, prim.Expression)
......
......@@ -243,12 +243,12 @@ def determine_accumulation_space(expr, number, measure, idims=None):
lfs = lfs_child(lfs, idims, shape=subel.value_shape(), symmetry=subel.symmetry())
subel = subel.sub_elements()[0]
iname = get_backend("lfs_inames", selector=option_switch("blockstructured"))(subel, ma.restriction, count=number)
iname, = get_backend("lfs_inames", selector=option_switch("blockstructured"))(subel, ma.restriction, count=number)
if get_option("blockstructured"):
from dune.perftool.blockstructured.tools import micro_index_to_macro_index
lfsi = micro_index_to_macro_index(iname)
else:
lfsi = Variable(iname[0])
lfsi = Variable(iname)
# If the LFS is not yet a pymbolic expression, make it one
......
......@@ -23,10 +23,8 @@ def get_pymbolic_indices(expr):
def ind(i):
if isinstance(i, int):
return (i,)
if isinstance(i, prim.Sum) or isinstance(i, prim.Product):
return sum((ind(child) for child in i.children if isinstance(child, prim.Expression)), ())
return (get_pymbolic_basename(i),)
return i
return get_pymbolic_basename(i)
if not isinstance(expr.index, tuple):
return (ind(expr.index),)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment