Skip to content
Snippets Groups Projects
Commit 17d694ca authored by Goffin Pierre-Yves's avatar Goffin Pierre-Yves
Browse files

Coupled .geo

parent 3577a2e5
No related branches found
No related tags found
2 merge requests!15Py bem,!14Py bem
scale = 1e-6;
bt = 1; // Beam tickness
bh = 1; // Beam height (position)
eL = 3; // Electrode base length
el = 2.5; // Electrode tip length
et = 1; // Electrode tickness
es = 1.5; // Electrode spacing
W = 3*es + 3*eL; // Width of the space
H = bh + 3*(bt+et); // Height of the space
nEs = 9;
nEL = 12;
nBt = 6;
nEt = 3;
nBh = 3;
bt = bt * scale;
bh = bh * scale;
eL = eL * scale;
el = el * scale;
et = et * scale;
es = es * scale;
W = W * scale;
H = H * scale;
e = (eL - el)/2; // Electrode constant
/*
* Points
*/
// Space
Point(1) = {0, 0, 0, 1};
Point(2) = {W, 0, 0, 1};
Point(3) = {W, H, 0, 1};
Point(4) = {0, H, 0, 1};
// Beam
Point(5) = {0, bh+bt, 0, 1};
Point(6) = {es, bh+bt, 0, 1};
Point(7) = {es+eL, bh+bt, 0, 1};
Point(8) = {2*es+eL, bh+bt, 0, 1};
Point(9) = {2*es+2*eL, bh+bt, 0, 1};
Point(10) = {3*es+2*eL, bh+bt, 0, 1};
Point(11) = {3*es+2*eL, bh, 0, 1};
Point(12) = {2*es+2*eL, bh, 0, 1};
Point(13) = {2*es+eL, bh, 0, 1};
Point(14) = {es+eL, bh, 0, 1};
Point(15) = {es, bh, 0, 1};
Point(16) = {0, bh, 0, 1};
// Electrodes
Point(17) = {es+e, bh+bt+et, 0, 1};
Point(18) = {es+e+el, bh+bt+et, 0, 1};
Point(19) = {2*es+eL+e, bh+bt+et, 0, 1};
Point(20) = {2*es+eL+e+el, bh+bt+et, 0, 1};
/*
* Lines
*/
// Space
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 5};
Line(5) = {16, 1};
// Beam
Line(6) = {5, 6};
Line(7) = {6, 7};
Line(8) = {7, 8};
Line(9) = {8, 9};
Line(10) = {9, 10};
Line(11) = {10, 11};
Line(12) = {11, 12};
Line(13) = {12, 13};
Line(14) = {13, 14};
Line(15) = {14, 15};
Line(16) = {15, 16};
Line(17) = {16, 5};
Line(18) = {15, 6};
Line(19) = {14, 7};
Line(20) = {13, 8};
Line(21) = {12, 9};
// Electrode
Line(22) = {6, 17};
Line(23) = {17, 18};
Line(24) = {18, 7};
Line(25) = {8, 19};
Line(26) = {19, 20};
Line(27) = {20, 9};
/*
* Transfinite Curve
*/
Transfinite Curve{6, 8, 10, 12, 14, 16} = nEs + 1 Using Progression 1;
Transfinite Curve{7, 9, 13, 15, 23, 26} = nEL + 1 Using Progression 1;
Transfinite Curve{17, 18, 19, 20, 21, 11} = nBt + 1 Using Progression 1;
Transfinite Curve{22, 24, 25, 27} = nEt + 1 Using Progression 1;
Transfinite Curve{1, 3} = 3*nEs + 3*nEL + 1 Using Progression 1;
Transfinite Curve(2) = 4 * (nBt + nEt) + 1 Using Progression 1;
Transfinite Curve(4) = 3*nEt + 2*nBt + 1 Using Progression 1;
Transfinite Curve(5) = nBh + 1 Using Progression 1;
/*
* Curve Loop
*/
Curve Loop(1) = {6, -18, 16, 17};
Curve Loop(2) = {7, -19, 15, 18};
Curve Loop(3) = {8, -20, 14, 19};
Curve Loop(4) = {9, -21, 13, 20};
Curve Loop(5) = {10, 11, 12, 21};
Curve Loop(6) = {23, 24, -7, 22};
Curve Loop(7) = {26, 27, -9, 25};
Curve Loop(8) = {1, 2, 3, 4, 6, 22, 23, 24, 8, 25, 26, 27, 10, 11, 12, 13, 14, 15, 16, 5};
/*
* Plane Surface
*/
Plane Surface(1) = {1};
Plane Surface(2) = {2};
Plane Surface(3) = {3};
Plane Surface(4) = {4};
Plane Surface(5) = {5};
Plane Surface(6) = {6};
Plane Surface(7) = {7};
Plane Surface(8) = {8};
/*
* Surface Manipulation
*/
Transfinite Surface{1, 2, 3, 4, 5, 6, 7};
Recombine Surface{1, 2, 3, 4, 5, 6, 7};
/*
* Physical Curve
*/
Physical Curve("clamp", 1) = {17};
Physical Curve("outside", 2) = {1, 2, 3, 4, 5};
Physical Curve("left_electrode", 3) = {24};
Physical Curve("right_electrode", 4) = {25};
Physical Curve("beam", 5) = {6, 22, 23, 8, 26, 10, 11, 12, 13, 14, 15, 16};
Physical Curve("BEM_FEM_boundary", 6) = {6, 22, 23, 24, 8, 26, 26, 27, 10, 11, 12, 13, 14, 15, 16};
/*
* Physical Surface
*/
Physical Surface("FEM_domain", 1) = {1, 2, 3, 4, 5, 6, 7};
Physical Surface("BEM_domain_1", 2) = {8};
/*
* Parameters
*/
SetNumber("Boundary Conditions/clamp/ux", 0.0);
SetNumber("Boundary Conditions/clamp/uy", 0.0);
SetNumber("Materials/FEM_domain/Young", 3e7);
SetNumber("Materials/FEM_domain/Poisson", 0.3);
SetNumber("Materials/FEM_domain/rho",7800);
SetNumber("Volumic Forces/FEM_domain/bx", 0.0);
SetNumber("Volumic Forces/FEM_domain/by", 0.0);
SetNumber("Boundary Conditions/outside/BEM_domain_1/neumann", 0);
SetNumber("Boundary Conditions/beam/BEM_domain_1/neumann", 0);
SetNumber("Boundary Conditions/left_electrode/BEM_domain_1/dirichlet", 200);
SetNumber("Boundary Conditions/right_electrode/BEM_domain_1/dirichlet", 0);
SetNumber("Materials/BEM_domain_1/Epsilon", 8.8541878128e-12); // dielectric permittivity
SetNumber("Non_linear_solver", 1);
\ No newline at end of file
scale = 1e-6;
bt = 0.2; // Beam tickness
bh = 1; // Beam height (position)
eL = 3; // Electrode base length
el = 1.5; // Electrode tip length
et = 1; // Electrode tickness
es = 1.5; // Electrode spacing
nbElectrodes = 10;
W = nbElectrodes*(es + eL) + 2*es; // Width of the space
H = W; // Height of the space
nEs = 12;
nEL = 15;
nBt = 8;
nEt = 3;
nBh = 3;
bt = bt * scale;
bh = bh * scale;
eL = eL * scale;
el = el * scale;
et = et * scale;
es = es * scale;
W = W * scale;
H = H * scale;
e = (eL - el)/2; // Electrode constant
/*
* Points
*/
For i In {0:nbElectrodes-1}
Point(6*i + 1) = {(es+eL)*i, bh + bt, 0, 1};
Point(6*i + 2) = {(es+eL)*i + es, bh + bt, 0, 1};
Point(6*i + 3) = {(es+eL)*i + es+e, bh + bt + et, 0, 1};
Point(6*i + 4) = {(es+eL)*i + es+e+el, bh + bt + et, 0, 1};
Point(6*i + 5) = {(es+eL)*i + es, bh, 0, 1};
Point(6*i + 6) = {(es+eL)*i, bh, 0, 1};
EndFor
Point(6*nbElectrodes + 1) = {nbElectrodes*(es+eL), bh + bt, 0, 1};
Point(6*nbElectrodes + 2) = {nbElectrodes*(es+eL) + es, bh + bt, 0, 1};
Point(6*nbElectrodes + 3) = {nbElectrodes*(es+eL) + es, bh, 0, 1};
Point(6*nbElectrodes + 6) = {nbElectrodes*(es+eL), bh, 0, 1};
Point(6*nbElectrodes + 4) = {0, 0, 0, 1};
Point(6*nbElectrodes + 5) = {W, 0, 0, 1};
Point(6*nbElectrodes + 7) = {W, H, 0, 1};
Point(6*nbElectrodes + 8) = {0, H, 0, 1};
For i In {0:nbElectrodes-1}
Point(6*nbElectrodes + 8 + i + 1) = {i*(es+eL) + es + e + el/2, bh + bt + et - ((el*e)/(2*et)), 0, 1};
EndFor
/*
* Lines
*/
Line(9*nbElectrodes + 1) = {6*nbElectrodes + 4, 6*nbElectrodes + 5};
Line(9*nbElectrodes + 2) = {6*nbElectrodes + 5, 6*nbElectrodes + 7};
Line(9*nbElectrodes + 3) = {6*nbElectrodes + 7, 6*nbElectrodes + 8};
Line(9*nbElectrodes + 4) = {6*nbElectrodes + 8, 1};
Line(9*nbElectrodes + 5) = {6, 6*nbElectrodes + 4};
Line(9*nbElectrodes + 9) = {6*nbElectrodes + 1, 6*nbElectrodes + 2};
Line(9*nbElectrodes + 6) = {6*nbElectrodes + 2, 6*nbElectrodes + 3};
Line(9*nbElectrodes + 7) = {6*nbElectrodes + 3, 6*nbElectrodes + 6};
Line(9*nbElectrodes + 8) = {6*nbElectrodes + 6, 6*nbElectrodes + 1};
For i In {0:nbElectrodes-1}
Line(9*i + 1) = {6*i + 1, 6*i + 2};
Line(9*i + 2) = {6*i + 2, 6*i + 3};
// Line(9*i + 3) = {6*i + 3, 6*i + 4};
Circle(9*i + 3) = {6*i + 3, 6*nbElectrodes + 8 + i + 1, 6*i + 4};
Line(9*i + 4) = {6*i + 4, 6*(i+1) + 1};
Line(9*i + 5) = {6*(i+1) + 1, 6*i + 2};
Line(9*i + 6) = {6*i + 5, 6*i + 2};
Line(9*i + 7) = {6*i + 5, 6*i + 6};
Line(9*i + 8) = {6*i + 6, 6*i + 1};
Line(9*i + 9) = {6*(i+1) + 6, 6*i + 5};
EndFor
/*
* Transfinite Curve
*/
For i In {0:nbElectrodes-1}
Transfinite Curve{9*i + 1, 9*i + 7} = nEs + 1 Using Progression 1;
Transfinite Curve{9*i + 2, 9*i + 4} = nEt + 1 Using Progression 1;
Transfinite Curve{9*i + 3, 9*i + 5, 9*i + 9} = nEL + 1 Using Progression 1;
Transfinite Curve{9*i + 8, 9*i + 6} = nBt + 1 Using Progression 1;
EndFor
Transfinite Curve{9*nbElectrodes + 9, 9*nbElectrodes + 7} = nEs + 1 Using Progression 1;
Transfinite Curve{9*nbElectrodes + 8, 9*nbElectrodes + 6} = nBt + 1 Using Progression 1;
Transfinite Curve{9*nbElectrodes + 1, 9*nbElectrodes + 3} = nbElectrodes*(nEs + nEL) + 2*nEs + 1 Using Progression 1;
Transfinite Curve{9*nbElectrodes + 2} = nbElectrodes*(nEs + nEL) + 2*nEs + 1 Using Progression 1;
Transfinite Curve{9*nbElectrodes + 4} = nbElectrodes*(nEs + nEL) + 2*nEs - nBt - nBh + 1 Using Progression 1;
Transfinite Curve{9*nbElectrodes + 5} = nBh + 1 Using Progression 1;
/*
* Curve Loop
*/
For i In {0:nbElectrodes-1}
Curve Loop(3*i + 1) = {9*i + 1, -(9*i + 6), 9*i + 7, 9*i + 8};
Curve Loop(3*i + 2) = {9*i + 3, 9*i + 4, 9*i + 5, 9*i + 2};
Curve Loop(3*i + 3) = {-(9*i + 5), -(9*(i+1) + 8), 9*i + 9, 9*i + 6};
EndFor
Curve Loop(3*nbElectrodes + 1) = {9*nbElectrodes + 9, 9*nbElectrodes + 6, 9*nbElectrodes + 7, 9*nbElectrodes + 8};
// How to fill a vector in a loop ?
bemBoundary = {1:6*nbElectrodes + 8};
bemBoundaryWithoutE = {1:4*nbElectrodes + 8 + 2};
bemFemBoundary = {1:6*nbElectrodes + 3};
leftElectrodesFull = {1:nbElectrodes};
rightElectrodesFull = {1:nbElectrodes};
For i In {0:nbElectrodes-1}
bemBoundary[4*i+1 - 1] = 9*i + 1;
bemBoundary[4*i+2 - 1] = 9*i + 2;
bemBoundary[4*i+3 - 1] = 9*i + 3;
bemBoundary[4*i+4 - 1] = 9*i + 4;
bemBoundaryWithoutE[2*i+1 - 1] = 9*i + 1;
bemBoundaryWithoutE[2*i+2 - 1] = 9*i + 3;
bemFemBoundary[4*i+1 - 1] = 9*i + 1;
bemFemBoundary[4*i+2 - 1] = 9*i + 2;
bemFemBoundary[4*i+3 - 1] = 9*i + 3;
bemFemBoundary[4*i+4 - 1] = 9*i + 4;
leftElectrodesFull[i] = 9*i + 2;
rightElectrodesFull[i] = 9*i + 4;
EndFor
bemBoundary[4*nbElectrodes + 1 - 1] = 9*nbElectrodes + 9;
bemBoundary[4*nbElectrodes + 2 - 1] = 9*nbElectrodes + 6;
bemBoundary[4*nbElectrodes + 3 - 1] = 9*nbElectrodes + 7;
bemBoundaryWithoutE[2*nbElectrodes + 1 - 1] = 9*nbElectrodes + 9;
bemBoundaryWithoutE[2*nbElectrodes + 2 - 1] = 9*nbElectrodes + 6;
bemBoundaryWithoutE[2*nbElectrodes + 3 - 1] = 9*nbElectrodes + 7;
bemFemBoundary[4*nbElectrodes + 1 - 1] = 9*nbElectrodes + 9;
bemFemBoundary[4*nbElectrodes + 2 - 1] = 9*nbElectrodes + 6;
bemFemBoundary[4*nbElectrodes + 3 - 1] = 9*nbElectrodes + 7;
For i In {0:nbElectrodes-1}
bemBoundary[4*nbElectrodes + 3 + 2*i+1 - 1] = 9*i + 9;
bemBoundary[4*nbElectrodes + 3 + 2*i+2 - 1] = 9*i + 7;
bemBoundaryWithoutE[2*nbElectrodes + 3 + 2*i+1 - 1] = 9*i + 9;
bemBoundaryWithoutE[2*nbElectrodes + 3 + 2*i+2 - 1] = 9*i + 7;
bemFemBoundary[4*nbElectrodes + 3 + 2*i+1 - 1] = 9*i + 9;
bemFemBoundary[4*nbElectrodes + 3 + 2*i+2 - 1] = 9*i + 7;
EndFor
bemBoundary[6*nbElectrodes + 4 - 1] = 9*nbElectrodes + 5;
bemBoundary[6*nbElectrodes + 5 - 1] = 9*nbElectrodes + 1;
bemBoundary[6*nbElectrodes + 6 - 1] = 9*nbElectrodes + 2;
bemBoundary[6*nbElectrodes + 7 - 1] = 9*nbElectrodes + 3;
bemBoundary[6*nbElectrodes + 8 - 1] = 9*nbElectrodes + 4;
bemBoundaryWithoutE[4*nbElectrodes + 4 - 1] = 9*nbElectrodes + 5;
bemBoundaryWithoutE[4*nbElectrodes + 5 - 1] = 9*nbElectrodes + 1;
bemBoundaryWithoutE[4*nbElectrodes + 6 - 1] = 9*nbElectrodes + 2;
bemBoundaryWithoutE[4*nbElectrodes + 7 - 1] = 9*nbElectrodes + 3;
bemBoundaryWithoutE[4*nbElectrodes + 8 - 1] = 9*nbElectrodes + 4;
leftElectrodes = {1:nbElectrodes-1};
rightElectrodes = {1:nbElectrodes-1};
For i In {0:nbElectrodes-2}
leftElectrodes[i] = leftElectrodesFull[i+1];
rightElectrodes[i] = rightElectrodesFull[i];
EndFor
bemBoundaryWithoutE[4*nbElectrodes + 9 - 1] = 2;
bemBoundaryWithoutE[4*nbElectrodes + 10 - 1] = 9*nbElectrodes + 4;
Curve Loop(3*nbElectrodes + 2) = bemBoundary[];
/*
* Plane Surface
*/
For i In {0:nbElectrodes-1}
Plane Surface(3*i + 1) = {3*i + 1};
Plane Surface(3*i + 2) = {3*i + 2};
Plane Surface(3*i + 3) = {3*i + 3};
EndFor
Plane Surface(3*nbElectrodes + 1) = {3*nbElectrodes + 1};
Plane Surface(3*nbElectrodes + 2) = {3*nbElectrodes + 2};
/*
* Surface operations
*/
Recombine Surface{1:3*nbElectrodes + 1};
Transfinite Surface{1:3*nbElectrodes + 1};
/*
* Physical Curve
*/
Physical Curve("clamp", 1) = {8};
Physical Curve("BEM_boundary_without_e", 2) = bemBoundaryWithoutE[];
Physical Curve("BEM_FEM_boundary", 3) = bemFemBoundary[];
Physical Curve("left_electrodes", 4) = leftElectrodes[];
Physical Curve("right_electrodes", 5) = rightElectrodes[];
Physical Surface("FEM_domain") = {1:3*nbElectrodes + 1};
Physical Surface("BEM_domain_1") = {3*nbElectrodes + 2};
/*
* Parameters
*/
SetNumber("Boundary Conditions/clamp/ux", 0.0);
SetNumber("Boundary Conditions/clamp/uy", 0.0);
// SetNumber("Materials/FEM_domain/Young", 3e7);
SetNumber("Materials/FEM_domain/Young", 1e7);
SetNumber("Materials/FEM_domain/Poisson", 0.3);
SetNumber("Materials/FEM_domain/rho",7800);
SetNumber("Volumic Forces/FEM_domain/bx", 0.0);
SetNumber("Volumic Forces/FEM_domain/by", 0.0);
SetNumber("Boundary Conditions/BEM_boundary_without_e/BEM_domain_1/neumann", 0);
SetNumber("Boundary Conditions/left_electrodes/BEM_domain_1/dirichlet", 8);
SetNumber("Boundary Conditions/right_electrodes/BEM_domain_1/dirichlet", 0);
SetNumber("Materials/BEM_domain_1/Epsilon", 8.8541878128e-12); // dielectric permittivity
SetNumber("Non_linear_solver", 1);
\ No newline at end of file
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment