<jupyter_start><jupyter_text>P1 discontinuous SWIPDG, stationary linear elliptic ESV2007 problem==================================This example is about approximating the solution $u$ of the elliptic problem$$\begin{align}  -\nabla\cdot( \kappa \nabla u ) &= f   &&\text{in } \Omega\\                                u &= g_D &&\text{on }\partial\Omega\end{align}$$with datafunction as defined in `dune/gdt/test/linearelliptic/problems/ESV2007.hh` (see below) using a piecewise linear SWIPDG scheme, as in `dune/gdt/test/linearelliptic/discretizers/ipdg.hh`.Note that the discretization below contains handling of arbitrary Dirichlet and Neumann boundary data, although the problem at hand contains only trivial Dirichlet data.<jupyter_code>from dune.xt import common, grid, functions, la
from dune import gdt
common.init_mpi()<jupyter_output><empty_output><jupyter_text>$$\begin{align}  \Omega &= [-1, 1]^2\\  \Gamma_D &= \partial\Omega\\  \Gamma_N &= \emptyset\end{align}$$<jupyter_code>g = grid.make_cube_grid__2d_simplex_aluconform(lower_left=[-1, -1],
                                               upper_right=[1, 1],
                                               num_elements=[4, 4],
                                               num_refinements=2,
                                               overlap_size=[0, 0])
#g.visualize('../swipdg_esv2007_grid')
boundary_info = grid.make_boundary_info_on_leaf_layer(g, 'xt.grid.boundaryinfo.alldirichlet')
apply_on_neumann_boundary = grid.make_apply_on_neumann_intersections_leaf_part(boundary_info)
apply_on_dirichlet_boundary = grid.make_apply_on_dirichlet_intersections_leaf_part(boundary_info)<jupyter_output><empty_output><jupyter_text>$$\begin{align}\kappa(x) &:= 1\\f(x) &:= \tfrac{1}{2} \pi^2 \cos(\tfrac{1}{2} \pi x_0) \cos(\tfrac{1}{2} \pi x_1)\\g_D(x) &:= 0\end{align}$$Note that the grid `g` is only provided to select the correct _type_ of function. These functions do not rely on the actual grid which is why we need to later on provide the grid again, i.e., for `visualize(g)`.<jupyter_code>kappa = functions.make_constant_function_1x1(g, 1.0, name='diffusion')
f = functions.make_expression_function_1x1(g,
                                           'x',
                                           '0.5*pi*pi*cos(0.5*pi*x[0])*cos(0.5*pi*x[1])',
                                           order=3,
                                           name='force')
g_D = functions.make_constant_function_1x1(g, 0.0, name='dirichlet')
g_N = functions.make_constant_function_1x1(g, 0.0, name='neumann')
#kappa.visualize(g, '../swipdg_esv2007_diffusion')
#f.visualize(g, '../swipdg_esv2007_force')
space = gdt.make_dg_leaf_part_to_1x1_fem_p1_space(g)
#space.visualize('../swipdg_esv2007_dg_space')

# There are two ways to create containers: 
# * manually create them and given them to the operators/functionals
# * let those create appropriate ones

# in the CG example we chose the latter, so here we do the former
system_matrix = la.IstlRowMajorSparseMatrixDouble(space.size(), space.size(), space.compute_pattern('face_and_volume'))
swipdg_operator = gdt.make_elliptic_swipdg_matrix_operator(kappa, boundary_info, system_matrix, space)

rhs_vector = la.IstlDenseVectorDouble(space.size(), 0.0)
l2_force_functional = gdt.make_l2_volume_vector_functional(f, rhs_vector, space)

# there are two equivalent ways to restrict the integration domain of the face functional:
# * provide an apply_on_... tag on construction
# * provide an apply_on_... tag when appending the functional to the system assembler

l2_neumann_functional = gdt.make_l2_face_vector_functional(g_N, rhs_vector, space, apply_on_neumann_boundary)

swipdg_dirichlet_functional = gdt.make_elliptic_swipdg_dirichlet_vector_functional(g_D, kappa, boundary_info, rhs_vector, space)

# compute everything in one grid walk
system_assembler = gdt.make_system_assembler(space)
system_assembler.append(swipdg_operator)
system_assembler.append(l2_force_functional)
system_assembler.append(l2_neumann_functional)
system_assembler.append(swipdg_dirichlet_functional)
system_assembler.assemble()

# solve the linear system
u_h = la.IstlDenseVectorDouble(space.size(), 0.0)
la.make_solver(system_matrix).apply(rhs_vector, u_h)

# visualize (this will write swipdg_esv2007_solution.vtu)
gdt.make_discrete_function(space, u_h, 'solution').visualize('../swipdg_esv2007_solution')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.11534103784273737, 'fraction_numerical': 0.02107409925220938, 'mean_word_length': 4.010215664018161, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 13, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 17, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '3', 'dir': 'jupyter-structured-clean-dedup', 'id': '447356', 'license': 'BSD-2-Clause', 'n_tokens_mistral': 1529, 'n_tokens_neox': 1419, 'n_words': 300, 'path': 'notebooks/linearelliptic_swipdg.ipynb', 'repo_name': 'dune-community/dune-gdt-pymor-interaction', 'size': '4413'}