I could not ifind a way to make cube “holes” in a cubes mesh. I thought bcube may
do that but looking here, plugin/seq/msh3.cpp, it looks like there may be
flag bits for that but none seem to be coded in. Is there some simple way to do this?
Before getting to the Dirac Eqn, I thought I would play with this, lol,
This code sort of runs except it seems to add the BC to each cube and solve inside
each one. And the adaptive thing does not work but I can probably figure that out.
Thanks.
load “msh3”
load “medit”
real h=.25;
real s1=.25;
real s2=.25;
real cen=1;
real first=1.5;
real d=.6;
mesh3 Th3 = cube(20,20,20,[2x,2y,2*z]);
mesh3 Th4 = bcube(10,10,10, [s1*x+first,s1*y+cen,s1*z+h]);
mesh3 Th5 = bcube(10,10,10, [s2*x+first-d,s2*y+cen,s2*z+h]);
mesh3 Thx=Th3+Th4+Th5;
mesh3 Tsb=buildBdMesh(Th3);
//meshS surface=Tsb.Gamma;
int[int] labs = [1];
//meshS Thl = extract(Th3.Gamma, label=labs);
meshS Thl = extract(Th3, label=labs);
//plot(Thx,wait=true);
real tol=.5;
real cut=.01;
for(int i=0; i<1; ++i)
{
fespace Xh(Thx,P13d);
Xh psi, wx;
macro grad(u) [dx(u),dy(u),dz(u)]// def of grad operator
// Solve
solve potential(psi, wx)
= int3d(Thx)(
grad(psi)'*grad(wx)) // scalar product
//+on(Thl,psi=1)
+on(2,psi=2)
+on(1,psi=1)
+on(3,psi=1)
+on(4,psi=1)
+on(5,psi=1)
+on(6,psi=1)
;
plot(psi, value=true, wait=1);
medit(“asdf”,Thx);
//Thx = adaptmesh(Thx, [dx(psi), dy(psi),dz(psi)], splitpbedge=1, nbvx=1e5, abserror=1, cutoff=cut, err=tol, inquire=0, ratio=1.5, hmin=1./100);
tol=tol/2;
cut=cut/2;
} // i