- Thank you received: 0
Algebraic Multigrid for Compressed VectorH1 space
3 years 8 months ago - 3 years 8 months ago #3566
by Ryanleaf
Replied by Ryanleaf on topic Algebraic Multigrid for Compressed VectorH1 space
Dear Christoph,
Thank you very much for fixing the bug. This update fixed the boundary issue.
I am sorry I think may meet another problem about multigrid.
Same script:
I noticed when I refine maxh to 1/8, I begin to get the following error:
Would you please let me know where to find a description about which version of algebraic multigrid is currently implemented?
Again thank you very much for your kind help.
Best,
Qi
Thank you very much for fixing the bug. This update fixed the boundary issue.
I am sorry I think may meet another problem about multigrid.
Same script:
Code:
from netgen.geom2d import SplineGeometry
from netgen.meshing import MeshingParameters
from xfem import *
from netgen.csg import CSGeometry, OrthoBrick, Pnt
import netgen.gui
from netgen.geom2d import unit_square
from ngsolve import *
from ngsolve.la import EigenValues_Preconditioner
from xfem.lsetcurv import *
def eps(u):
return 0.5 * (Grad(u) + Grad(u).trans)
def Setup():
cube = CSGeometry()
cube.Add(OrthoBrick(Pnt(-1.5, -1.5, -1.5), Pnt(1.5, 1.5, 1.5)))
mesh = Mesh(cube.GenerateMesh(maxh=1/8, quad_dominated=False))
lsetmeshadap = LevelSetMeshAdaptation(mesh, order=3, threshold=1000, discontinuous_qn=True)
phi = Norm(CoefficientFunction((x, y, z))) - 1
deformation = lsetmeshadap.CalcDeformation(phi)
lsetp1 = lsetmeshadap.lset_p1
mesh.SetDeformation(deformation)
ci = CutInfo(mesh, lsetp1)
Vhbase = VectorH1(mesh, order=2, dirichlet=[1,2,3,4,5,6])
Vhneg = Restrict(Vhbase, ci.GetElementsOfType(HASNEG))
Vhpos = Restrict(Vhbase, ci.GetElementsOfType(HASPOS))
WhGV = FESpace([Vhneg, Vhpos], flags={"dgjumps": True})
lset_neg = {"levelset": lsetp1, "domain_type": NEG}
lset_pos = {"levelset": lsetp1, "domain_type": POS}
lset_doms = [lset_neg, lset_pos]
elements_doms = [HASNEG, HASPOS]
u, v = WhGV.TnT()
a = BilinearForm(WhGV)
for i in [0, 1]: # loop over domains
a += SymbolicBFI(lset_doms[i], form=2* InnerProduct(eps(u[i]), eps(v[i])),
definedonelements=ci.GetElementsOfType(elements_doms[i]))
a += SymbolicBFI(lset_doms[i], form=u[i] * v[i], definedonelements=ci.GetElementsOfType(elements_doms[i]))
return mesh, WhGV, a
mesh, WhGV, a = Setup()
gfu = GridFunction(WhGV)
gfu.components[1].Set(CoefficientFunction((x,y,z)),BND)
mg = Preconditioner(a, "h1amg") # register mg to a
a.Assemble()
print('NDOF = ', WhGV.ndof)
preI = Projector(mask=WhGV.FreeDofs(), range=True)
lams = EigenValues_Preconditioner(mat=a.mat, pre=mg.mat)
print("condition number", max(lams)/min(lams))
Am I having inappropriate settings in this script to causing this problem?Traceback (most recent call last):
File "/tmp/iterative_solver_krylov_hybrid.py/main.py", line 42, in <module>
a.Assemble()
netgen.libngpy._meshing.NgException: Inverse matrix: Matrix singularin Assemble BilinearForm 'biform_from_py'
Would you please let me know where to find a description about which version of algebraic multigrid is currently implemented?
Again thank you very much for your kind help.
Best,
Qi
Time to create page: 0.115 seconds