NGSolve coding issue


I am trying to code a problem that involves two domains connected by an interface. There is pressure in the domains and at the interface. Therefore, I have a 2D and a 1D space. I tried to develop these pieces of code, but I am getting an anumerical factorization error( The matrix is singular). Any help would be greatly appreciated.

V   = VectorL2(mesh, order=order)
Qi = H1(mesh, order=order, orderinner=0, definedon=mesh.Boundaries("gammai"), dirichlet="gamma1|gamma2")
Q1 = L2(mesh,order=order,definedon=mesh.Materials("Domain1"))
Q2 = L2(mesh,order=order,definedon=mesh.Materials("Domain2"))
Qbar1 = FacetFESpace(mesh, order=order,definedon=mesh.Materials("Domain1"), dirichlet="gamma1")
Qbar2 = FacetFESpace(mesh, order=order, definedon=mesh.Materials("Domain2"), dirichlet="gamma2")

dx1 = dx(definedon=mesh.Materials("Domain1"))
ds1 = dx(element_boundary = True, definedon=mesh.Materials("Domain1"))

dx2 = dx(definedon=mesh.Materials("Domain2"))
ds2 = dx(element_boundary = True, definedon=mesh.Materials("Domain2"))

dsi = ds(skeleton = True, definedon=mesh.Boundaries("gammai"))
a += u*v*dx
a += p1*q1*dx1
a += p2*q2*dx2
a += pbar1*qbar1*ds1
a += pbar2*qbar2*ds2
a += p_i*qi*dsi
res = f.vec.CreateVector() = f.vec - a.mat * gfu.vec,inverse="umfpack")*res

Hi Dori,

if your mesh is not too fine, you can check the stiffness matrix by hand.

The piqidsi integral does not need the skeleton=True argument - it needs only the pressure at the interface.


Hi Joakim,
Thank you for your reply.I believe the error is on this line: Qi = H1(mesh, order=order, orderinner=0, definedon=mesh.Boundaries(“gammaf”), dirichlet=“gamma1|gamma2”).In the interface, we set the Dirichlet condition on the boundary (the two interface endpoints), as shown above. I am unsure if this is the correct approach. Any suggestions you may have would be greatly appreciated.
Thank you,

I am grateful for your assistance. However, I am experiencing another problem. Within my model, there is an equation in the interface that should be calculated in a linear form, but the code does not seem to be executing this correctly. No matter what coefficient function I apply to this right-hand side, the output remains consistent.
Any assistance with resolving this issue would be greatly appreciated.