normalized units : 2*pi*c=a so we find exactly

as in the convergence section of the paper

NonLinearEVP/NonLinearEVP_data.geo

@@ -3,10 +3,6 @@
 ////   _data.geo               ////
 ///////////////////////////////////
 
-nm       = 1.;
-epsilon0 = 8.854187817e-3*nm;
-mu0      = 400.*Pi*nm;
-cel      = 1.0/(Sqrt[epsilon0 * mu0]);
 deg2rad  = Pi/180;
 
 pp0        = "1Geometry/0";
@@ -27,7 +23,6 @@ colorpp5   = "Ivory";
 DefineConstant[ a_lat = {50     , Name StrCat[pp0  , "1grating period d [nm]"]          , Highlight Str[colorpp0]  , Closed close_menu} ];
 
 // normalization factor
-norm = a_lat/(2.*Pi*cel);
 
 DefineConstant[
   d_sq           = {0.806 , Name StrCat[pp0 , "2sq [d]"] , Highlight Str[colorpp0]  , Closed close_menu} , 
@@ -58,6 +53,15 @@ DefineConstant[
   flag_outEigvec = {1     , Name StrCat[pp4, "output eigenvector?"], Choices{0,1}}
 ];
 
+
+cel      = a_lat/(2*Pi);
+epsf     = 8.854187817e-3;
+muf      = 400.*Pi;
+nm       = 2*Pi/(a_lat*Sqrt[epsf*muf]);
+epsilon0 = epsf*nm;
+mu0      = muf*nm;
+norm     = a_lat/(2.*Pi*cel);
+
 // Normalization
 d_sq           = d_sq          * a_lat;
 space2pml      = space2pml     * a_lat;