From 2c456f7bc87860313300b8b7ab458ffba8e00a7f Mon Sep 17 00:00:00 2001
From: Guillaume Demesy <guillaume.demesy@fresnel.fr>
Date: Mon, 18 Nov 2019 18:30:44 +0100
Subject: [PATCH] try bi-sin surface

---
 DiffractionGratings/grating3D.geo             | 58 +++++++++++++++++--
 DiffractionGratings/grating3D.pro             |  4 +-
 DiffractionGratings/grating3D_data_torus.geo  | 38 ++++++------
 .../gratings_common_materials.pro             |  2 +-
 4 files changed, 76 insertions(+), 26 deletions(-)

diff --git a/DiffractionGratings/grating3D.geo b/DiffractionGratings/grating3D.geo
index 6a9ade2..26a7400 100644
--- a/DiffractionGratings/grating3D.geo
+++ b/DiffractionGratings/grating3D.geo
@@ -1,7 +1,8 @@
-Include "grating3D_data_torus.geo";
+// Include "grating3D_data_torus.geo";
 // Include "grating3D_data_hole.geo";
 // Include "grating3D_data_pyramid.geo";
 // Include "grating3D_data_halfellipsoid.geo";
+Include StrCat["grating3D_data_",test_case,".geo"];
 
 SetFactory("OpenCASCADE");
 // Geometry.LineWidth = 5;
@@ -41,7 +42,9 @@ Box(1) = {-period_x/2,-period_y/2, PML_bot_hh,period_x,period_y, PML_bot  };
 Box(2) = {-period_x/2,-period_y/2,     hh_L_6,period_x,period_y, thick_L_6};
 Box(3) = {-period_x/2,-period_y/2,     hh_L_5,period_x,period_y, thick_L_5};
 Box(4) = {-period_x/2,-period_y/2,     hh_L_4,period_x,period_y, thick_L_4};
-Box(5) = {-period_x/2,-period_y/2,     hh_L_3,period_x,period_y, thick_L_3};
+If (tag_geom!=6)
+    Box(5) = {-period_x/2,-period_y/2,     hh_L_3,period_x,period_y, thick_L_3};
+EndIf
 Box(6) = {-period_x/2,-period_y/2,     hh_L_2,period_x,period_y, thick_L_2};
 Box(7) = {-period_x/2,-period_y/2,     hh_L_1,period_x,period_y, thick_L_1};
 Box(8) = {-period_x/2,-period_y/2, PML_top_hh,period_x,period_y, PML_top  };
@@ -87,9 +90,55 @@ If (tag_geom==4)
     BooleanDifference{ Volume{9}; Delete; }{ Volume{4}; }
 EndIf
 
+If (tag_geom==5)
+    Box(9) = {-rx/2,-ry/2, hh_L_2-rz, rx, ry, rz};
+    Rotate {{0, 0, 1}, {0,0,0}, Pi/4} {Volume{9};}
+EndIf
 
-Coherence;
+If (tag_geom==6)
+    nsin = 2;
+    pp=newp;
+    ll=newl;
+    ss=news;
+    For i In {0:nsin}
+        xx~{i} = -period_x/2+i*period_x/(nsin);
+        yy~{i} = -period_y/2+i*period_y/(nsin);
+    EndFor
+    lss={};
+    For i In {0:nsin-1}
+        For j In {0:nsin-1}
+            pp=newp;ll=newl;lll=newll;
+            ssij=ss+i*nsin+j;
+            Point(pp+0) = {xx~{i},yy~{j},rz/4*(Cos(2*Pi*xx~{i}/period_x)+Cos(2*Pi*yy~{j}/period_y))+hh_L_3+thick_L_3/2};
+            Point(pp+1) = {xx~{i+1},yy~{j},rz/4*(Cos(2*Pi*xx~{i+1}/period_x)+Cos(2*Pi*yy~{j}/period_y))+hh_L_3+thick_L_3/2};
+            Point(pp+2) = {xx~{i+1},yy~{j+1},rz/4*(Cos(2*Pi*xx~{i+1}/period_x)+Cos(2*Pi*yy~{j+1}/period_y))+hh_L_3+thick_L_3/2};
+            Point(pp+3) = {xx~{i},yy~{j+1},rz/4*(Cos(2*Pi*xx~{i}/period_x)+Cos(2*Pi*yy~{j+1}/period_y))+hh_L_3+thick_L_3/2};
+            Line(ll+0)  = {pp,pp+1};
+            Line(ll+1)  = {pp+1,pp+2};
+            Line(ll+2)  = {pp+2,pp+3};
+            Line(ll+3)  = {pp+3,pp};
+            // If(i==0)
+            //     Line{ll} In Surface{25};
+            // EndIf            
+            Line Loop(lll) = {ll,ll+1,ll+2,ll+3};
+            Surface(ssij) = {lll};
+            lss(i*nsin+j) = ssij;
+            // Surface{ss} In Volume{5};
+        EndFor
+    EndFor
+    BooleanUnion { Surface{lss(0)}; Delete; }{ Surface{lss({1:#lss()-1})}; Delete; }
+    
+    // Line(305) = {25, 57};
+    // Line(306) = {57, 34};
+    // Line(307) = {29, 167};
+    // Line(308) = {167, 38};
+    // Line(309) = {31, 177};
+    // Line(310) = {177, 40};
+    // Line(311) = {36, 78};
+    // Line(312) = {78, 27};
 
+EndIf
+Coherence;
 // // // Physical Surface("ENDPML",770) = {69,37};
 
 Call SetPBCs;
@@ -141,4 +190,5 @@ Characteristic Length{pts_ROD()}      =  scat_lc;
 If (tag_geom==3)
     Mesh.Algorithm = 6;
 EndIf
-Mesh.Lines = 0;
+Mesh.SurfaceEdges = 0;
+Mesh.VolumeEdges = 0;
diff --git a/DiffractionGratings/grating3D.pro b/DiffractionGratings/grating3D.pro
index 09e9235..f6e0ad8 100644
--- a/DiffractionGratings/grating3D.pro
+++ b/DiffractionGratings/grating3D.pro
@@ -1,5 +1,5 @@
 myDir = "run_results3D/";
-Include "grating3D_data_torus.geo";
+Include StrCat["grating3D_data_",test_case,".geo"];
 // Include "grating3D_data_hole.geo";
 // Include "grating3D_data_pyramid.geo";
 // Include "grating3D_data_halfellipsoid.geo";
@@ -393,7 +393,7 @@ PostOperation {
       EndIf
       // // For DEBUG
       Print [ epsr_xx , OnElementsOf Omega, File StrCat[myDir,"epsr_xx.pos"]];
-      // Print [ epsr_xx , OnElementsOf Scat, File StrCat[myDir,"epsr_xx.pos"]];
+      // Print [ epsr_0xx , OnElementsOf Scat, File StrCat[myDir,"epsr_xx.pos"]];
       // Print [ Ecm , OnElementsOf Omega_plot, File StrCat[myDir,"Ecm.pos"]];
       // Print [ E1, OnElementsOf Omega_plot, File StrCat[myDir,"E1.pos"]];
       // Print [ Damp_pml_top, OnElementsOf PMLtop, File StrCat[myDir,"Damp_pml_top.pos"]];
diff --git a/DiffractionGratings/grating3D_data_torus.geo b/DiffractionGratings/grating3D_data_torus.geo
index 69f8f8e..a5216cd 100644
--- a/DiffractionGratings/grating3D_data_torus.geo
+++ b/DiffractionGratings/grating3D_data_torus.geo
@@ -6,26 +6,26 @@ pp4 = "4Layer Materials";
 pp5 = "5Computational Paramameters";
 pp6 = "6Output";
 DefineConstant[
-    lambda0       = {1000   , Name StrCat[pp1,"/1lambda0 [nm]"]},
+    lambda0       = {1000  , Name StrCat[pp1,"/1lambda0 [nm]"]},
     thetadeg      = {0     , Name StrCat[pp1,"/2theta0 [deg]"]},
     phideg        = {0     , Name StrCat[pp1,"/3phi0 [deg]"]},
     psideg        = {0     , Name StrCat[pp1,"/4psi0 [deg]"]},
     period_x      = { 300  , Name StrCat[pp2,"/1X period [nm]"]},
     period_y      = { 300  , Name StrCat[pp2,"/2Y period [nm]"]},
     thick_L_1     = {50    , Name StrCat[pp2,"/3thickness layer 1 [nm] (superstrate)"]},
-    thick_L_2     = {200    , Name StrCat[pp2,"/4thickness layer 2 [nm]"]},
+    thick_L_2     = {200   , Name StrCat[pp2,"/4thickness layer 2 [nm]"]},
     thick_L_3     = {100   , Name StrCat[pp2,"/5thickness layer 3 [nm]"]},
-    thick_L_4     = {200    , Name StrCat[pp2,"/6thickness layer 4 [nm]"]},
+    thick_L_4     = {200   , Name StrCat[pp2,"/6thickness layer 4 [nm]"]},
     thick_L_5     = {50    , Name StrCat[pp2,"/7thickness layer 5 [nm]"]},
     thick_L_6     = {50    , Name StrCat[pp2,"/8thickness layer 6 [nm] (substrate)"]},
     
-    tag_geom      = {  3   , Name StrCat[pp3,"/0Shape"], Choices {1="Pyramid",2="Cylindrical Hole",3="Torus",4="HalfEllipspoid"}},
+    tag_geom      = {  3   , Name StrCat[pp3,"/0Shape"], Choices {1="Pyramid",2="Cylindrical Hole",3="Torus",4="HalfEllipspoid",5="Checkerboard",6="bi-sinusoidal"}},
     rx            = {150/2 , Name StrCat[pp3,"/1rx"]},
     ry            = { 50   , Name StrCat[pp3,"/2ry"]},
     rz            = { 25   , Name StrCat[pp3,"/3rz"]},
     flag_mat_scat = { 0    , Name StrCat[pp3,"/4Scatterer permittivity model"], Choices {0="Custom (Value Below)",1="SiO2",2="Ag (palik)",3="Al (palik)",4="Au (johnson)",5="Nb2O5",6="ZnSe",7="MgF2",8="TiO2",9="PMMA",10="Si",11="ITO",12="Cu (palik)"} },
-    eps_re_Scat   = {-21   , Name StrCat[pp3,"/5Custom eps_re_Scat"]},
-    eps_im_Scat   = { 20   , Name StrCat[pp3,"/6Custom eps_im_Scat"]},
+    eps_re_Scat   = {-21   , Name StrCat[pp3,"/5Custom real part of relative permittivity"]},
+    eps_im_Scat   = { 20   , Name StrCat[pp3,"/6Custom real part of relative permittivity"]},
 
     flag_mat_1    = { 0    , Name StrCat[pp4,"/1Layer 1"], Choices {0="Custom (Value Below)",1="SiO2",2="Ag (palik)",3="Al (palik)",4="Au (johnson)",5="Nb2O5",6="ZnSe",7="MgF2",8="TiO2",9="PMMA",10="Si",11="ITO",12="Cu (palik)"} },
     flag_mat_2    = { 0    , Name StrCat[pp4,"/2Layer 2"], Choices {0="Custom (Value Below)",1="SiO2",2="Ag (palik)",3="Al (palik)",4="Au (johnson)",5="Nb2O5",6="ZnSe",7="MgF2",8="TiO2",9="PMMA",10="Si",11="ITO",12="Cu (palik)"} },
@@ -46,17 +46,17 @@ DefineConstant[
     eps_re_L_6    = {2.25  , Name StrCat[pp4,"/7Custom Values/layer 6: real part of relative permittivity"]},
     eps_im_L_6    = {0     , Name StrCat[pp4,"/7Custom Values/layer 6: imag part of relative permittivity"]},
     
-    paramaille    = {10         , Name StrCat[pp5,"/1Number of mesh elements per wavelength [-]"]},
-    scat_lc       = {20         , Name StrCat[pp5,"/2Scatterer absolute mesh size [nm]"]},
-    PML_top       = {lambda0*1.5, Name StrCat[pp5,"/4PML top thickness [nm]"]},
-    PML_bot       = {lambda0*1.5, Name StrCat[pp5,"/5PML bot thickness [nm]"]},
-    Nmax          = {1          , Name StrCat[pp5,"/6Max number of non specular order to output [-]"]},
-    refine_mesh_L1= {1          , Name StrCat[pp5,"/7refine layers/1refine mesh layer 1 [-]"]},
-    refine_mesh_L2= {1          , Name StrCat[pp5,"/7refine layers/2refine mesh layer 2 [-]"]},
-    refine_mesh_L3= {1          , Name StrCat[pp5,"/7refine layers/3refine mesh layer 3 [-]"]},
-    refine_mesh_L4= {1          , Name StrCat[pp5,"/7refine layers/4refine mesh layer 4 [-]"]},
-    refine_mesh_L5= {1          , Name StrCat[pp5,"/7refine layers/5refine mesh layer 5 [-]"]},
-    refine_mesh_L6= {1          , Name StrCat[pp5,"/7refine layers/6refine mesh layer 6 [-]"]},
+    paramaille    = {10      , Name StrCat[pp5,"/1Number of mesh elements per wavelength [-]"]},
+    scat_lc       = {20      , Name StrCat[pp5,"/2Scatterer absolute mesh size [nm]"]},
+    PML_top       = {lambda0 , Name StrCat[pp5,"/4PML top thickness [nm]"]},
+    PML_bot       = {lambda0 , Name StrCat[pp5,"/5PML bot thickness [nm]"]},
+    Nmax          = {1       , Name StrCat[pp5,"/6Max number of non specular order to output [-]"]},
+    refine_mesh_L1= {1       , Name StrCat[pp5,"/7refine layers/1refine mesh layer 1 [-]"]},
+    refine_mesh_L2= {1       , Name StrCat[pp5,"/7refine layers/2refine mesh layer 2 [-]"]},
+    refine_mesh_L3= {1       , Name StrCat[pp5,"/7refine layers/3refine mesh layer 3 [-]"]},
+    refine_mesh_L4= {1       , Name StrCat[pp5,"/7refine layers/4refine mesh layer 4 [-]"]},
+    refine_mesh_L5= {1       , Name StrCat[pp5,"/7refine layers/5refine mesh layer 5 [-]"]},
+    refine_mesh_L6= {1       , Name StrCat[pp5,"/7refine layers/6refine mesh layer 6 [-]"]},
     
     InterpSampling     = { 30   , Name StrCat[pp6,"/0Interpolation grid step [nm]"]},
     Flag_interp_cubic  = { 0    , Name StrCat[pp6,"/1Interpolate on cubic grid?"], Choices {0,1} },
@@ -78,8 +78,8 @@ PML_bot_hh = hh_L_6-PML_bot;
 PML_top_hh = hh_L_1+thick_L_1;
 
 theta0 = thetadeg*Pi/180;
-phi0 = phideg*Pi/180;
-psi0 = psideg*Pi/180;
+phi0   = phideg*Pi/180;
+psi0   = psideg*Pi/180;
 
 DomainZsizeSca  = PML_top_hh+PML_bot-(hh_L_6-PML_bot);
 DomainZsizeTot  = PML_top_hh-hh_L_6;
diff --git a/DiffractionGratings/gratings_common_materials.pro b/DiffractionGratings/gratings_common_materials.pro
index 2ecd9b2..df21bcf 100644
--- a/DiffractionGratings/gratings_common_materials.pro
+++ b/DiffractionGratings/gratings_common_materials.pro
@@ -31,7 +31,7 @@
 // 11_Si (odp)
 // 12_ITO
 // 13_Cu
-nb_available_materials = 13;
+nb_available_materials = 12;
 lambdamat_1    = {1.984000e-07,2.066000e-07,2.139000e-07,2.144000e-07,2.267000e-07,2.302000e-07,2.378000e-07,2.399000e-07,2.483000e-07,2.652000e-07,2.699000e-07,2.753000e-07,2.803000e-07,2.894000e-07,2.967000e-07,3.021000e-07,3.303000e-07,3.341000e-07,3.404000e-07,3.466000e-07,3.611000e-07,3.650000e-07,4.047000e-07,4.358000e-07,4.678000e-07,4.861000e-07,5.086000e-07,5.461000e-07,5.770000e-07,5.791000e-07,5.876000e-07,5.893000e-07,6.438000e-07,6.563000e-07,6.678000e-07,7.065000e-07,8.521000e-07,8.943000e-07,1.014000e-06,1.083000e-06,1.128700e-06,1.362200e-06,1.395100e-06,1.469500e-06,1.529500e-06,1.660600e-06,1.681000e-06,1.693200e-06,1.709100e-06,1.813100e-06,1.970100e-06,2.058100e-06,2.152600e-06,2.325400e-06,2.437400e-06,3.243900e-06};
 epsr_mat_re_1  = {2.414916e+00,2.380849e+00,2.354076e+00,2.352236e+00,2.318920e+00,2.310704e+00,2.294316e+00,2.290380e+00,2.275271e+00,2.250000e+00,2.244304e+00,2.237717e+00,2.232036e+00,2.223081e+00,2.216228e+00,2.211764e+00,2.191880e+00,2.189808e+00,2.186258e+00,2.183006e+00,2.175920e+00,2.174150e+00,2.159724e+00,2.151209e+00,2.144174e+00,2.140662e+00,2.137152e+00,2.131892e+00,2.128389e+00,2.128097e+00,2.127222e+00,2.126931e+00,2.121975e+00,2.121101e+00,2.120227e+00,2.117316e+00,2.109756e+00,2.108014e+00,2.103370e+00,2.100760e+00,2.099311e+00,2.091494e+00,2.090338e+00,2.088025e+00,2.086002e+00,2.081383e+00,2.080518e+00,2.080229e+00,2.079652e+00,2.075616e+00,2.069282e+00,2.065544e+00,2.061522e+00,2.053202e+00,2.047475e+00,1.996852e+00};
 epsr_mat_im_1  = {-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00,-0.000000e+00};
-- 
GitLab