Regarding computation of the error norm

General
1

Hello everyone,
I wish to compute the error norm, based on the solution of Stokes equation and the steady state solution of the Navier-Stokes (NS) equation.

Following is my code where u is the velocity field obtained by solving the NS equations.

from dolfin import *
from mshr import *
import matplotlib.pyplot as plt

mesh = Mesh('../../res/mesh_generation/annulus_refined.xml')

V = VectorFunctionSpace(mesh, "CG", 4)
W = FunctionSpace(mesh, "CG", 4)

def curl2D(w):
    return as_vector([w.dx(1),-w.dx(0)])

def laplacian(w):
    return as_vector([w[0].dx(0).dx(0)+w[0].dx(1).dx(1),w[1].dx(0).dx(0)+w[1].dx(1).dx(1)])

mu = 1.0

wi = Expression("sin(x[1])*cos(x[0])",degree=6)
ti = Expression("x[1]",degree=6)

w = project(wi,W)
t = project(ti,W)
u_bar = curl2D(w)
f     = -mu*laplacian(u_bar) + grad(t)

errornorm(u_bar, u, 'l2', 4)

The error I receive

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
/home/art/thesis/source/Re_1_ST/plotData.ipynb Cell 5 in <cell line: 35>()
     32 u_bar = curl2D(w)
     33 f     = -mu*laplacian(u_bar) + grad(t)
---> 35 errornorm(u_bar, u, 'l2', 4)

File /usr/lib/petsc/lib/python3/dist-packages/dolfin/fem/norms.py:246, in errornorm(u, uh, norm_type, degree_rise, mesh)
    243 degree = uh.ufl_element().degree() + degree_rise
    245 # Check degree of 'exact' solution u
--> 246 degree_u = u.ufl_element().degree()
    247 if degree_u is not None and degree_u < degree:
    248     cpp.warning("Degree of exact solution may be inadequate for accurate result in errornorm.")

AttributeError: 'ListTensor' object has no attribute 'ufl_element'

Following is my mesh :

annulus
annulus1920Ă—920 77 KB

If I am not wrong, the function errornorm works if the arguements to it are of type dolfin objects. However u_bar is an ufl-element (tensor). Therefore is it possible to convert it to dolfin object?

Any leads would highly be appreciated.

2

I don’t have an install of old dolfin on any of my machines but perhaps you can try computing the error manually since errornorm expects an Expression for the analytical solution. Something like:

error_L2 = assemble((u - ubar)**2 * dx)**0.5
print(error_L2)

If I recall correctly old dolfin does this assembly collectively over all processes so should give you the total over all processes.

3

@nate , I am thankful to you for the prompt response. I am using the fenics 2019.1.0

When I compute it manually, I get the following error:

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
/home/art/thesis/source/Re_1_ST/plotData.ipynb Cell 5 in <cell line: 35>()
     32 u_bar = curl2D(w)
     33 f     = -mu*laplacian(u_bar) + grad(t)
---> 35 error_l2 = assemble(((u_bar-u)**2 )* dx)**0.5

File /usr/lib/python3/dist-packages/ufl/measure.py:428, in Measure.__rmul__(self, integrand)
    426 domain = self.ufl_domain()
    427 if domain is None:
--> 428     domains = extract_domains(integrand)
    429     if len(domains) == 1:
    430         domain, = domains

File /usr/lib/python3/dist-packages/ufl/domain.py:346, in extract_domains(expr)
    344 for t in traverse_unique_terminals(expr):
    345     domainlist.extend(t.ufl_domains())
--> 346 return sorted(join_domains(domainlist))

TypeError: '<' not supported between instances of 'Mesh' and 'Mesh'

Furthermore, when I use the expanded expression suggested on the page for errornorm , :

norm_l2 = inner(u,u) + inner(u_bar,u_bar) - 2*inner(u,u_bar)
error_norm = assemble(norm_l2*dx)**0.5

I get the same error. What exactly is the error trying to suggest?

4

It doesn’t look like your code shows us what u is. The error indicates that it lives on a different mesh than u_bar, so you’ll have to come up with a clever way to get them on the same mesh, and carefully consider what that will mean in terms of the error.

1 Like
5

@nate , you were correct the problem was with the mesh.