Skip to content
Snippets Groups Projects
Commit d081faae authored by René Heß's avatar René Heß
Browse files

[skip ci] Very close to working poisson example

[ ] Go through system part -> GFS
[ ] Export boundary grid function
parent 4f1eaa30
No related branches found
No related tags found
No related merge requests found
......@@ -8,22 +8,22 @@ from dune.codegen.generation import (class_basename,
from dune.codegen.pdelab.driver import (get_form_ident,
name_initree,
)
from dune.codegen.pdelab.driver.gridfunctionspace import (name_leafview,
type_leafview,
)
@template_parameter(classtag="driver_block")
def driver_block_template_parameter():
from dune.codegen.pdelab.driver.gridfunctionspace import type_leafview
return type_leafview()
@class_member(classtag="driver_block")
def driver_block_grid_view():
from dune.codegen.pdelab.driver.gridfunctionspace import name_leafview, type_leafview
return "{} {};".format(type_leafview(), name_leafview())
def init_driver_block():
from dune.codegen.pdelab.driver.gridfunctionspace import name_leafview, type_leafview
driver_block_template_parameter()
gridview_argument = "_{}".format(name_leafview())
constructor_parameter("{}&".format(type_leafview()), gridview_argument, classtag="driver_block")
......@@ -39,6 +39,8 @@ def driver_block_basename(name):
def type_driver_block():
from dune.codegen.pdelab.driver.gridfunctionspace import type_leafview
# palpo TODO: driver_ident!
form_ident = get_form_ident()
......@@ -47,6 +49,7 @@ def type_driver_block():
@preamble(section="solver", kernel="main")
def define_driver_block(name):
from dune.codegen.pdelab.driver.gridfunctionspace import name_leafview
driver_block_type = type_driver_block()
return "{} {}({}, {});".format(driver_block_type, name, name_leafview(), name_initree())
......
......@@ -9,8 +9,8 @@ from dune.codegen.pdelab.driver import (get_form_ident,
get_trial_element,
preprocess_leaf_data,
)
from dune.codegen.pdelab.driver.gridfunctionspace import (name_leafview,
name_trial_gfs,
from dune.codegen.pdelab.driver.gridfunctionspace import (main_type_trial_gfs,
name_leafview,
name_trial_subgfs,
type_range,
)
......@@ -19,7 +19,8 @@ from dune.codegen.pdelab.driver.interpolate import (interpolate_vector,
)
from dune.codegen.pdelab.driver.solve import (define_vector,
dune_solve,
name_vector,
main_name_vector,
main_type_vector,
)
from ufl import MixedElement, TensorElement, VectorElement
......@@ -52,8 +53,11 @@ def type_discrete_grid_function(gfs):
@preamble(section="error", kernel="main")
def define_discrete_grid_function(gfs, vector_name, dgf_name):
dgf_type = type_discrete_grid_function(gfs)
return ["using {} = Dune::PDELab::DiscreteGridFunction<decltype({}),decltype({})>;".format(dgf_type, gfs, vector_name),
"{} {}({},{});".format(dgf_type, dgf_name, gfs, vector_name)]
gfs_type = main_type_trial_gfs()
form_ident = get_form_ident()
vector_type = main_type_vector(form_ident)
return ["using {} = Dune::PDELab::DiscreteGridFunction<{}, {}>;".format(dgf_type, gfs_type, vector_type),
"{} {}(*{},*{});".format(dgf_type, dgf_name, gfs, vector_name)]
def name_discrete_grid_function(gfs, vector_name):
......@@ -65,11 +69,11 @@ def name_discrete_grid_function(gfs, vector_name):
@preamble(section="error", kernel="main")
def typedef_difference_squared_adapter(name, treepath):
sol = name_exact_solution_gridfunction(treepath)
vector = name_vector(get_form_ident())
vector = main_name_vector(get_form_ident())
gfs = name_trial_subgfs(treepath)
dgf = name_discrete_grid_function(gfs, vector)
return 'using {} = Dune::PDELab::DifferenceSquaredAdapter<decltype({}), decltype({})>;'.format(name, sol, dgf)
dgf_type = type_discrete_grid_function(gfs)
return 'using {} = Dune::PDELab::DifferenceSquaredAdapter<decltype({}), {}>;'.format(name, sol, dgf_type)
def type_difference_squared_adapter(treepath):
......@@ -82,7 +86,7 @@ def type_difference_squared_adapter(treepath):
def define_difference_squared_adapter(name, treepath):
t = type_difference_squared_adapter(treepath)
sol = name_exact_solution_gridfunction(treepath)
vector = name_vector(get_form_ident())
vector = main_name_vector(get_form_ident())
gfs = name_trial_subgfs(treepath)
dgf = name_discrete_grid_function(gfs, vector)
......
......@@ -17,6 +17,8 @@ from dune.codegen.pdelab.driver import (FEM_name_mangling,
name_initree,
preprocess_leaf_data,
)
from dune.codegen.pdelab.driver.driverblock import (name_driver_block,
type_driver_block,)
from dune.codegen.loopy.target import type_floatingpoint
from ufl import FiniteElement, MixedElement, TensorElement, VectorElement, TensorProductElement, TensorProductCell
......@@ -254,32 +256,65 @@ def name_trial_gfs():
return name_gfs(element, is_dirichlet)
def main_name_trial_gfs():
element = get_trial_element()
is_dirichlet = preprocess_leaf_data(element, "is_dirichlet")
return name_gfs(element, is_dirichlet, main=True)
def name_test_gfs():
element = get_test_element()
is_dirichlet = preprocess_leaf_data(element, "is_dirichlet")
return name_gfs(element, is_dirichlet)
def name_gfs(element, is_dirichlet, treepath=(), root=True):
def name_gfs(element, is_dirichlet, treepath=(), root=True, main=False):
"""Generate name of grid function space
This function will call itself recursively to build the grid function space
tree.
Parameters
----------
element : UFL FiniteElement
is_dirichlet : Tuple
Which parts of the treepath use dirichlet constraints?
treepath :
Treepath for the grid function space tree
root : bool
Called for the root of the tree?
main : bool
Can be called for driver block generation or for generation from the
main of the program. Instead of copying the code to a new method
main_name_gfs we use switches to avoid code duplication.
"""
if isinstance(element, (VectorElement, TensorElement)):
subel = element.sub_elements()[0]
subgfs = name_gfs(subel, is_dirichlet[:subel.value_size()], treepath=treepath + (0,), root=False)
subgfs = name_gfs(subel, is_dirichlet[:subel.value_size()], treepath=treepath + (0,), root=False, main=main)
name = "{}_pow{}gfs_{}".format(subgfs,
element.num_sub_elements(),
"_".join(str(t) for t in treepath))
define_power_gfs(element, is_dirichlet, name, subgfs, root)
if main:
# TODO
assert False
else:
define_power_gfs(element, is_dirichlet, name, subgfs, root)
return name
elif isinstance(element, MixedElement):
k = 0
subgfs = []
for i, subel in enumerate(element.sub_elements()):
subgfs.append(name_gfs(subel, is_dirichlet[k:k + subel.value_size()], treepath=treepath + (i,), root=False))
subgfs.append(name_gfs(subel, is_dirichlet[k:k + subel.value_size()], treepath=treepath + (i,), root=False, main=main))
k = k + subel.value_size()
name = "_".join(subgfs)
if len(subgfs) == 1:
name = "{}_dummy".format(name)
name = "{}_{}".format(name, "_".join(str(t) for t in treepath))
define_composite_gfs(element, is_dirichlet, name, tuple(subgfs), root)
if main:
# TODO
assert False
else:
define_composite_gfs(element, is_dirichlet, name, tuple(subgfs), root)
return name
else:
assert isinstance(element, (FiniteElement, TensorProductElement))
......@@ -287,7 +322,10 @@ def name_gfs(element, is_dirichlet, treepath=(), root=True):
"_dirichlet" if is_dirichlet[0] else "",
"_".join(str(t) for t in treepath))
define_gfs(element, is_dirichlet, name, root)
if main:
main_define_gfs(element, is_dirichlet, name, root)
else:
define_gfs(element, is_dirichlet, name, root)
return name
......@@ -460,11 +498,17 @@ def type_constraintsassembler(is_dirichlet):
return name
def name_trial_subgfs(treepath):
if len(treepath) == 0:
return name_trial_gfs()
else:
return name_subgfs(treepath)
def type_subgfs(tree_path):
include_file('dune/pdelab/gridfunctionspace/subspace.hh', filetag='driver')
gfs = type_trial_gfs()
return "Dune::PDELab::GridFunctionSubSpace<{}, Dune::TypeTree::TreePath<{}> >".format(gfs, ', '.join(str(t) for t in tree_path))
@preamble(section="vtk", kernel="driver_block")
def define_subgfs(name, treepath):
t = type_subgfs(treepath)
gfs = name_trial_gfs()
return "{} {}({});".format(t, name, gfs)
def name_subgfs(treepath):
......@@ -474,14 +518,48 @@ def name_subgfs(treepath):
return name
@preamble(section="vtk", kernel="driver_block")
def define_subgfs(name, treepath):
t = type_subgfs(treepath)
gfs = name_trial_gfs()
return "{} {}({});".format(t, name, gfs)
def name_trial_subgfs(treepath):
if len(treepath) == 0:
return name_trial_gfs()
else:
return name_subgfs(treepath)
def type_subgfs(tree_path):
include_file('dune/pdelab/gridfunctionspace/subspace.hh', filetag='driver')
gfs = type_trial_gfs()
return "Dune::PDELab::GridFunctionSubSpace<{}, Dune::TypeTree::TreePath<{}> >".format(gfs, ', '.join(str(t) for t in tree_path))
def main_name_trial_subgfs(treepath):
if len(treepath) == 0:
return name_trial_gfs()
else:
# TODO
assert False
# return name_subgfs(treepath)
@class_member(classtag="driver_block")
def driver_block_get_gridfunctionsspace(element, is_dirichlet, root):
gfs_type = type_gfs(element, is_dirichlet, root=root)
name = name_gfs(element, is_dirichlet, root=root)
return ["std::shared_ptr<{}> get_gridfunctionsspace(){{".format(gfs_type),
" return {};".format(name),
"}"]
@preamble(section="postprocessing", kernel="main")
def main_typedef_trial_gfs(name, element, is_dirichlet):
driver_block_type = type_driver_block()
db_gfs_type = type_gfs(element, is_dirichlet)
return "using {} = {}::{};".format(name, driver_block_type, db_gfs_type)
def main_type_trial_gfs():
element = get_trial_element()
is_dirichlet = preprocess_leaf_data(element, "is_dirichlet")
gfs_type = type_gfs(element, is_dirichlet)
main_typedef_trial_gfs(gfs_type, element, is_dirichlet)
return gfs_type
@preamble(section="postprocessing", kernel="main")
def main_define_gfs(element, is_dirichlet, name, root):
driver_block_name = name_driver_block()
driver_block_get_gridfunctionsspace(element, is_dirichlet, root)
return "auto {} = {}.get_gridfunctionsspace();".format(name, driver_block_name)
......@@ -27,10 +27,10 @@ def interpolate_dirichlet_data(name):
@preamble(section="vector", kernel="driver_block")
def interpolate_vector(func, gfs, name):
return "Dune::PDELab::interpolate({}, {}, {});".format(func,
*gfs,
*name,
)
return "Dune::PDELab::interpolate({}, *{}, *{});".format(func,
gfs,
name,
)
@cached
......
......@@ -8,14 +8,13 @@ from dune.codegen.pdelab.driver import (get_form_ident,
preprocess_leaf_data,
)
from dune.codegen.pdelab.driver.gridfunctionspace import (name_leafview,
name_gfs,
name_trial_gfs,
main_name_trial_gfs,
type_leafview,
)
from dune.codegen.pdelab.driver.solve import name_vector
from dune.codegen.pdelab.driver.solve import main_name_vector
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def define_vtkfile(name):
ini = name_initree()
include_file("string", filetag="driver")
......@@ -27,7 +26,7 @@ def name_vtkfile():
return "vtkfile"
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def typedef_vtkwriter(name):
include_file("dune/grid/io/file/vtk/subsamplingvtkwriter.hh", filetag="driver")
gv = type_leafview()
......@@ -39,7 +38,7 @@ def type_vtkwriter():
return "VTKWriter"
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def define_subsamplinglevel(name):
ini = name_initree()
degree = get_trial_element().degree()
......@@ -55,7 +54,7 @@ def name_subsamplingintervals():
return "subint"
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def define_vtkwriter(name):
_type = type_vtkwriter()
gv = name_leafview()
......@@ -68,16 +67,16 @@ def name_vtkwriter():
return "vtkwriter"
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def vtkoutput():
include_file("dune/pdelab/gridfunctionspace/vtk.hh", filetag="driver")
vtkwriter = name_vtkwriter()
gfs = name_trial_gfs()
gfs = main_name_trial_gfs()
vtkfile = name_vtkfile()
predicate = name_predicate()
vec = name_vector(get_form_ident())
vec = main_name_vector(get_form_ident())
return ["Dune::PDELab::addSolutionToVTKWriter({}, {}, {}, Dune::PDELab::vtk::defaultNameScheme(), {});".format(vtkwriter, gfs, vec, predicate),
return ["Dune::PDELab::addSolutionToVTKWriter({}, *{}, *{}, Dune::PDELab::vtk::defaultNameScheme(), {});".format(vtkwriter, gfs, vec, predicate),
"{}.write({}, Dune::VTK::ascii);".format(vtkwriter, vtkfile)]
......@@ -86,7 +85,7 @@ def type_predicate():
return "CuttingPredicate"
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def define_predicate(name):
t = type_predicate()
return "{} {};".format(t, name)
......@@ -97,7 +96,7 @@ def name_predicate():
return "predicate"
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def typedef_vtk_sequence_writer(name):
include_file("dune/grid/io/file/vtk/vtksequencewriter.hh", filetag="driver")
gv_type = type_leafview()
......@@ -109,7 +108,7 @@ def type_vtk_sequence_writer():
return "VTKSW"
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def define_vtk_sequence_writer(name):
vtksw_type = type_vtk_sequence_writer()
vtkw_type = type_vtkwriter()
......@@ -123,15 +122,15 @@ def name_vtk_sequence_writer():
return "vtkSequenceWriter"
@preamble(section="vtk")
@preamble(section="vtk", kernel="main")
def visualize_initial_condition():
include_file("dune/pdelab/gridfunctionspace/vtk.hh", filetag="driver")
vtkwriter = name_vtk_sequence_writer()
element = get_trial_element()
gfs = name_trial_gfs()
vector = name_vector(get_form_ident())
gfs = main_name_trial_gfs()
vector = main_name_vector(get_form_ident())
predicate = name_predicate()
from dune.codegen.pdelab.driver.instationary import name_time
time = name_time()
return ["Dune::PDELab::addSolutionToVTKWriter({}, {}, {}, Dune::PDELab::vtk::defaultNameScheme(), {});".format(vtkwriter, gfs, vector, predicate),
return ["Dune::PDELab::addSolutionToVTKWriter({}, *{}, *{}, Dune::PDELab::vtk::defaultNameScheme(), {});".format(vtkwriter, gfs, vector, predicate),
"{}.write({}, Dune::VTK::appendedraw);".format(vtkwriter, time)]
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