Hello,
I am simulating a capacitive displacement sensor. To simulate the influence of a human, the object “detect” is introduced. A very high permiattivity is used for testing.
3 simulations are performed.
1: “detect” far away
2: “detect” close
3: “detect” close and permiattivity of “detect” set to the ambient value.
Results:
2.3658178708067787e-11
2.4048430781629297e-11
2.403768025497859e-11
The results for the first two simulations seem fine, the capacitance increases as expected. In the third simulation, 2.3637829146278974e-11 should come out, as the result should be identical to the value if “detect” was not present, but a much larger value comes out.
Does anyone have an idea where my error lies.
NGSolve-6.2.2302 on Win10 machine.
I am new to the forum, so I am not allowed to upload anything. Therefore the source code follows in the text.
from netgen.occ import *
from ngsolve import *
from collections import defaultdict
from ngsolve.internal import visoptions, viewoptions, Rotate, Zoom
import netgen.gui
ngsglobals.msg_level = 0
def run(value1,value2):
outside= (Cylinder(Pnt(0,0,00), Z, r=20, h=0.5) - Cylinder(Pnt(0,0,00), Z, r=8, h=0.5)).bc("plus").mat("egal")
inside= Cylinder(Pnt(0,0,00), Z, r=7.9, h=0.5).bc("minus").mat("egal")
detect = Box(Pnt(-5,-5,value1),Pnt(5,5,35)).mat("human")
air = Sphere( Pnt(0,0,0), 40 ) - outside - inside - detect - Box(Pnt(-40,-40,-40),Pnt(40,40,0))
air= air.mat("air")
geo = OCCGeometry(Glue([air, outside, inside, detect]))
mesh = Mesh(geo.GenerateMesh (maxh=5))
mesh.Curve(3)
diel_perm =defaultdict(lambda: -1)
diel_perm["air"]=1.00059
diel_perm["human"]=value2
diel_perm_cf = CoefficientFunction([ diel_perm[mat] for mat in mesh.GetMaterials() ] )
voltage = 4
potential =defaultdict(lambda: 0)
potential["plus"]=voltage/2
potential["minus"]=-voltage/2
potential_cf = [ potential[bnd] for bnd in mesh.GetBoundaries() ]
potential_str="plus|minus"
dirichlet_str="air|plus|minus|human"
eps0=8.8541878128E-12
scale=0.01 #=>geometry values in cm
###################
fes = H1(mesh, order=3, dirichlet=dirichlet_str)
(ut,vt) = fes.TnT()
a = BilinearForm(fes)
a += diel_perm_cf*grad(ut)*grad(vt)*dx
a.Assemble()
u = GridFunction(fes, name="Potential")
u.Set(potential_cf, definedon=mesh.Boundaries(potential_str))
f = u.vec.CreateVector()
f.data = a.mat * u.vec
u.vec.data -= a.mat.Inverse(fes.FreeDofs(), inverse="sparsecholesky") * f
E = -Grad(u)/scale
#cap = eps0*Integrate(E*E*diel_perm_cf, mesh, VOL, definedon=mesh.Boundaries("plus")) / voltage**2 *scale*scale*scale
cap = eps0*Integrate(E*E*diel_perm_cf, mesh, VOL) / voltage**2 *scale*scale*scale
####################
Draw(u, mesh,"U")
Rotate(0,180)
Zoom(50)
viewoptions.clipping.enable = 1
visoptions.clipsolution = "scal"
visoptions.usetexture = 1
visoptions.lineartexture = 1
viewoptions.drawnetgenlogo=0
#netgen.Redraw(blocking=True, fr=1e8)
netgen.libngguipy.Snapshot(10,10)
netgen.gui.Snapshot(w=800,h=800, filename=str(value1)+"_"+str(value2)+".jpg")
###################
print(value1,value2,cap)
run(34.9,8000)
run(15,8000)
run(15,1.00059)
Many greetings
Rüdiger