Hi all. I want to know that why L2 error increases with increase in mesh size if I am using Lagrange P2 element . Whether I am doing something wrong. The code which is given below I am using this one and get 2.952056 error with mesh size 1/2 and 7.1425555 with mesh size 1/16. Can someone please tell me why I am getting this much error.
from __future__ import print_function
from fenics import *
import matplotlib.pyplot as plt
# Create mesh and define function space
mesh = UnitSquareMesh(2,2)
V = FunctionSpace(mesh, 'P', 2)
# Define boundary condition
u_D = Expression('1 + x[0]*x[0] + 2*x[1]*x[1]', degree=2)
def boundary(x, on_boundary):
return on_boundary
bc = DirichletBC(V, u_D, boundary)
# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Constant(-6.0)
a = dot(grad(u), grad(v))*dx
L = f*v*dx
# Compute solution
u = Function(V)
solve(a == L, u, bc)
# Plot solution and mesh
plot(u)
plot(mesh)
# Save solution to file in VTK format
vtkfile = File('poisson/solution.pvd')
vtkfile << u
# Compute error in L2 norm
error_L2 = errornorm(u_D, u, 'L2')
# Compute maximum error at vertices
vertex_values_u_D = u_D.compute_vertex_values(mesh)
vertex_values_u = u.compute_vertex_values(mesh)
import numpy as np
error_max = np.max(np.abs(vertex_values_u_D - vertex_values_u))
# Print errors
print('error_L2 =', error_L2)
print('error_max =', error_max)
# Hold plot
plt.show()
I can not reproduce this behaviour. When executing your code, I get values for error_L2 somewhere around 1e-15, which is practically zero, no matter the mesh size. This is expected with your chosen exact solution, as it is contained in your finite element function space.
and I did the following and run and got error_L2 = 0.14162139723192527
with mesh size 1/2 and error_L2 = 0.48687396343176065 with mesh size 1/8. Thats what I was saying L2 error increases with finer and finer mesh. Kindly help me regarding the same. Thanku.
from __future__ import print_function
from fenics import *
import matplotlib.pyplot as plt
# Create mesh and define function space
mesh = UnitSquareMesh(4,4)
V = FunctionSpace(mesh, 'Lagrange', 2)
# Define boundary condition
u_D = Expression('sin(pi*x[0])*sin(pi*x[1])*sin(pi*x[0])*sin(pi*x[1]) + cos(pi*x[0])*cos(pi*x[1])', degree=2)
def boundary(x, on_boundary):
return on_boundary
bc = DirichletBC(V, u_D, boundary)
# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Expression('-8*pi*pi*pi*pi*sin(pi*x[1])*sin(pi*x[1])*cos(2*pi*x[0]) - 8*pi*pi*pi*pi*sin(pi*x[0])*sin(pi*x[0])*cos(2*pi*x[1]) + 8*pi*pi*pi*pi*cos(2*pi*x[0])*cos(2*pi*x[1]) \
+ 4 *pi*pi*pi*pi*cos(pi*x[0])*cos(pi*x[1])', degree=0)
n = FacetNormal(mesh)
h = 2.0*Circumradius(mesh)
h_avg = (h('+') + h('-'))/2
# Parameters
alpha = 6
# Bilinear form
a = dot(div(grad(u)), div(grad(v)))*dx \
- dot(avg(div(grad(u))), jump(grad(v), n))*dS \
- dot(jump(grad(u), n), avg(div(grad(v))))*dS \
+ alpha/h_avg*dot(jump(grad(u), n), jump(grad(v),n))*dS
#print('a=',a)
L = f*v*dx
# Compute solution
u = Function(V)
solve(a == L, u, bc)
# Plot solution and mesh
plot(u)
plot(mesh)
# Save solution to file in VTK format
vtkfile = File('biharmonic/solution.pvd')
vtkfile << u
# Compute error in L2 norm
error_L2 = errornorm(u_D, u, 'L2')
# Compute maximum error at vertices
vertex_values_u_D = u_D.compute_vertex_values(mesh)
vertex_values_u = u.compute_vertex_values(mesh)
import numpy as np
error_max = np.max(np.abs(vertex_values_u_D - vertex_values_u))
# Print errors
print('error_L2 =', error_L2)
print('error_max =', error_max)
# Hold plot
plt.show()
I again did for poisson problem but got the same . Kindly tell me why this is happend . Error must be reduced by doing finner mesh but I got opposite of it. Kindly help me. Thanku.
ravina@ravina:~/Downloads$ python3 poisson2.py
Solving linear variational problem.
*** Warning: Degree of exact solution may be inadequate for accurate result in errornorm.
error_L2 = 2.9520565726363387e-15(with mesh size(1/2))
error_max = 6.661338147750939e-16
ravina@ravina:~/Downloads$ python3 poisson2.py
Solving linear variational problem.
*** Warning: Degree of exact solution may be inadequate for accurate result in errornorm.
error_L2 = 5.603873381595536e-15(with mesh size 1/4)
error_max = 5.329070518200751e-15
from __future__ import print_function
from fenics import *
import matplotlib.pyplot as plt
# Create mesh and define function space
mesh = UnitSquareMesh(4,4)
V = FunctionSpace(mesh, 'P', 2)
# Define boundary condition
u_D = Expression('1 + x[0]*x[0] + 2*x[1]*x[1]', degree=2)
def boundary(x, on_boundary):
return on_boundary
bc = DirichletBC(V, u_D, boundary)
# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Constant(-6.0)
a = dot(grad(u), grad(v))*dx
L = f*v*dx
# Compute solution
u = Function(V)
solve(a == L, u, bc)
# Plot solution and mesh
plot(u)
plot(mesh)
# Save solution to file in VTK format
vtkfile = File('poisson2/solution.pvd')
vtkfile << u
# Compute error in L2 norm
error_L2 = errornorm(u_D, u, 'L2')
# Compute maximum error at vertices
vertex_values_u_D = u_D.compute_vertex_values(mesh)
vertex_values_u = u.compute_vertex_values(mesh)
import numpy as np
error_max = np.max(np.abs(vertex_values_u_D - vertex_values_u))
# Print errors
print('error_L2 =', error_L2)
print('error_max =', error_max)
# Hold plot
plt.show()
In poisson It is not showing clearly that error increases in some sense. but we can easily see in case of biharmonic.
ravina@ravina:~/Downloads$ python3 biharmonic.py
Solving linear variational problem.
*** Warning: Degree of exact solution may be inadequate for accurate result in errornorm.
error_L2 = 0.14162139723192527 (with mesh size 1/2)
error_max = 0.3623117469334274
ravina@ravina:~/Downloads$ python3 biharmonic.py
Solving linear variational problem.
*** Warning: Degree of exact solution may be inadequate for accurate result in errornorm.
error_L2 = 0.370765075499436 (with mesh size 1/4)
error_max = 0.6318280069181568
Here is me code file:
from __future__ import print_function
from fenics import *
import matplotlib.pyplot as plt
# Create mesh and define function space
mesh = UnitSquareMesh(4,4)
V = FunctionSpace(mesh, 'Lagrange', 2)
# Define boundary condition
u_D = Expression('sin(pi*x[0])*sin(pi*x[1])*sin(pi*x[0])*sin(pi*x[1]) + cos(pi*x[0])*cos(pi*x[1])', degree=2)
def boundary(x, on_boundary):
return on_boundary
bc = DirichletBC(V, u_D, boundary)
# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Expression('-8*pi*pi*pi*pi*sin(pi*x[1])*sin(pi*x[1])*cos(2*pi*x[0]) - 8*pi*pi*pi*pi*sin(pi*x[0])*sin(pi*x[0])*cos(2*pi*x[1]) + 8*pi*pi*pi*pi*cos(2*pi*x[0])*cos(2*pi*x[1]) \
+ 4 *pi*pi*pi*pi*cos(pi*x[0])*cos(pi*x[1])', degree=0)
n = FacetNormal(mesh)
h = 2.0*Circumradius(mesh)
h_avg = (h('+') + h('-'))/2
# Parameters
alpha = 6
# Bilinear form
a = dot(div(grad(u)), div(grad(v)))*dx \
- dot(avg(div(grad(u))), jump(grad(v), n))*dS \
- dot(jump(grad(u), n), avg(div(grad(v))))*dS \
+ alpha/h_avg*dot(jump(grad(u), n), jump(grad(v),n))*dS
#print('a=',a)
L = f*v*dx
# Compute solution
u = Function(V)
solve(a == L, u, bc)
# Plot solution and mesh
plot(u)
plot(mesh)
# Save solution to file in VTK format
vtkfile = File('biharmonic/solution.pvd')
vtkfile << u
# Compute error in L2 norm
error_L2 = errornorm(u_D, u, 'L2')
# Compute maximum error at vertices
vertex_values_u_D = u_D.compute_vertex_values(mesh)
vertex_values_u = u.compute_vertex_values(mesh)
import numpy as np
error_max = np.max(np.abs(vertex_values_u_D - vertex_values_u))
# Print errors
print('error_L2 =', error_L2)
print('error_max =', error_max)
# Hold plot
