i am trying to calculate error and convergence in terms of mesh. how do i do it, any help, how to refine the mesh is the problem in the for loop. This the mesh
mesh = RectangleMesh(Point(0.0, 0.0),Point(2.0math.pi, 2.0math.pi),2**(6), 2**(6))
You can use mesh=refine(mesh) inside the loop
Does that include changing the dimension of x and y or i have to define that as well
For a convergence study, it is common to have a mesh of the same area, and increasing the number of elements describing this area. Refining the mesh will keep the same outer boundaries, But increase the «resolution» of your discretization.
ok, i do understand but it gives error and as well, with the mesh define above, i would like to know how to just change the values of the mesh without using the mesh=refine (mesh)
What do you mean by changing the values of the mesh? A mesh connsists of points, and cells (connectivity between points). You can either move the points, (without changing connectivity), add more points and change connectivity with either doing refine or remeshing
i want to know how to do mesh convergence, i am trying to write a code in spatial discretization.
the mesh i used was rectangle, which is
mesh = RectangleMesh(Point(0.0, 0.0),Point(2.0math.pi, 2.0math.pi),2**(6), 2**(6)),
want to know how to change the point ((2**(6), 2**(6)) to say [4, 8, 16, 32, 64, 128, 264] in the for loop
That is trivial:
for i in range(2,8):
mesh = RectangleMesh(Point(0.0, 0.0),Point(2.0math.pi, 2.0math.pi),2**(i), 2**(i)),
how to i then compute the error, since the mesh are not equal. how do i interpolate the mesh to be able to compute the error and rate of convergence.
It Depends on what error you are measuring. If you are measuring against an exact solution, you can do
e_i = assemble((uh-u_ex)**2*dx(domain=mesh_i))
Where the convergence rate is
r=log(e_i/e_{i-1})/log(h_i/h_{i-1}) Where h_i is the largest celldiameter in mesh_i.
This is the normal way of computing convergence rates.
To interpolate data over non-matching meshes, see:
Maybe you have a function named ‘ch’ elsewhere in your code? What happens if you try different name? And for interpolating last solution obtained, try looking at the previous post on projecting data from one mesh to other. Once projected, do as you would do with exact solution.
1 Like
how do you interpolate unto same space
As i said earlier, use the LagrangeInterpolator:
ok and how do i compute the error and convergence if i an not measuring with exact solution, as well the domain=mesh_i gives error, saying mesh_i not define
Given a function f_i on mesh_i, you would have to interpolate f_i on to mesh_{i+1} using the Langrange interpolator (Lets call the interpolated result U_i. Then you can measure error
e_{i+1}= numpy.sqrt(assemble(inner(u_{i+1}-U_i, u_{i+1}-U_i)*dx))
The rate of convergence would then be:
r_{i+1} = numpy.log(e_{i+1}/e_{i})/numpy.log(h_{i+1}/h_{i})
where h_i denotes the mesh size for the i-th mesh.
As when it comes to errors, we can’t really help you as you have not supplied a minimal example illustrating the error. And when I say minimal, I mean minimal, not your full script, but a self contained code that illustrates the error in as few lines as possible.
1 Like
As I said, you should NOT supply the full script, but a minimal example. Additionally, format code by enclosing it in ``` on each side
ok, will do that. so do i give an example of the error analysis and convergence. because i do not see how to input the Langrange interpolator
the mesh_k gives error and as well says not define from using the mesh define earlier
From the link that I shared with you several times now, which you clearly havent looked at, you can do the following:
class Quadratic2D(UserExpression):
def eval(self, values, x):
values[0] = x[0]*x[0] + x[1]*x[1] + 1.0
def test_functional2D():
"""Test integration of function interpolated in non-matching meshes"""
f = Quadratic2D(degree=2)
# Interpolate quadratic function on course mesh
mesh0 = UnitSquareMesh(8, 8)
V0 = FunctionSpace(mesh0, "Lagrange", 2)
u0 = Function(V0)
LagrangeInterpolator.interpolate(u0, f)
# Interpolate FE function on finer mesh
mesh1 = UnitSquareMesh(31, 31)
V1 = FunctionSpace(mesh1, "Lagrange", 2)
u1 = Function(V1)
LagrangeInterpolator.interpolate(u1, u0)
assert round(assemble(u0*dx) - assemble(u1*dx), 10) == 0
mesh1 = UnitSquareMesh(30, 30)
V1 = FunctionSpace(mesh1, "Lagrange", 2)
u1 = Function(V1)
LagrangeInterpolator.interpolate(u1, u0)
assert round(assemble(u0*dx) - assemble(u1*dx), 10) == 0
The first operation interpolates f onto u0.
The next interpolation operator interpolates u0 onto the finer space u1.
Therefore, if you have another function u2=Function(V1) on the finer mesh, you can now do
assemble(inner(u2-u1, u2-u1)*dx).
1 Like
I am going to close this topic as you are clearly not reading my replies, only hoping that i will write out the specific code for you. I have given several examples now on how to use the LagrangeInterpolator.