/////////////////////
/// GMSH GEOMETRY ///
/////////////////////

/// FILE: permanentMagnet.geo
/// AUTHOR: Frederic Trillaud <ftrillaudp@gmail.com>
/// UNIVERSITY: UNAM
/// INSTITUTE: IINGEN
/// DATE (english): 02/23/2013
/// UNITS: SI otherwise indicated

/// INTRODUCTION: automatic generation of the geometry to be meshed. Case of a brushless motor with permanent magnets. The user can use saliency id required.

/// The mesh is generated from the graphical interface or the command prompt, taken parameters from the geometry file (see manual):
// > gmsh filename.geo -2 -o fileName.msh
// -2: 2D mesh
// -o: output to "fileName.msh" file


//////////////////////
/// INCLUDED FILES ///
//////////////////////

/// Include parameter file to define parameters of the geometry ///
Include "./permanentMagnet.par";


////////////////
/// MESHING ///
///////////////


/// To copy all the properties of the geometries during transformation for meshing ///
//Geometry.CopyMeshingMethod = 1;

/// list of target element size (see parameter file) ///
n[] = {n0, n1, n2, n3, n4, n5, n6, n7, n8, n9, n10, n11, n12, n13, n14, n15, n16, n17, n18, n19, n20, n21, n22, n23, n24, n25, n26, n27, n28, n29, n30, n31, n32};


//////////////////
/// GEOMETRIES ///
//////////////////

/// Origin ///
x000 = 0; y000 = 0; z000 = 0;

/// Coordinates ///
/// Rotor
x1 = radiusRotor;                      y1 = 0.;
x2 = radiusRotor*Cos[anglePM/2.];      y2 = radiusRotor*Sin[anglePM/2.];
x3 = 0.;                               y3 = radiusRotor;
x4 = (-1)*radiusRotor*Cos[anglePM/2.]; y4 = radiusRotor*Sin[anglePM/2.];
x5 = (-1)*radiusRotor;                 y5 = 0.;
x6 = (-1)*radiusRotor*Cos[anglePM/2.]; y6 = (-1)*radiusRotor*Sin[anglePM/2.];
x7 = 0.;                               y7 = (-1)*radiusRotor;
x8 = radiusRotor*Cos[anglePM/2.];      y8 = (-1)*radiusRotor*Sin[anglePM/2.];
/// Permanent magnet 1:
x9 = radiusRotor*Cos[anglePM/2.];                y9 = (-1)*radiusRotor*Sin[anglePM/2.];
x10 = (radiusRotor+thicknessPM)*Cos[anglePM/2.]; y10 = (-1)*(radiusRotor+thicknessPM)*Sin[anglePM/2.];
x11 = radiusRotor+thicknessPM;                   y11 = 0.;
x12 = (radiusRotor+thicknessPM)*Cos[anglePM/2.]; y12 = (radiusRotor+thicknessPM)*Sin[anglePM/2.];
x13 = radiusRotor*Cos[anglePM/2.];               y13 = radiusRotor*Sin[anglePM/2.];
x14 = radiusRotor;                               y14 = 0.;
/// Permanent magnet 2:
x15 = (-1)*radiusRotor;                               y15 = 0.;
x16 = (-1)*radiusRotor*Cos[anglePM/2.];               y16 = radiusRotor*Sin[anglePM/2.];
x17 = (-1)*(radiusRotor+thicknessPM)*Cos[anglePM/2.]; y17 = (radiusRotor+thicknessPM)*Sin[anglePM/2.];
x18 = (-1)*(radiusRotor+thicknessPM);                 y18 = 0.;
x19 = (-1)*(radiusRotor+thicknessPM)*Cos[anglePM/2.]; y19 = (-1)*(radiusRotor+thicknessPM)*Sin[anglePM/2.];
x20 = (-1)*radiusRotor*Cos[anglePM/2.];               y20 = (-1)*radiusRotor*Sin[anglePM/2.];
// Stator:
x21 = innerRadiusStator;      y21 = 0.;
x22 = 0.;                     y22 = innerRadiusStator;
x23 = (-1)*innerRadiusStator; y23 = 0.;
x24 = 0.;                     y24 = (-1)*innerRadiusStator;
x25 = outerRadiusStator;      y25 = 0.;
x26 = 0.;                     y26 = outerRadiusStator;
x27 = (-1)*outerRadiusStator; y27 = 0.;
x28 = 0.;                     y28 = (-1)*outerRadiusStator;
// Infinite boundary space:
x29 = radiusInfiniteBoundarySpace;      y29 = 0.;
x30 = 0.;                               y30 = radiusInfiniteBoundarySpace;
x31 = (-1)*radiusInfiniteBoundarySpace; y31 = 0.;
x32 = 0.;                               y32 = (-1)*radiusInfiniteBoundarySpace;
/// List of coordinates
x[] = {x000, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x26, x27, x28, x29, x30, x31, x32};
y[] = {y000, y1, y2, y3, y4, y5, y6, y7, y8, y9, y10, y11, y12, y13, y14, y15, y16, y17, y18, y19, y20, y21, y22, y23, y24, y25, y26, y27, y28, y29, y30, y31, y32};
z = 0.; // 2D geometry

/// Points ///
For i In {0:32}
 Point(i) = {x[i], y[i], z, n[i]};
EndFor

/// Rotor ///
For i In {1:7}
 Circle(i) = {i, 0, i+1};
EndFor
Circle(8) = {8, 0, 1};
For i In {1:8}
 listRotor[] += {i};
EndFor
Line Loop(1) = listRotor[];
Plane Surface(1) = {1};

/// Permanent magnets ///
/// PM1:
Line(9) = {9, 10};
Circle(10) = {10, 0, 11};
Circle(11) = {11, 0, 12};
Line(12) = {12, 13};
Circle(13) = {13, 0, 14};
Circle(14) = {14, 0, 9};
For i In {9:14}
 listPM1[] += {i};
EndFor
Line Loop(2) = listPM1[];
Plane Surface(2) = {2};
/// PM2:
Circle(15) = {15, 0, 16};
Line(16) = {16, 17};
Circle(17) = {17, 0, 18};
Circle(18) = {18, 0, 19};
Line(19) = {19, 20};
Circle(20) = {20, 0, 15};
For i In {15:20}
 listPM2[] += {i};
EndFor
Line Loop(3) = listPM2[];
Plane Surface(3) = {3};

/// Air gap ///
Circle(21) = {11, 0, 12};
Line(22) = {12, 13};
Circle(23) = {13, 0, 3};
Circle(24) = {3, 0, 16};
Line(25) = {16, 17};
Circle(26) = {17, 0, 18};
Circle(27) = {18, 0, 19};
Line(28) = {19, 20};
Circle(29) = {20, 0, 7};
Circle(30) = {7, 0, 8};
Line(31) = {8, 10};
Circle(32) = {10, 0, 11};
For i In {21:32}
 listInnerAirGap[] += {i};
EndFor
Line Loop(4) = listInnerAirGap[];
For i In {1:3}
Circle(32+i) = {20+i, 0, 21+i};
EndFor
Circle(36) = {24, 0, 21};
For i In {33:36}
 listOuterAirGap[] += {i};
EndFor
Line Loop(5) = listOuterAirGap[];
Plane Surface(4) = {5, 4};

/// Stator ///
For i In {1:3}
 Circle(36+i) = {24+i, 0, 25+i};
EndFor
Circle(40) = {28, 0, 25};
For i In {37:40}
 listStator[] += {i};
EndFor
Line Loop(6) = listStator[];
Plane Surface(5) = {6, 5};

/// Infinite boundary space ///
For i In {1:3}
 Circle(40+i) = {28+i, 0, 29+i};
EndFor
Circle(44) = {32, 0, 29};
For i In {41:44}
 listInfinite[] += {i};
EndFor
Line Loop(7) = listInfinite[];
Plane Surface(6) = {7, 6};

/// Renove duplicate nodes of the mesh ///
Coherence Mesh;

///////////////////////
/// PHYSICAL GROUPS ///
///////////////////////
/// NB: define groups of similar properties or to be considered as one by the solver

/// Rotor ///
Physical Surface(rotorID) = {1};
/// Permanent magnets ///
/// PM1:
Physical Surface(PMID) = {2};
/// PM2:
Physical Surface(PMID+1) = {3};
/// Air gap ///
Physical Surface(airGapID) = {4};
/// Stator ///
Physical Surface(statorID) = {5};
/// Infinite boundary space ///
Physical Surface(infiniteBoundarySpaceID) = {6};

/// Boundary line ///
Physical Line(boundaryLineID) = listInfinite[];
