FE Function from Numerical Data

General
1

Hi everyone,

To start, I’m using FEniCS 2019.1.0. So, suppose I’m working with a differential equation \partial_tu + Lu = f where u: [0,T] \times [0,1]^2 \to \mathbb{R}, the analytical solution of which (I’ll call this u_e) is given by an integral over time. Suppose also that I have created a mesh and function space V_h for the domain [0,1]^2, and that I can compute values of u_e at the nodal points of the grid using some numerical integration technique. If I want to compute the (say, L_2) error between my finite element solution and the analytical solution, how would I then interpolate/project u_e onto V_h?

2

Without having any specific formula for u_e, it is quite hard to give you a general rule that would work for your example. Please provide a minimal toy example that illustrates what u_e could be, i.e. what is u_e(x,y,t)?

3

Hi Dokken, in this case, u_e = \frac{1}{4\pi}\int_{0}^t \frac{1}{t - \tau} \exp\biggl( -\frac{\mathbf{r}^2}{4(t - \tau)} \biggr) d\tau .

4

I would do something along the lines of:

from dolfin import *
parameters["form_compiler"]["representation"] = "quadrature"
mesh = UnitSquareMesh(4, 4)

x = SpatialCoordinate(mesh)

# Create a 3rd order quadrature space
el = VectorElement("Quadrature", mesh.ufl_cell(), 3, quad_scheme="default")
V = FunctionSpace(mesh, el)
u = TrialFunction(V)
v = TestFunction(V)
a = inner(u, v)*dx(degree=3)
L = inner(x, v)*dx(degree=3)

x_c = Function(V)
solve(a == L, x_c)
x, y = x_c.sub(0, deepcopy=True).vector().get_local(), x_c.sub(1, deepcopy=True).vector().get_local()
import numpy as np
r = np.sqrt(x**2 + y**2)

import sympy
r_, t_, tau_ = sympy.symbols("r, t, tau")
ue_symbolic = sympy.integrate(1/4*sympy.pi*1/(t_-tau_)*sympy.exp(-r_*r_)/(4*(t_-tau_)), (tau_, 0, t_))
for t_i in [0, 0.1]:
    for r_i in r:
        print(ue_symbolic.subs({r_: r_i, t_: t_i}))

which would give you the value of u_ex at every point.

However, your integral does not converge, so you only get nan-s. See for instance: https://www.wolframalpha.com/input?i=integrate+1%2F(t-x)*exp(-1%2F(4*(t-x)))+%2C+x+from+0+to+t where i’ve chosen r=1