Toy Maxwell equation code

General Discussion
1

Hi, I made a FreeFem code to test a toy 2D MaxWell equation model which is proposed as a benchmark in the presentation by Monique Dauge & Martin Costabel
( https://perso.univ-rennes1.fr/monique.dauge/publis/MaxwellMD.pdf )

Correct eigenvalues seem independent of the div . div regularization term s while spurious eigenvalues disperse linearly - which is the expected behavior from the presentation.

image
image1918×966 33.1 KB

The FreeFem++ script is listed below.

/////////////////////////////////////

int Left=99;
int Right=98;
int Top=96;
int Bottom=97;

border D1(t=0., pi){x=t; y=0; label=Bottom;}
border D2(t=0, pi){x=pi; y=t; label=Right;}
border D3(t=pi, 0){x=t; y=pi; label=Top;}
border D4(t=pi, 0){x=0; y=t; label=Left;}

int n = 5;
mesh Th=buildmesh( D1(10n) + D2(10n) + D3(10n) + D4(10n) );

fespace Vh(Th, [P2, P2]); // change to P2 if convergence problems

Vh [Ex, Ey];
Vh [Ux, Uy];
real sigma = 1.;
real s;

ofstream fout(“evalout.dat”);

for (s = 0; s <= 4.01; s += 0.05) {

varf TE([Ex, Ey], [Ux, Uy])= int2d(Th)( (dx(Uy) - dy(Ux)) * (dx(Ey) - dy(Ex) ) )

  • int2d(Th)( s * (dx(Ux) + dy(Uy)) * (dx(Ex) + dy(Ey) ) )
    • int2d(Th)( sigma * ( Ex * Ux + Ey * Uy ) )
      +on(Left,Ey=0) + on(Right, Ey=0) + on(Top, Ex=0) + on(Bottom, Ex = 0)
      ;

varf TEid([Ex, Ey], [Ux, Uy])=int2d(Th)( Ex * Ux + Ey * Uy);

int nev = 30;
real[int] ev(nev);

matrix A = TE(Vh, Vh,solver=UMFPACK);
matrix B = TEid(Vh, Vh, solver=UMFPACK);

// Vh[int] eVx, eVy;
// int k = EigenValue(A, B, sym=true, sigma=sigma, value=ev, vector=eVx, tol=1e-10, maxit=3, ncv=0);

int k = EigenValue(A, B, sym=true, sigma=sigma, value=ev, tol=1e-10);

for (int i = 0; i < k; i++){
fout << s << " " << i << " " << ev[i] << endl;
}
fout << endl;
fout.flush;

//plot( eVx[25], fill=true, value=true, wait= 1);
//plot( eVy[25], fill=true, value=true, wait= 1);

}