Skip to content
Snippets Groups Projects
Select Git revision
  • b43da7130702982e1f4f7baba94f5084191d9742
  • master default protected
  • albertpiwonski-master-patch-57409
  • quadspheres
  • fix_Tmatrix_code_epsr_background
  • albertpiwonski-master-patch-12427
  • cavity
  • c1
8 results

NonLinearEVP.geo

Blame
  • NonLinearEVP.geo 15.59 KiB
    // Copyright (C) 2018 Guillaume Demésy
    //
    // This file is part of the model NonLinearEVP.pro.
    // 
    // This program is free software: you can redistribute it and/or modify
    // it under the terms of the GNU General Public License as published by
    // the Free Software Foundation, either version 3 of the License, or
    // (at your option) any later version.
    // 
    // This program is distributed in the hope that it will be useful,
    // but WITHOUT ANY WARRANTY; without even the implied warranty of
    // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    // GNU General Public License for more details.
    // 
    // You should have received a copy of the GNU General Public License
    // along with This program. If not, see <https://www.gnu.org/licenses/>.
    
    Include "NonLinearEVP_data.geo";
    lc_cell = a_lat/paramaille;
    lc_sq   = lc_cell/4;
    lc_pml  = lc_cell*1.2;
    lc_sqa  = lc_sq;
    lc_sq_inside = lc_sq*1.4;
    epsc         = lc_sq*0.8;
    
    
    If (flag_Tmesh==0)
      Point(1)  = {-a_lat/2. ,-d_sq/2.-space2pml, 0. , lc_cell};
      Point(2)  = { a_lat/2. ,-d_sq/2.-space2pml, 0. , lc_cell};
      Point(3)  = { a_lat/2. , d_sq/2.+space2pml, 0. , lc_cell};
      Point(4)  = {-a_lat/2. , d_sq/2.+space2pml, 0. , lc_cell};
      Point(9)  = {-a_lat/2. ,-d_sq/2.-space2pml-pmlsize, 0. , lc_pml};
      Point(10) = { a_lat/2. ,-d_sq/2.-space2pml-pmlsize, 0. , lc_pml};
      Point(11) = { a_lat/2. , d_sq/2.+space2pml+pmlsize, 0. , lc_pml};
      Point(12) = {-a_lat/2. , d_sq/2.+space2pml+pmlsize, 0. , lc_pml};
      Point(13) = {-a_lat/2. , a_lat/2., 0. , lc_sqa};
      Point(14) = {-a_lat/2. ,-a_lat/2., 0. , lc_sqa};
      Point(15) = { a_lat/2. , a_lat/2., 0. , lc_sqa};
      Point(16) = { a_lat/2. ,-a_lat/2., 0. , lc_sqa};
      If (flag_rounding==0)
        Point(5)  = {-d_sq/2. ,-d_sq/2., 0. , lc_sq};
        Point(6)  = { d_sq/2. ,-d_sq/2., 0. , lc_sq};
        Point(7)  = { d_sq/2. , d_sq/2., 0. , lc_sq};
        Point(8)  = {-d_sq/2. , d_sq/2., 0. , lc_sq};
        Line(1) = {1, 2};
        Line(3) = {3, 4};
        Line(5) = {5, 6};
        Line(6) = {6, 7};
        Line(7) = {7, 8};
        Line(8) = {8, 5};
        Line(15) = {1, 14};
        Line(16) = {14, 13};
        Line(17) = {13, 4};
        Line(18) = {2, 16};
        Line(19) = {16, 15};
        Line(20) = {15, 3};
        Line(9) = {4, 12};
        Line(10) = {12, 11};
        Line(11) = {11, 3};
        Line(12) = {1, 9};
        Line(13) = {9, 10};
        Line(14) = {10, 2};
        Line Loop(1) = {12, 13, 14, -1};
        Plane Surface(1) = {1};
        Line Loop(2) = {5, 6, 7, 8};
        Plane Surface(2) = {2};
        Line Loop(3) = {17, -3, -20, -19, -18, -1, 15, 16};
        Plane Surface(3) = {2, -3};
        Line Loop(4) = {9, 10, 11, 3};
        Plane Surface(4) = {-4};
        Periodic Line { 14,18,19,20,11 } = {12,15,16,17,9 } Translate {a_lat,0,0} ;    
        Physical Line("SCATBOUND",1005)   = {5,6,7,8}; // bound for lag
      Else
        // // with a circle
        // Point(20)  = {0.,0., 0. , lc_sq};
        // Point(21)  = {-d_sq/2. ,0, 0. , lc_sq};
        // Point(22)  = {0        ,-d_sq/2., 0. , lc_sq};
        // Point(23)  = {d_sq/2. ,0, 0. , lc_sq};
        // Point(24)  = {0        ,d_sq/2., 0. , lc_sq};
        // Line(1) = {1, 2};
        // Line(3) = {3, 4};
        // Line(15) = {1, 14};
        // Line(16) = {14, 13};
        // Line(17) = {13, 4};
        // Line(18) = {2, 16};
        // Line(19) = {16, 15};
        // Line(20) = {15, 3};
        // Line(9) = {4, 12};
        // Line(10) = {12, 11};
        // Line(11) = {11, 3};
        // Line(12) = {1, 9};
        // Line(13) = {9, 10};
        // Line(14) = {10, 2};
        // Circle(21) = {21, 20, 22};
        // Circle(22) = {22, 20, 23};
        // Circle(23) = {23, 20, 24};
        // Circle(24) = {24, 20, 21};
        // Line Loop(1) = {12, 13, 14, -1};
        // Plane Surface(1) = {1};
        // Curve Loop(2) = {21,22,23,24};
        // Plane Surface(2) = {2};
        // Curve Loop(3) = {20, 3, -17, -16, -15, 1, 18, 19};
        // Plane Surface(3) = {2, 3};
        // Line Loop(4) = {9, 10, 11, 3};
        // Plane Surface(4) = {-4};
        // // Rotate {{0, 0, 1}, {0, 0, 0}, 2.*Pi/180.} { Surface{ 2 } ; }
        // Periodic Line { 14,18,19,20,11 } = {12,15,16,17,9 } Translate {a_lat,0,0} ;
        // Physical Line("SCATBOUND",1005)    = {21,22,23,24}; // bound for lag
    
        Point(20)  = {-d_sq/2.+corner_rad ,-d_sq/2., 0. , lc_sq};
        Point(21)  = {-d_sq/2. ,-d_sq/2.+corner_rad, 0. , lc_sq};
        Point(22)  = {-d_sq/2.+corner_rad ,-d_sq/2.+corner_rad, 0. , lc_sq};
        Point(23)  = { d_sq/2.-corner_rad ,-d_sq/2., 0. , lc_sq};
        Point(24)  = { d_sq/2. ,-d_sq/2.+corner_rad, 0. , lc_sq};
        Point(25)  = { d_sq/2.-corner_rad ,-d_sq/2.+corner_rad, 0. , lc_sq};
        Point(26)  = { d_sq/2.-corner_rad , d_sq/2., 0. , lc_sq};
        Point(27)  = { d_sq/2. , d_sq/2.-corner_rad, 0. , lc_sq};
        Point(28)  = { d_sq/2.-corner_rad , d_sq/2.-corner_rad, 0. , lc_sq};
        Point(29)  = {-d_sq/2.+corner_rad , d_sq/2., 0. , lc_sq};
        Point(30)  = {-d_sq/2. , d_sq/2.-corner_rad, 0. , lc_sq};
        Point(31)  = {-d_sq/2.+corner_rad , d_sq/2.-corner_rad, 0. , lc_sq};
    
        Line(1) = {1, 2};
        Line(3) = {3, 4};
        Line(15) = {1, 14};
        Line(16) = {14, 13};
        Line(17) = {13, 4};
        Line(18) = {2, 16};
        Line(19) = {16, 15};
        Line(20) = {15, 3};
        Line(9) = {4, 12};
        Line(10) = {12, 11};
        Line(11) = {11, 3};
        Line(12) = {1, 9};
        Line(13) = {9, 10};
        Line(14) = {10, 2};
        Circle(21) = {30, 31, 29};
        Circle(22) = {26, 28, 27};
        Circle(23) = {24, 25, 23};
        Circle(24) = {21, 22, 20};
        Line(25) = {29, 26};
        Line(26) = {27, 24};
        Line(27) = {23, 20};
        Line(28) = {21, 30};    
        Line Loop(1) = {12, 13, 14, -1};
        Plane Surface(1) = {1};
        Curve Loop(2) = {25, 22, 26, 23, 27, -24, 28, 21};
        Plane Surface(2) = {2};
        Curve Loop(3) = {20, 3, -17, -16, -15, 1, 18, 19};
        Plane Surface(3) = {2, 3};
        Line Loop(4) = {9, 10, 11, 3};
        Plane Surface(4) = {-4};
        // Rotate {{0, 0, 1}, {0, 0, 0}, 2.*Pi/180.} { Surface{ 2 } ; }
        Periodic Line { 14,18,19,20,11 } = {12,15,16,17,9 } Translate {a_lat,0,0} ;
        Physical Line("SCATBOUND",1005)    = {21,22,23,24,25,26,27,28}; // bound for lag
      EndIf
      Physical Surface("SCAT",100)       = {2};   // 1 dom in
      Physical Surface("OUT",101)        = {3};  // 2 dom out
      Physical Surface("PMLBOT",102)     = {1};  // PML bot
      Physical Surface("PMLTOP",103)     = {4};  // PML top
      Physical Line("BLOCHXL",1001)      = {12,15,16,17,9}; // bloch x left
      Physical Line("BLOCHXR",1002)      = {14,18,19,20,11}; // bloch x right
      Physical Line("TOP",1003)          = {10}; // top bound
      Physical Line("BOT",1004)          = {13}; // bot bound
      Physical Point("PRINTPOINT",10000) = {1};   // Printpoint    
    EndIf
    
    If (flag_Tmesh==1)
      Point(1)  = {-a_lat/2. ,-d_sq/2.-space2pml, 0. , lc_cell};
      Point(2)  = { a_lat/2. ,-d_sq/2.-space2pml, 0. , lc_cell};
      Point(3)  = { a_lat/2. , d_sq/2.+space2pml, 0. , lc_cell};
      Point(4)  = {-a_lat/2. , d_sq/2.+space2pml, 0. , lc_cell};
    
      Point(5)  = {-d_sq/2. ,-d_sq/2., 0. , lc_sq*15};
      Point(6)  = { d_sq/2. ,-d_sq/2., 0. , lc_sq*15};
      Point(7)  = { d_sq/2. , d_sq/2., 0. , lc_sq*15};
      Point(8)  = {-d_sq/2. , d_sq/2., 0. , lc_sq*15};
    
      Point(9)  = {-a_lat/2. ,-d_sq/2.-space2pml-pmlsize, 0. , lc_pml};
      Point(10) = { a_lat/2. ,-d_sq/2.-space2pml-pmlsize, 0. , lc_pml};
      Point(11) = { a_lat/2. , d_sq/2.+space2pml+pmlsize, 0. , lc_pml};
      Point(12) = {-a_lat/2. , d_sq/2.+space2pml+pmlsize, 0. , lc_pml};
    
      Point(13) = {-d_sq/2.-epsc , d_sq/2.+epsc, 0. , lc_sqa};
      Point(14) = {-d_sq/2.      , d_sq/2.+epsc, 0. , lc_sqa};
      Point(15) = {-d_sq/2.+epsc , d_sq/2.+epsc, 0. , lc_sqa};
      Point(16) = { d_sq/2.-epsc , d_sq/2.+epsc, 0. , lc_sqa};
      Point(17) = { d_sq/2.      , d_sq/2.+epsc, 0. , lc_sqa};
      Point(18) = { d_sq/2.+epsc , d_sq/2.+epsc, 0. , lc_sqa};
      Point(19) = { d_sq/2.+epsc , d_sq/2.     , 0. , lc_sqa};
      Point(20) = { d_sq/2.+epsc , d_sq/2.-epsc, 0. , lc_sqa};
      Point(21) = { d_sq/2.+epsc ,-d_sq/2.+epsc, 0. , lc_sqa};
      Point(22) = { d_sq/2.+epsc ,-d_sq/2.     , 0. , lc_sqa};
      Point(23) = { d_sq/2.+epsc ,-d_sq/2.-epsc, 0. , lc_sqa};
      Point(24) = { d_sq/2.      ,-d_sq/2.-epsc, 0. , lc_sqa};
      Point(25) = { d_sq/2.-epsc ,-d_sq/2.-epsc, 0. , lc_sqa};
      Point(26) = {-d_sq/2.+epsc ,-d_sq/2.-epsc, 0. , lc_sqa};
      Point(27) = {-d_sq/2.      ,-d_sq/2.-epsc, 0. , lc_sqa};
      Point(28) = {-d_sq/2.-epsc ,-d_sq/2.-epsc, 0. , lc_sqa};
      Point(29) = {-d_sq/2.-epsc ,-d_sq/2.     , 0. , lc_sqa};
      Point(30) = {-d_sq/2.-epsc ,-d_sq/2.+epsc, 0. , lc_sqa};
      Point(31) = {-d_sq/2.-epsc , d_sq/2.-epsc, 0. , lc_sqa};
      Point(32) = {-d_sq/2.-epsc , d_sq/2.     , 0. , lc_sqa};
    
      Point(33)  = {-d_sq/2.+epsc ,-d_sq/2.+epsc, 0. , lc_sqa};
      Point(34)  = { d_sq/2.-epsc ,-d_sq/2.+epsc, 0. , lc_sqa};
      Point(35)  = { d_sq/2.-epsc , d_sq/2.-epsc, 0. , lc_sqa};
      Point(36)  = {-d_sq/2.+epsc , d_sq/2.-epsc, 0. , lc_sqa};
    
      Point(37)  = {-d_sq/2.+epsc ,-d_sq/2., 0. , lc_sqa};
      Point(38)  = { d_sq/2.-epsc ,-d_sq/2., 0. , lc_sqa};
      Point(39)  = { d_sq/2.-epsc , d_sq/2., 0. , lc_sqa};
      Point(45)  = {-d_sq/2.+epsc , d_sq/2., 0. , lc_sqa};
      Point(46)  = {-d_sq/2. ,-d_sq/2.+epsc, 0. , lc_sqa};
      Point(47)  = { d_sq/2. ,-d_sq/2.+epsc, 0. , lc_sqa};
      Point(48)  = { d_sq/2. , d_sq/2.-epsc, 0. , lc_sqa};
      Point(49)  = {-d_sq/2. , d_sq/2.-epsc, 0. , lc_sqa};
    
      Point(54) = {-a_lat/2. , a_lat/2., 0. , lc_sqa};
      Point(55) = {-a_lat/2. ,-a_lat/2., 0. , lc_sqa};
      Point(56) = { a_lat/2. , a_lat/2., 0. , lc_sqa};
      Point(58) = { a_lat/2. ,-a_lat/2., 0. , lc_sqa};
    
      Line(1) = {1, 2};
      Line(3) = {3, 4};
      Line(9) = {4, 12};
      Line(10) = {12, 11};
      Line(11) = {11, 3};
      Line(12) = {1, 9};
      Line(13) = {9, 10};
      Line(14) = {10, 2};
    
      Line(35) = {8, 45};
      Line(36) = {45, 45};
      Line(37) = {39, 39};
      Line(38) = {45, 39};
      Line(39) = {39, 7};
      Line(40) = {7, 48};
      Line(41) = {48, 47};
      Line(42) = {47, 6};
      Line(43) = {6, 38};
      Line(44) = {38, 37};
      Line(45) = {37, 5};
      Line(46) = {5, 46};
      Line(47) = {46, 49};
      Line(48) = {49, 8};
      Line(49) = {31, 49};
      Line(50) = {49, 36};
      Line(51) = {36, 45};
      Line(52) = {45, 15};
      Line(53) = {32, 8};
      Line(54) = {8, 14};
      Line(55) = {13, 8};
      Line(56) = {8, 36};
      Line(57) = {31, 8};
      Line(58) = {8, 15};
      Line(59) = {16, 39};
      Line(60) = {39, 35};
      Line(61) = {35, 48};
      Line(62) = {48, 20};
      Line(63) = {19, 7};
      Line(64) = {7, 17};
      Line(65) = {16, 7};
      Line(66) = {7, 20};
      Line(67) = {35, 7};
      Line(68) = {7, 18};
      Line(69) = {34, 47};
      Line(70) = {47, 21};
      Line(71) = {34, 38};
      Line(72) = {38, 25};
      Line(73) = {25, 6};
      Line(74) = {6, 24};
      Line(75) = {6, 22};
      Line(76) = {6, 21};
      Line(77) = {23, 6};
      Line(78) = {6, 34};
      Line(79) = {30, 46};
      Line(80) = {46, 33};
      Line(81) = {33, 37};
      Line(82) = {37, 26};
      Line(83) = {26, 5};
      Line(84) = {5, 27};
      Line(85) = {5, 28};
      Line(86) = {5, 29};
      Line(87) = {5, 30};
      Line(88) = {5, 33};
      Line(89) = {33, 36};
      Line(90) = {36, 35};
      Line(91) = {35, 34};
      Line(92) = {34, 33};
      Line(93) = {13, 14};
      Line(94) = {14, 15};
      Line(95) = {15, 16};
      Line(96) = {16, 17};
      Line(97) = {17, 18};
      Line(98) = {18, 19};
      Line(99) = {19, 20};
      Line(100) = {20, 21};
      Line(101) = {21, 22};
      Line(102) = {22, 23};
      Line(103) = {23, 24};
      Line(104) = {24, 25};
      Line(105) = {25, 26};
      Line(106) = {26, 27};
      Line(107) = {27, 28};
      Line(108) = {28, 29};
      Line(109) = {29, 30};
      Line(110) = {30, 31};
      Line(111) = {31, 32};
      Line(112) = {32, 13};
      Line(113) = {1, 55};
      Line(114) = {55, 54};
      Line(115) = {54, 4};
      Line(116) = {2, 58};
      Line(117) = {58, 56};
      Line(118) = {56, 3};
    
      Line Loop(1) = {53, -55, -112};
      Plane Surface(1) = {1};
      Line Loop(2) = {55, 54, -93};
      Plane Surface(2) = {2};
      Line Loop(3) = {54, 94, -58};
      Plane Surface(3) = {-3};
      Line Loop(4) = {35, 52, -58};
      Plane Surface(4) = {4};
      Line Loop(5) = {111, 53, -57};
      Plane Surface(5) = {-5};
      Line Loop(6) = {49, 48, -57};
      Plane Surface(6) = {6};
      Line Loop(7) = {48, 56, -50};
      Plane Surface(7) = {-7};
      Line Loop(8) = {56, 51, -35};
      Plane Surface(8) = {8};
      Line Loop(9) = {51, 38, 60, -90};
      Plane Surface(9) = {-9};
      Line Loop(10) = {38, -59, -95, -52};
      Plane Surface(10) = {10};
      Line Loop(11) = {59, 39, -65};
      Plane Surface(11) = {11};
      Line Loop(12) = {65, 64, -96};
      Plane Surface(12) = {12};
      Line Loop(13) = {64, 97, -68};
      Plane Surface(13) = {-13};
      Line Loop(14) = {68, 98, 63};
      Plane Surface(14) = {-14};
      Line Loop(15) = {63, 66, -99};
      Plane Surface(15) = {15};
      Line Loop(16) = {66, -62, -40};
      Plane Surface(16) = {-16};
      Line Loop(17) = {40, -61, 67};
      Plane Surface(17) = {-17};
      Line Loop(18) = {67, -39, 60};
      Plane Surface(18) = {18};
      Line Loop(19) = {91, 69, -41, -61};
      Plane Surface(19) = {19};
      Line Loop(20) = {70, -100, -62, 41};
      Plane Surface(20) = {20};
      Line Loop(21) = {42, 76, -70};
      Plane Surface(21) = {21};
      Line Loop(22) = {76, 101, -75};
      Plane Surface(22) = {-22};
      Line Loop(23) = {42, 78, 69};
      Plane Surface(23) = {-23};
      Line Loop(24) = {71, -43, 78};
      Plane Surface(24) = {24};
      Line Loop(25) = {72, 73, 43};
      Plane Surface(25) = {25};
      Line Loop(26) = {73, 74, 104};
      Plane Surface(26) = {-26};
      Line Loop(27) = {74, -103, 77};
      Plane Surface(27) = {27};
      Line Loop(28) = {77, 75, 102};
      Plane Surface(28) = {-28};
      Line Loop(29) = {92, 81, -44, -71};
      Plane Surface(29) = {29};
      Line Loop(30) = {82, -105, -72, 44};
      Plane Surface(30) = {30};
      Line Loop(31) = {83, -45, 82};
      Plane Surface(31) = {-31};
      Line Loop(32) = {83, 84, -106};
      Plane Surface(32) = {32};
      Line Loop(33) = {84, 107, -85};
      Plane Surface(33) = {-33};
      Line Loop(34) = {85, 108, -86};
      Plane Surface(34) = {-34};
      Line Loop(35) = {86, 109, -87};
      Plane Surface(35) = {-35};
      Line Loop(36) = {87, 79, -46};
      Plane Surface(36) = {-36};
      Line Loop(37) = {46, 80, -88};
      Plane Surface(37) = {-37};
      Line Loop(38) = {88, 81, 45};
      Plane Surface(38) = {-38};
      Line Loop(39) = {110, 49, -47, -79};
      Plane Surface(39) = {-39};
      Line Loop(40) = {50, -89, -80, 47};
      Plane Surface(40) = {-40};
      Line Loop(43) = {3, 9, 10, 11};
      Plane Surface(43) = {-43};
      Line Loop(44) = {1, -14, -13, -12};
      Plane Surface(44) = {-44};
      Line Loop(45) = {113, 114, 115, -3, -118, -117, -116, -1};
      Line Loop(46) = {106, 107, 108, 109, 110, 111, 112, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105};
      Plane Surface(45) = {-45,-46};
      Line Loop(47) = {89, 90, 91, 92};
      Plane Surface(46) = {-47};
    
    
      Periodic Line { 14,116,117,118,11 } = {12,113, 114, 115,9 } Translate {a_lat,0,0} ;
    
      Periodic Line {47}  = {110} Translate {epsc,0,0} ;
      Periodic Line {89}  = {47}  Translate {epsc,0,0} ;
      Periodic Line {41}  = {91}  Translate {epsc,0,0} ;
      Periodic Line {100} = {41}  Translate {epsc,0,0} ;
    
      Periodic Line {44} = {105} Translate {0,epsc,0} ;
      Periodic Line {92} = {44}  Translate {0,epsc,0} ;
      Periodic Line {38} = {90}  Translate {0,epsc,0} ;
      Periodic Line {95} = {38}  Translate {0,epsc,0} ;
    
      // Transfinite Surface { expression-list } | "*" < = { expression-list } > < Left | Right | Alternate | AlternateRight | AlternateLeft > ;
      Transfinite Surface { 9 } Left;
      Transfinite Surface { 10 } Right;
    
      Transfinite Surface { 19 } Left;
      Transfinite Surface { 20 } Left;
    
      Transfinite Surface { 29 } Left;
      Transfinite Surface { 30 } Left;
    
      Transfinite Surface { 39 } Left;
      Transfinite Surface { 40 } Right;
    
    
      Physical Surface(100) = {7,8,9,17,18,19,23,24,29,37,38,40,46};   // 1 dom in
      Physical Surface(101) = {45, 1, 2, 3, 4, 10, 5, 6, 39, 36, 35, 34, 33, 32, 31, 30, 25, 26, 27, 28, 22, 21, 20, 16, 15, 14, 13, 12, 11}; // out
      Physical Surface(102) = {43};  // PML top
      Physical Surface(103) = {44};  // PML bot
      Physical Line(1001)   = {12,113,114,115,9}; // bloch x left
      Physical Line(1002)   = {14,116,117,118,11}; // bloch x right
      Physical Line(1003)   = {10}; // top bound
      Physical Line(1004)   = {13}; // bot bound
      Physical Line(1005)   = {38,39,40,41,42,43,44,45,46,47,48,35}; // bound for lag
      Physical Point(10000) = {1};   // Printpoint
    EndIf
    
    If (flag_o2g==1)
      Mesh.ElementOrder = 2;
    Else
      Mesh.ElementOrder = 1;
    EndIf