Hello, I am new to Fenics and I was wondering if there is a way to evaluate the source and boundary condition functions exactly at the quadrature points rather than interpolating from the nodes as the ‘Expression’ functionality does. I thought that quadrature elements might be the solution, but I read that is just an optimization strategy whereby the values are interpolated and stored at the quadrature points.
For source terms, you can define them directly in UFL using SpatialCoordinate, e.g.,
from dolfin import *
N = 100
mesh = UnitIntervalMesh(N)
# Spatial coordinate accessible in UFL,
# without using Expression:
x = SpatialCoordinate(mesh)
# Source term using this coordinate:
f = sin(pi*x[0])
# Can be used directly in variational problems:
V = FunctionSpace(mesh,"CG",1)
u = TrialFunction(V)
v = TestFunction(V)
uh = Function(V)
solve(u*v*dx==f*v*dx,uh)
SpatialCoordinate can also be used similarly for weakly-enforced BCs built into the variational problem (e.g., Neumann BCs or Dirichlet BCs using Nitsche’s method). For strongly-enforced Dirichlet BCs (applied using a DirichletBC object), the data needs to be interpolated at the nodes of the space they are applied to anyway, so there is no “loss” associated with using an Expression interpolated into the same space.
Thank you for your feedback. However, I tried defining my exact solution in UFL, but the errornorm function complains. Instead, I used the norm_L2 = sqrt( assemble( ( u_exact - u_numerical )**2 * dx ) ), but I am concerned about floating point errors. Is there a way to avoid interpolating the exact solution without committing floating point errors? Thank you.
Anecdotally, I mostly use the “naive” error norm evaluation, and I’ve never run into floating point precision problems. If you’re seeing the convergence rates you expect when assembling the error directly, it’s probably not worth worrying about.