do script " cd '.'; /usr/local/ff++/mpich3/bin/FreeFem++ 'test.edp' -glut '/usr/local/ff++/mpich3/bin/ffglut' ;"
-- FreeFem++ v4.9 (Mer 21 avr 2021 16:09:39 CEST - git v4.9)
 Load: lg_fem lg_mesh lg_mesh3 eigenvalue 
    1 : load "MUMPS"init MUMPS_SEQ: MPI_Init

    2 : load "iovtk" load: iovtk 

    3 : load "medit"
    4 : load "msh3"
    5 : load "Element_Mixte3d"
    6 : 
    7 : bool debug = true;
    8 : 
    9 : real r    = 0.00127;
   10 : real d    = r;
   11 : int  nx   = 20;
   12 : int  nz   = 20;
   13 : real zmin = 0.0;
   14 : real zmax = 0.01;
   15 : 
   16 : // build the mesh
   17 : border C(t=0.0, 2.0*pi) {
   18 :     x=r*cos(t);
   19 :     y=r*sin(t);
   20 :     label=1;
   21 : };
   22 : mesh Thx = buildmesh(C(nx));
   23 : 
   24 : int [int] rup = [0, 1], rdown = [0, 2], rmid = [1, 3];
   25 : int guide = 3, in = 1, out = 2;
   26 : mesh3 Th = buildlayers(Thx, nz, zbound=[zmin, zmax]
   27 :         , coef = 1.
   28 :         , labelup = rup
   29 :         , labelmid = rmid
   30 :         , labeldown = rdown
   31 :         );
   32 : //medit("Cyl", Th, wait=debug);
   33 : //plot(Th, wait=debug);
   34 : 
   35 : cout << "labels(Th) = " << labels(Th) << endl;
   36 : 
   37 : real f  = 94*10^9;
   38 : real er = 1;
   39 : real c0 = 299792458;
   40 : real k  = 2*pi*f*sqrt(er) / c0;
   41 : int m = 1, n = 1;
   42 : real beta = sqrt(k^2 - (m*pi/d)^2 - (n*pi/d)^2);
   43 : real Z = sqrt(er) * 120 * pi;
   44 : 
   45 : real ukb = 1 / (k^2 - beta^2);
   46 : func expbz = exp(-1i*beta*z);
   47 : func ExTE = (1i*k*Z) * ukb * (n*pi) / d * cos(m*pi*x/d) * sin(n*pi*y/d) * expbz;
   48 : func EyTE = -(1i*k*Z) * ukb * (n*pi) / d * sin(m*pi*x/d) * cos(n*pi*y/d) * expbz;
   49 : 
   50 : func Gix = +2*1i*beta*ExTE;
   51 : func Giy = +2*1i*beta*EyTE;
   52 : func Giz = 0;
   53 : 
   54 : 
   55 : fespace Nh(Th, Edge13d);
   56 : Nh<complex> [Ex, Ey, Ez], [vx, vy, vz];
   57 : 
   58 : // Macros
   59 : macro Curl(ux,uy,uz) [dy(uz)-dz(uy),dz(ux)-dx(uz),dx(uy)-dy(ux)]  )  // EOM
   60 : macro Nvec(ux,uy,uz) [uy*N.z-uz*N.y,uz*N.x-ux*N.z,ux*N.y-uy*N.x]  )  // EOM //uxN 
   61 : macro Curlabs(ux,uy,uz) [abs(dy(uz)-dz(uy)),abs(dz(ux)-dx(uz)),abs(dx(uy)-dy(ux))]  )  //EOM
   62 : 
   63 : problem waveguide([Ex, Ey, Ez], [vx, vy, vz], solver=sparsesolver) =
   64 :       int3d(Th)(Curl(Ex,Ey,Ez)   [dy(Ez)-dz(Ey),dz(Ex)-dx(Ez),dx(Ey)-dy(Ex)] ' * Curl(vx,vy,vz)   [dy(vz)-dz(vy),dz(vx)-dx(vz),dx(vy)-dy(vx)] )
   65 :     - int3d(Th)(k^2*[Ex,Ey,Ez]' * [vx, vy, vz])
   66 :     + int2d(Th,in,out)(1i * beta * Nvec(Ex,Ey,Ez)   [Ey*N.z-Ez*N.y,Ez*N.x-Ex*N.z,Ex*N.y-Ey*N.x] ' * Nvec(vx,vy,vz)   [vy*N.z-vz*N.y,vz*N.x-vx*N.z,vx*N.y-vy*N.x] )
   67 :     - int2d(Th,in)([vx,vy,vz]' * [Gix,Giy,Giz])
   68 :     + on(guide,Ex=0,Ey=0,Ez=0);
   69 : 
   70 : waveguide;
   71 : 
   72 : cout << "r = " << r << " Integrated Ex = " << int3d(Th)(abs(Ex)^2) << endl;
   73 : medit("real", Th, [real(Ex),real(Ey),real(Ez)]);
   74 : //cout << (abs(Ex)^2) << end;
   75 : //cout << (abs(Ey)^2) << end;
   76 : //cout << (abs(Ez)^2) << end;
   77 : 
   78 : //int[int] Order = [1];
   79 : //savevtk("result_x.vtu", Th, abs(real(Ex)), dataname="RealEx", order=Order);
   80 : //savevtk("result_y.vtu", Th, abs(real(Ey)), dataname="RealEy", order=Order);
   81 : //savevtk("result_z.vtu", Th, abs(real(Ez)), dataname="RealEz", order=Order);
   82 : //medit("real", Th, [real(Ex),real(Ey),real(Ez)]);
   83 :  sizestack + 1024 =9384  ( 8360 )

  --  mesh:  Nb of Triangles =     68, Nb of Vertices 45
labels(Th) = 3	
	  1	  2	  3	
  -- Build Nodes/DF on mesh :   n.v. 945, n. elmt. 4080, n b. elmt. 936
     nb of Nodes 14120    nb of DoF   28240  DFon=0220
  -- FESpace: Nb of Nodes 14120 Nb of DoF 28240
  -- Solve : 
          min (0,0)  max (0,0)
  -- Medit,  Release 3.0a (Nov. 30, 2007)
     Copyright (c) LJLL, 1999-2007.
     compiled:  Mer 20 mai 2020 16:08:51 CEST (with ff++ 4.6).

medit with binary version of popen : mesh(es) and solution(s) 
mesh_name= real
 Loading data file(s)
End of mesh
.sol: Dimension 3 (mesh)3 (lecture)1 
  Input seconds:     0.00

 medit1() 

 Building scene(s)
  Creating scene 1
   Loading default options
  Scene seconds:     0.03

 Rendering scene(s)

 Total running seconds:  0.40
 Thank you for using Medit.
r = 0.00127 Integrated Ex = nan
version de medit ffmedit -popen -filebin -addsol 1  real
times: compile 0.022873s, execution 8.47454s,  mpirank:0
 ######## We forget of deleting   6 Nb pointer,   0Bytes  ,  mpirank 0, memory leak =44677024
 CodeAlloc : nb ptr  4277,  size :496712 mpirank: 0
Ok: Normal End
close  MUMPS_SEQ: MPI_Finalize