Line integral and Refine the mesh

General
1

Hello,

I have two questions:
First, I am trying to understand why the following code is not calculating the line integral correctly: I get the same values for all three integrals. Is it because of the integrate function or because of the ds1, ds2 definition:

from ngsolve import *
from netgen.geom2d import SplineGeometry
from math import pi
ngsglobals.msg_level = 1
from sympy import Symbol
import sys

N = 8
hmax = 0.4
order = 1

def topBottom():
    geometry = SplineGeometry()
    pnts = [(0, 0), (0.5, 0), (1, 0), (1, 1), (0.5, 1), (0, 1)]
    pnums = [geometry.AppendPoint(*p) for p in pnts]
    lines = [(pnums[0], pnums[1], "gamma1", 1, 0),
             (pnums[1], pnums[2], "gamma2", 2, 0),
             (pnums[2], pnums[3], "gamma2", 2, 0),
             (pnums[3], pnums[4], "gamma2", 2, 0),
             (pnums[4], pnums[5], "gamma1", 1, 0),
             (pnums[5], pnums[0], "gamma1", 1, 0),
             (pnums[1], pnums[4], "gammaf", 1, 2)]
    for p1, p2, bc, left, right in lines:
        geometry.Append(["line", p1, p2], bc=bc, leftdomain=left, rightdomain=right)
    geometry.SetMaterial(1, "Domain1")
    geometry.SetMaterial(2, "Domain2")
    return geometry

mesh = Mesh(topBottom().GenerateMesh(maxh=hmax))

n = specialcf.normal(mesh.dim)
h = specialcf.mesh_size

Draw(mesh)
Qbar = H1(mesh, order=order, orderinner=0, dirichlet="gamma1|gamma2")
fes = H1(mesh, order=1)

myfunc = CoefficientFunction(100 * (4-x+y)**3)

ds1 = ds(element_boundary=True, definedon=mesh.Materials("Domain1"))
ds2 = ds(element_boundary=True, definedon=mesh.Materials("Domain2"))
ds12 = ds(element_boundary=True, definedon=mesh.Materials("Domain1|Domain2"))

val1 = Integrate(myfunc * ds1, mesh)
val2 = Integrate(myfunc * ds2, mesh)
val12 = Integrate(myfunc * ds12, mesh)

print(val1, val2, val12)

My second question concerns a previous post to which I still didn’t get any reply (Zeromatrix) post.

Let

Q1 = H1(mesh, order=order,dirichlet="gamma1|gamma2")
Q = Compress(Q1, active_dofs=Q1.GetDofs(mesh.Boundaries("gammaf")))

To refine the mesh I am the following solve function:

def Solve():
        Q1.Update()
        Q.Update()
        gfu.Update()
        a.Assemble()
        f.Assemble()
....
for i in range(4):
      mesh.Refine()
     Solve()

However, it’s not functioning as intended.

I would greatly appreciate any assistance you can provide.