Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- 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
Advertisement
Add Comment
Please, Sign In to add comment