L2 projection from corse mesh FESpace to refined mesh FESpace

More
2 years 4 months ago #4438 by Wenzheng Kuang
Dear developers,

Assume Q0 is a lowest order L2 FESpace on the current coarse mesh. After hierarchical refinement of the mesh, I construct a new L2 FESpace Q on the refined mesh, without updating the "old" Q0 on the coarser mesh. I am wondering if NGSolve supports calculating the projection from L2 space Q0 on the coarser mesh to L2 space Q on the finer mesh, or vice versa? I try to do this by using BilinearForm(trialspace=Q0, testspace=Q), but this throws an error as expected. Is there any way to do this?
More
2 years 4 months ago #4439 by mneunteufel
Hi,

it is possible to interpolate a GridFunction to another one living on a different mesh.

gf1 = GridFunction(L2(mesh1,order=0))
gf1.Set(…)

gf2 = GridFunction(L2(mesh2,order=0))
gf2.Set(gf1)

This works in both directions. The Set method does an L2-projection.

Best,
Michael
The following user(s) said Thank You: Wenzheng Kuang
More
2 years 4 months ago #4441 by hvwahl
Hi,

you can also do the L2 projection of a FE function from one mesh to the other explicitly

from ngsolve import *
from ngsolve.meshes import MakeStructured2DMesh
from ngsolve.webgui import Draw

mesh1 = MakeStructured2DMesh(nx=10, ny=10)
Q0 = L2(mesh1, order=0)
gfu_q0 = GridFunction(Q0)
gfu_q0.Set(y)

Draw(gfu_q0, mesh1, 'q0')

mesh2 = MakeStructured2DMesh(nx=40, ny=40)
Q = L2(mesh2, order=0)

f = LinearForm(Q)
f += gfu_q0 * Q.TestFunction() * dx
f.Assemble()

gfu_q = GridFunction(Q)
gfu_q.vec.data = Q.InvM() * f.vec

Draw(gfu_q, mesh2, 'q')

Best wishes,
Henry
The following user(s) said Thank You: Wenzheng Kuang
More
2 years 4 months ago #4442 by joachim
Hi,

I guess you want to utilize the hierarchy of spaces, and not considering coarse and fine spaces as unrelated finite element spaces.

You can assemble matrices only on the current (fine) mesh. However, some spaces provide grid-transfer operators (=prolongation matrix). 

An examples is the (P1-iso-P2)-P1 element pairing for Stokes:
docu.ngsolve.org/latest/i-tutorials/unit....html#(P1-iso-P2)-P1

or multigrid-algorithms in the i-FEM tutorials:
github.com/JSchoeberl/iFEM/blob/master/m...rid/algorithms.ipynb

Should work for the embedding of the coarse-L2(order=0) into a finer one.

Joachim
The following user(s) said Thank You: Wenzheng Kuang
More
1 year 10 months ago #4599 by Wenzheng Kuang
Hi Dr. Schoberl,

Thanks for your previous reply.
However, when I follow the link  github.com/JSchoeberl/iFEM/blob/master/m...rid/algorithms.ipynb  and just change the H1 space to L2 (order=0) space, there is error in "fes.Prolongtaion().LevelDofs(level-1)" saying "Illegal level 0 num levels = 0". It seems that for now there is no LevelDofs() updating for the L2 space when the mesh is refined, right?
More
1 year 10 months ago - 1 year 10 months ago #4600 by str
Dear All,
I've been also curious about the L2 projection operator from the fine to the coarse levels - I guess , this is what is meant in the topic title). I've been testing an implementation of it a little (see other recent posts of mine) and found that MinRes, for example, throws an error when the prolongation matrix is computed between L2 spaces:
netgen.libngpy._meshing.NgException: Prolongation::GetNDofLevel not overloaded
This behaviour doesn't manifest for prolongations between H1 spaces. I attach the file where the only thing that happens is the multiplication of a course vector with something that involves prolongation matrix. In fact, it should be the coarse mass matrix computed from the fine level mass matrix. Any help is appreciated!

from netgen.csg import CSGeometry, OrthoBrick, Pnt
from ngsolve import *

# Geometry
cube = CSGeometry()
cube.Add(OrthoBrick(Pnt(-2.0, -2.0, -2.0), Pnt(2.0, 2.0, 2.0)))

# Mesh sizes
N=4
maxh = 4.0/N #to be refined later
maxH = 4.0/N

#Mesh generation
mesh   = Mesh(cube.GenerateMesh(maxh=maxh, quad_dominated=False))
mesh_c = Mesh(cube.GenerateMesh(maxh=maxH, quad_dominated=False))

#forcing
pi=3.141592
coef=CoefficientFunction(((3 * pi * pi + 1)*sin(pi * x)*sin(pi * y)*sin(pi * z)))

#FE Spaces
#changing all the L2 spaces to H1 resolves the bug
Vh = L2(mesh, order=1) #fine level space
u, v = Vh.TnT()

Rh = L2(mesh_c, order=1)#coarse level space

Rh_avatar= L2(mesh, order=1)#used only to compute prolongation operator from coarse to fine

#Mesh refinement from level=0 to level=1
mesh.ngmesh.Refine()
Vh.Update()
Rh_avatar.Update() #Rh remains to be coarse of level=0

# Bilinear forms:
mass = BilinearForm(Vh, symmetric=True) #fine level mass matrix M_h
mass += (u * v) * dx
mass.Assemble()

#Grid function
gfz=GridFunction(Rh)
gfz.Set(coef);

# Prolongation operator from mesh-level 0 to level 1:
prol = Rh_avatar.Prolongation().Operator(1)
coarse_to_fine_to_coarse_mass   = prol.T @ mass.mat @ prol#=coarse mass matrix

#computing coarse mass_matrix * coarse vector
diff = gfz.vec.CreateVector()
diff +=  coarse_to_fine_to_coarse_mass * gfz.vec
#!!!error: netgen.libngpy._meshing.NgException: Prolongation::GetNDofLevel not overloaded
exit(0);
Last edit: 1 year 10 months ago by str.
Time to create page: 0.120 seconds