# Region and label

**URL:** <https://community.freefem.org/t/region-and-label/2217>\
**Category:** General Discussion\
**Created:** [January 13, 2023, 3:18pm UTC](https://community.freefem.org/t/region-and-label/2217 "2023-01-13T15:18:29Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![marchywka](https://avatars.discourse-cdn.com/v4/letter/m/ee59a6/32.png) [@marchywka](https://community.freefem.org/u/marchywka)\
**Post date:** [January 13, 2023, 3:18pm UTC](https://community.freefem.org/t/region-and-label/2217/1 "2023-01-13T15:18:29Z")

</div>

In movemesh23, the code I copied had a “region” parameter but I also added  
a “label” parameter that appears to work with the boundary conditions. Do these  
conflict or are they different? The resulting code seems to work as I could include  
neumann conditions on the X planes and verify in R the solution was mostly planar  
although the slope at the one end was closer to .1 than .2 . Did I setup the  
weak formulation Neumann conditions right? It looks like you get zero by default.

Also, when dumping the psi.txt, I used code that I copied from somewhere else  
and it looks a bit odd and is probably dumping the same point multiple times.  
Is the mesh indexing discussed somewhere?

This code is supposed to be two cubes ( cars ) in an airstream ( copied the code  
from the airfoil example in the manual ). I was debating about trying a full Navier  
Stokes with this example or moving on the Dirac or Landau-Lifschitz- Gilbert. 🙂  
Has anyone encountered issues with weak formulations that may be in the examples?  
Thanks, so far so good.

///home/ubuntu/dev/freefem/FreeFem-sources-master/examples/3d$ vi Poisson-cube-ballon.edp  
//marchywka@happy:/home/ubuntu/dev/freefem/FreeFem-sources-master/examples/3d$

load “msh3”  
load “tetgen”  
load “medit”

//  
real volumetet; // use in tetg.  
meshS Thx0,Thx1, Thy0, Thy1, Thz0, Thz1;  
int nx=15,ny=15,nz=15;  
func meshS mycube ( real x0, real x1, real y0,real y1, real z0, real z1, int orient, int lb)  
{  
meshS ThHex;  
// first build the 6 faces of the hex.  
//real x0=-1,x1=1;  
//real y0=-1.1,y1=1.1;  
//real z0=-1.2,z1=1.2;

// a volume of on tet.  
volumetet= (x1-x0)_(y1-y0)_(z1-z0)/ (nx_ny_ny) /6.;

mesh Thx = square(ny,nz,[y0+(y1-y0)\*x,z0+(z1-z0)\*y]);  
mesh Thy = square(nx,nz,[x0+(x1-x0)\*x,z0+(z1-z0)\*y]);  
mesh Thz = square(nx,ny,[x0+(x1-x0)\*x,y0+(y1-y0)\*y]);

int[int] refz=[0,5]; // bas  
int[int] refZ=[0,6]; // haut  
int[int] refy=[0,3]; // devant  
int[int] refY=[0,4]; // derriere  
int[int] refx=[0,1]; // gauche  
int[int] refX=[0,2]; // droite  
int f1=-1;  
int f2=1;  
if (orient==1) { f1=1; f2=-1; }  
int[int] l0=[0,lb+0];  
//meshS  
Thx0 = movemesh23(Thx,transfo=[x0,x,y],label=l0,orientation=f1,region=refx,removeduplicate=false);  
++l0[1];  
//meshS  
Thx1 = movemesh23(Thx,transfo=[x1,x,y],label=l0,orientation=f2,region=refX,removeduplicate=false);  
++l0[1];  
//meshS  
Thy0 = movemesh23(Thy,transfo=[x,y0,y],label=l0,orientation=f2,region=refy,removeduplicate=false);  
++l0[1];  
//meshS  
Thy1 = movemesh23(Thy,transfo=[x,y1,y],label=l0,orientation=f1,region=refY,removeduplicate=false);  
++l0[1];  
//meshS  
Thz0 = movemesh23(Thz,transfo=[x,y,z0],label=l0,orientation=f1,region=refz,removeduplicate=false);  
++l0[1];  
//meshS  
Thz1 = movemesh23(Thz,transfo=[x,y,z1],label=l0,orientation=f2,region=refZ,removeduplicate=false);

cout\<\<“l0=”\<\<l0\<\<endl;

//medit(" — ", Thx0,Thx1,Thy0,Thy1,Thz0,Thz1);  
ThHex = Thx0+Thx1+Thy0+Thy1+Thz0+Thz1;  
return ThHex;  
}

real x0=-1,x1=1; real yy0=-1.1,yy1=.1; real z0=-1.2,z1=.2;  
real xa0=-.5; real ya0=-.8; real za0=-1;  
real xb0=-.1; real yb0=-.8; real zb0=-1;  
real sz=.25;

meshS ThHex = mycube(x0,x1,yy0,yy1,z0,z1,1,0);  
real vol=volumetet;  
meshS Thsrc=Thx0;  
nx=10; ny=10; nz=10;  
meshS Tha = mycube(xa0,xa0+sz,ya0,ya0+sz,za0,za0+sz,-1,10);  
meshS Thb = mycube(xb0,xb0+sz,yb0,yb0+sz,zb0,zb0+sz,-1,20);  
meshS Thc= Tha+Thb+ThHex;  
//medit(" asdf ",Thc);

real[int] domaine = [-.9,-1,-1.1,.5,vol];  
mesh3 Th = tetg(Thc,switch=“pqaAAYYQ”,nbofregions=1,regionlist=domaine);  
// Tetrahelize the interior of the cube with tetgen  
//medit(“tetg”,Th,wait=1);

fespace Xh(Th,P23d);  
Xh psi, wx,dxpsi,dypsi;  
macro grad(u) [dx(u),dy(u),dz(u)]// def of grad operator  
// Solve  
solve potential(psi, wx)  
= int3d(Th)(  
grad(psi)'_grad(wx)) // scalar product  
//+int2d(Thsrc) (2_wx)  
+int2d(Th,1) (.200_wx)  
+int2d(Th,2) (-.200_wx)  
//+on(Thl,psi=1)  
//+on(102,psi=1)  
+on(1,psi=0)  
//+on(1,psi=0)  
//+on(1,psi=1)  
//+on(3,dx(psi)=0)  
//+on(4,psi=1)  
//+on(5,psi=1)  
//+on(6,psi=1)  
;

dxpsi=dx(psi);  
dypsi=dy(psi);  
{  
ofstream ff(“psi.txt”);  
for (int ix = 0; ix \< Th.nv; ix++)  
{  
// this does produce dups, maybe confusing everything…  
for (int j = 0; j \< 3; j++) // foo[ix]+=Thx[ix][j]\*Thx[ix][j];  
//ff \<\< Th[ix][j].x \<\< " “\<\< Th[ix][j].y \<\< " " \<\< psi[Xh(ix,j)] \<\< endl;  
ff \<\< Th[ix][j].x \<\< " “\<\< Th[ix][j].y \<\< " "  
\<\<Th[ix][j].z\<\<” "  
\<\< dxpsi[Xh(ix,j)]\<\<” “\<\<dypsi[Xh(ix,j)]\<\<” "\<\<psi[Xh(ix,j)] \<\< endl;  
//ff \<\< Th[ix][0].x \<\< " " \<\< Th[ix][0].y \<\< " " \<\< psi[Xh(ix,0)] \<\< “\n\n\n”;  
}  
} // scoping

//real[int] viso=[1.2,1.15,1.1,1.02,1.01];  
//real[int] viso=[2,1.1,1.05,1.01,1,.999,.99,.9];  
//real[int] viso=[1.7,1.1,1,.999,.99,.9];  
real[int] viso=[0,.01,.05,.1,.2,.3];  
plot(psi, viso=viso, value=true, wait=1);  
//plot(psi, value=true, wait=1);  
//real[int] xxx=[dx(psi),dy(psi)];  
//plot(xxx, value=true, wait=1);

R-code ,

df\<-read.table(“psi.txt”,header=F)  
library(“rgl”)  
x=which((df$V3\> -.9)&(df$V3\< -.8))  
plot3d(df$V1,df$V2,df$V4)  
plot3d(df$V1,df$V2,df$V5,add=T, col=“red”)  
plot3d(df$V1,df$V2,df$V6)  
plot3d(df$V1,df$V2,df$V6)  
history()  
?history  
savehistory(file=“xxx.hist”)

---

<div class="post-metadata">

**Author:** ![marchywka](https://avatars.discourse-cdn.com/v4/letter/m/ee59a6/32.png) [@marchywka](https://community.freefem.org/u/marchywka)\
**Post date:** [January 13, 2023, 11:59pm UTC](https://community.freefem.org/t/region-and-label/2217/2 "2023-01-13T23:59:40Z")

</div>

This is probably behaving as described in the manual 🙂  
I can probably figure out the variable dump format …

---

<div class="post-metadata">

**Author:** ![julienG](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/julieng/32/222_2.png) [@julienG](https://community.freefem.org/u/julienG)\
**Post date:** [January 14, 2023, 10:48am UTC](https://community.freefem.org/t/region-and-label/2217/3 "2023-01-14T10:48:35Z")

</div>

> [@marchywka](#):
>
> moving on the Dirac or Landau-Lifschitz- Gilbert. 🙂  
> Has anyone encountered issues with weak formulations that may be in the examples?

No big problem with the weak forms. But for LLG you’ll have hard times to figure out how to work with a vector field with constant norm. At least I did not find a convenient way to deal wit this…

---

<div class="post-metadata">

**Author:** ![marchywka](https://avatars.discourse-cdn.com/v4/letter/m/ee59a6/32.png) [@marchywka](https://community.freefem.org/u/marchywka)\
**Post date:** [January 14, 2023, 10:59am UTC](https://community.freefem.org/t/region-and-label/2217/4 "2023-01-14T10:59:11Z")

</div>

angles? Thanks. There does seem to be a lot of literature on numerical  
analysis of LLG. Since you seem familiar with it, curious what manifestations  
there may be of spin inertia on linger time scales.

One other thing I was considering was modelling coagulation or maybe some  
classical electrochemistry ( drift and diffusion and reactions etc ). Just offhand  
I did not expect any issues, just a learning curve.

---

<div class="post-metadata">

**Author:** ![julienG](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/julieng/32/222_2.png) [@julienG](https://community.freefem.org/u/julienG)\
**Post date:** [January 16, 2023, 9:53am UTC](https://community.freefem.org/t/region-and-label/2217/5 "2023-01-16T09:53:44Z")

</div>

Well, I think angles is a complicated option because angle are compact variables. So this means you can have discontinuity of you angular variable due to the identification between zero and 2\pi…

What people usually do is to work in Cartesian coordinates for the vector field. More precisely, the magnetization say \boldsymbol{m} that has value on S^2, that is \boldsymbol{m}\cdot\boldsymbol{m}=1 everywhere. In general, algorithms will not preserve the norm of \boldsymbol{m}, so typically after each update, you have to project back \boldsymbol{m} on S^2. Computation-wise, this is quite time consuming I think…

That said, whenever you drop the constraint of normalization, everything is fine 🙂

---

<div class="post-metadata">

**Author:** ![marchywka](https://avatars.discourse-cdn.com/v4/letter/m/ee59a6/32.png) [@marchywka](https://community.freefem.org/u/marchywka)\
**Post date:** [January 16, 2023, 11:04am UTC](https://community.freefem.org/t/region-and-label/2217/6 "2023-01-16T11:04:48Z")

</div>

Thanks. I’ll have to look at it some more but does theta actually appear anywhere  
or just cos and sin? Overall tracking rotations will likely come up somewhere 🙂  
Is norm conservation fundamentally different from energy conservation?

It will probably be a while before I get there,. Right now still playing with navier  
stokes, thinking about coagulation or some eletrochemistry or a simple  
configuration of metal ( like LC tank or antenna ) and solving for vector potential  
A and rho. Previously I had been playing with libmesh but this  
seems simpler for starting out and I was also looking at  
some front end code for FEM and this looks like a good example  
of that if nothing else.  
Thanks.
