//**** Definition of points

p1 = newp;
lref= 0.100;
Point(p1+1) ={   0.5000,   0.0000,  -0.1500, lref};
Point(p1+2) ={   0.5000,   0.1000,  -0.1500, lref};
Point(p1+3) ={   0.3500,   0.1000,  -0.1500, lref};
Point(p1+4) ={   0.5000,  -0.0000,  -0.0500, lref};
Point(p1+5) ={   0.3500,  -0.0000,  -0.0500, lref};
Point(p1+6) ={   0.5000,  -0.1000,  -0.1500, lref};
Point(p1+7) ={   0.3500,  -0.1000,  -0.1500, lref};
Point(p1+8) ={   0.5000,   0.0000,  -0.2500, lref};
Point(p1+9) ={   0.3500,   0.0000,  -0.2500, lref};
Point(p1+10) ={   0.3500,   0.0000,  -0.1500, lref};

lref= 0.010;
Point(p1+11) ={   0.5000,   0.0210,  -0.1500, lref};
Point(p1+12) ={   0.5000,   0.0310,  -0.1500, lref};
Point(p1+13) ={   0.5000,  -0.0000,  -0.1290, lref};
Point(p1+14) ={   0.5000,  -0.0000,  -0.1190, lref};
Point(p1+15) ={   0.5000,  -0.0210,  -0.1500, lref};
Point(p1+16) ={   0.5000,  -0.0310,  -0.1500, lref};
Point(p1+17) ={   0.5000,   0.0000,  -0.1710, lref};
Point(p1+18) ={   0.5000,   0.0000,  -0.1810, lref};

lref= 0.521;
Point(p1+20) ={   1.6654,   0.0000,  -0.1500, lref};
Point(p1+21) ={   1.5500,   0.5206,  -0.1500, lref};
Point(p1+22) ={   1.5500,  -0.0000,   0.3706, lref};
Point(p1+23) ={   1.5500,  -0.5206,  -0.1500, lref};
Point(p1+24) ={   1.5500,   0.0000,  -0.6706, lref};

Point(p1+25) ={   0.4335,   0.0000,  -0.1500, lref};

Point(p1+26) ={   0.3500,   1.0000,  -0.1500, lref};
Point(p1+27) ={   1.7000,   1.0000,  -0.1500, lref};
Point(p1+28) ={   0.3500,  -0.0000,   0.8500, lref};
Point(p1+29) ={   1.7000,  -0.0000,   0.8500, lref};
Point(p1+30) ={   0.3500,  -1.0000,  -0.1500, lref};
Point(p1+31) ={   1.7000,  -1.0000,  -0.1500, lref};
Point(p1+32) ={   0.3500,   0.0000,  -1.1500, lref};
Point(p1+33) ={   1.7000,   0.0000,  -1.1500, lref};
Point(p1+34) ={   1.7000,   0.0000,  -0.1500, lref};

//**** points to embed in auxiliary surfaces

Point(p1+35) ={   0.9276,   0.0180,  -0.1500,   0.0283};
Point(p1+36) ={   0.9267,   0.0359,  -0.1500,   0.0283};
Point(p1+37) ={   0.9250,   0.0538,  -0.1500,   0.0283};
Point(p1+38) ={   0.9227,   0.0717,  -0.1500,   0.0283};
Point(p1+39) ={   0.9198,   0.0894,  -0.1500,   0.0283};
Point(p1+40) ={   0.9162,   0.1070,  -0.1500,   0.0283};
Point(p1+41) ={   0.9120,   0.1245,  -0.1500,   0.0283};
Point(p1+42) ={   0.9072,   0.1418,  -0.1500,   0.0283};
Point(p1+43) ={   0.9017,   0.1589,  -0.1500,   0.0283};
Point(p1+44) ={   0.8956,   0.1758,  -0.1500,   0.0283};
Point(p1+45) ={   0.8889,   0.1925,  -0.1500,   0.0283};
Point(p1+46) ={   1.2863,   0.0310,  -0.1500,   0.0565};
Point(p1+47) ={   1.2846,   0.0620,  -0.1500,   0.0565};
Point(p1+48) ={   1.2818,   0.0929,  -0.1500,   0.0565};
Point(p1+49) ={   1.2779,   0.1237,  -0.1500,   0.0565};
Point(p1+50) ={   1.2728,   0.1543,  -0.1500,   0.0565};
Point(p1+51) ={   1.2667,   0.1847,  -0.1500,   0.0565};
Point(p1+52) ={   1.2594,   0.2149,  -0.1500,   0.0565};
Point(p1+53) ={   1.2510,   0.2447,  -0.1500,   0.0565};
Point(p1+54) ={   1.2416,   0.2743,  -0.1500,   0.0565};
Point(p1+55) ={   1.2311,   0.3035,  -0.1500,   0.0565};
Point(p1+56) ={   1.2195,   0.3323,  -0.1500,   0.0565};
Point(p1+57) ={   0.9276,  -0.0000,  -0.1320,   0.0283};
Point(p1+58) ={   0.9267,  -0.0000,  -0.1141,   0.0283};
Point(p1+59) ={   0.9250,  -0.0000,  -0.0962,   0.0283};
Point(p1+60) ={   0.9227,  -0.0000,  -0.0783,   0.0283};
Point(p1+61) ={   0.9198,  -0.0000,  -0.0606,   0.0283};
Point(p1+62) ={   0.9162,  -0.0000,  -0.0430,   0.0283};
Point(p1+63) ={   0.9120,  -0.0000,  -0.0255,   0.0283};
Point(p1+64) ={   0.9072,  -0.0000,  -0.0082,   0.0283};
Point(p1+65) ={   0.9017,  -0.0000,   0.0089,   0.0283};
Point(p1+66) ={   0.8956,  -0.0000,   0.0258,   0.0283};
Point(p1+67) ={   0.8889,  -0.0000,   0.0425,   0.0283};
Point(p1+68) ={   1.2863,  -0.0000,  -0.1190,   0.0565};
Point(p1+69) ={   1.2846,  -0.0000,  -0.0880,   0.0565};
Point(p1+70) ={   1.2818,  -0.0000,  -0.0571,   0.0565};
Point(p1+71) ={   1.2779,  -0.0000,  -0.0263,   0.0565};
Point(p1+72) ={   1.2728,  -0.0000,   0.0043,   0.0565};
Point(p1+73) ={   1.2667,  -0.0000,   0.0347,   0.0565};
Point(p1+74) ={   1.2594,  -0.0000,   0.0649,   0.0565};
Point(p1+75) ={   1.2510,  -0.0000,   0.0947,   0.0565};
Point(p1+76) ={   1.2416,  -0.0000,   0.1243,   0.0565};
Point(p1+77) ={   1.2311,  -0.0000,   0.1535,   0.0565};
Point(p1+78) ={   1.2195,  -0.0000,   0.1823,   0.0565};
Point(p1+79) ={   0.9276,  -0.0180,  -0.1500,   0.0283};
Point(p1+80) ={   0.9267,  -0.0359,  -0.1500,   0.0283};
Point(p1+81) ={   0.9250,  -0.0538,  -0.1500,   0.0283};
Point(p1+82) ={   0.9227,  -0.0717,  -0.1500,   0.0283};
Point(p1+83) ={   0.9198,  -0.0894,  -0.1500,   0.0283};
Point(p1+84) ={   0.9162,  -0.1070,  -0.1500,   0.0283};
Point(p1+85) ={   0.9120,  -0.1245,  -0.1500,   0.0283};
Point(p1+86) ={   0.9072,  -0.1418,  -0.1500,   0.0283};
Point(p1+87) ={   0.9017,  -0.1589,  -0.1500,   0.0283};
Point(p1+88) ={   0.8956,  -0.1758,  -0.1500,   0.0283};
Point(p1+89) ={   0.8889,  -0.1925,  -0.1500,   0.0283};
Point(p1+90) ={   1.2863,  -0.0310,  -0.1500,   0.0565};
Point(p1+91) ={   1.2846,  -0.0620,  -0.1500,   0.0565};
Point(p1+92) ={   1.2818,  -0.0929,  -0.1500,   0.0565};
Point(p1+93) ={   1.2779,  -0.1237,  -0.1500,   0.0565};
Point(p1+94) ={   1.2728,  -0.1543,  -0.1500,   0.0565};
Point(p1+95) ={   1.2667,  -0.1847,  -0.1500,   0.0565};
Point(p1+96) ={   1.2594,  -0.2149,  -0.1500,   0.0565};
Point(p1+97) ={   1.2510,  -0.2447,  -0.1500,   0.0565};
Point(p1+98) ={   1.2416,  -0.2743,  -0.1500,   0.0565};
Point(p1+99) ={   1.2311,  -0.3035,  -0.1500,   0.0565};
Point(p1+100) ={   1.2195,  -0.3323,  -0.1500,   0.0565};
Point(p1+101) ={   0.9276,   0.0000,  -0.1680,   0.0283};
Point(p1+102) ={   0.9267,   0.0000,  -0.1859,   0.0283};
Point(p1+103) ={   0.9250,   0.0000,  -0.2038,   0.0283};
Point(p1+104) ={   0.9227,   0.0000,  -0.2217,   0.0283};
Point(p1+105) ={   0.9198,   0.0000,  -0.2394,   0.0283};
Point(p1+106) ={   0.9162,   0.0000,  -0.2570,   0.0283};
Point(p1+107) ={   0.9120,   0.0000,  -0.2745,   0.0283};
Point(p1+108) ={   0.9072,   0.0000,  -0.2918,   0.0283};
Point(p1+109) ={   0.9017,   0.0000,  -0.3089,   0.0283};
Point(p1+110) ={   0.8956,   0.0000,  -0.3258,   0.0283};
Point(p1+111) ={   0.8889,   0.0000,  -0.3425,   0.0283};
Point(p1+112) ={   1.2863,   0.0000,  -0.1810,   0.0565};
Point(p1+113) ={   1.2846,   0.0000,  -0.2120,   0.0565};
Point(p1+114) ={   1.2818,   0.0000,  -0.2429,   0.0565};
Point(p1+115) ={   1.2779,   0.0000,  -0.2737,   0.0565};
Point(p1+116) ={   1.2728,   0.0000,  -0.3043,   0.0565};
Point(p1+117) ={   1.2667,   0.0000,  -0.3347,   0.0565};
Point(p1+118) ={   1.2594,   0.0000,  -0.3649,   0.0565};
Point(p1+119) ={   1.2510,   0.0000,  -0.3947,   0.0565};
Point(p1+120) ={   1.2416,   0.0000,  -0.4243,   0.0565};
Point(p1+121) ={   1.2311,   0.0000,  -0.4535,   0.0565};
Point(p1+122) ={   1.2195,   0.0000,  -0.4823,   0.0565};

//**** Definition of lines


line1 = newreg;


// archi thruster box

Circle(line1+1) ={p1+  3,p1+ 10,p1+  5};
Circle(line1+2) ={p1+  5,p1+ 10,p1+  7};
Circle(line1+3) ={p1+  7,p1+ 10,p1+  9};
Circle(line1+4) ={p1+  9,p1+ 10,p1+  3};

Circle(line1+5) ={p1+  2,p1+  1,p1+  4};
Circle(line1+6) ={p1+  4,p1+  1,p1+  6};
Circle(line1+7) ={p1+  6,p1+  1,p1+  8};
Circle(line1+8) ={p1+  8,p1+  1,p1+  2};

// archi thruster exit plane

// internal radius

Circle(line1+9) ={p1+ 11,p1+  1,p1+ 13};
Circle(line1+10) ={p1+ 13,p1+  1,p1+ 15};
Circle(line1+11) ={p1+ 15,p1+  1,p1+ 17};
Circle(line1+12) ={p1+ 17,p1+  1,p1+ 11};

// external radius

Circle(line1+13) ={p1+ 12,p1+  1,p1+ 14};
Circle(line1+14) ={p1+ 14,p1+  1,p1+ 16};
Circle(line1+15) ={p1+ 16,p1+  1,p1+ 18};
Circle(line1+16) ={p1+ 18,p1+  1,p1+ 12};

// radial lines thruster exit plane

Line(line1+17) ={p1+  1,p1+ 11};
Line(line1+18) ={p1+ 11,p1+ 12};
Line(line1+19) ={p1+  1,p1+ 13};
Line(line1+20) ={p1+ 13,p1+ 14};
Line(line1+21) ={p1+  1,p1+ 15};
Line(line1+22) ={p1+ 15,p1+ 16};
Line(line1+23) ={p1+  1,p1+ 17};
Line(line1+24) ={p1+ 17,p1+ 18};

// radial lines aux surfaces

Circle(line1+25) ={p1+ 20,p1+ 25,p1+ 21};
Circle(line1+26) ={p1+ 20,p1+ 25,p1+ 22};
Circle(line1+27) ={p1+ 20,p1+ 25,p1+ 23};
Circle(line1+28) ={p1+ 20,p1+ 25,p1+ 24};

// longitudinal lines thruster box

Line(line1+29) ={p1+  2,p1+  3};
Line(line1+30) ={p1+  4,p1+  5};
Line(line1+31) ={p1+  6,p1+  7};
Line(line1+32) ={p1+  8,p1+  9};

// thruster axis

Line(line1+33) ={p1+  1,p1+ 20};

// longitudinal lines aux surfaces

Line(line1+34) ={p1+ 12,p1+ 21};
Line(line1+35) ={p1+ 14,p1+ 22};
Line(line1+36) ={p1+ 16,p1+ 23};
Line(line1+37) ={p1+ 18,p1+ 24};

// far field lines

Circle(line1+38) ={p1+ 26,p1+ 10,p1+ 28};
Circle(line1+39) ={p1+ 28,p1+ 10,p1+ 30};
Circle(line1+40) ={p1+ 30,p1+ 10,p1+ 32};
Circle(line1+41) ={p1+ 32,p1+ 10,p1+ 26};
Circle(line1+42) ={p1+ 27,p1+ 34,p1+ 29};
Circle(line1+43) ={p1+ 29,p1+ 34,p1+ 31};
Circle(line1+44) ={p1+ 31,p1+ 34,p1+ 33};
Circle(line1+45) ={p1+ 33,p1+ 34,p1+ 27};
Line(line1+46) ={p1+ 26,p1+ 27};
Line(line1+47) ={p1+ 28,p1+ 29};
Line(line1+48) ={p1+ 30,p1+ 31};
Line(line1+49) ={p1+ 32,p1+ 33};

//**** Definition of surfaces

surf1 = newreg;

Line Loop(surf1+1) ={ line1+29, line1 +1,-line1-30,-line1 -5};
Ruled Surface(surf1+5) ={surf1+  1};
Line Loop(surf1+2) ={ line1+30, line1 +2,-line1-31,-line1 -6};
Ruled Surface(surf1+6) ={surf1+  2};
Line Loop(surf1+3) ={ line1+31, line1 +3,-line1-32,-line1 -7};
Ruled Surface(surf1+7) ={surf1+  3};
Line Loop(surf1+4) ={ line1+32, line1 +4,-line1-29,-line1 -8};
Ruled Surface(surf1+8) ={surf1+  4};

Line Loop(surf1+8) ={ line1 +5, line1 +6, line1 +7, line1 +8};
Line Loop(surf1+9) ={line1+13,line1+14,line1+15,line1+16};
Plane Surface(surf1+10) ={surf1+  8,surf1+  9};
Line Loop(surf1+11) ={line1 +9,line1+10,line1+11,line1+12};
Plane Surface(surf1+12) ={surf1+  9,surf1+ 11};
Plane Surface(surf1+13) ={surf1+ 11};

//**** Auxiliary surfaces

Line Loop(surf1+14) ={ line1+17, line1+18, line1+34,-line1-25,-line1-33};
Plane Surface(surf1+18) ={surf1+ 14};
Line Loop(surf1+15) ={ line1+19, line1+20, line1+35,-line1-26,-line1-33};
Plane Surface(surf1+19) ={surf1+ 15};
Line Loop(surf1+16) ={ line1+21, line1+22, line1+36,-line1-27,-line1-33};
Plane Surface(surf1+20) ={surf1+ 16};
Line Loop(surf1+17) ={ line1+23, line1+24, line1+37,-line1-28,-line1-33};
Plane Surface(surf1+21) ={surf1+ 17};
Line Loop(surf1+22) ={line1+38,line1+39,line1+40,line1+41};
Line Loop(surf1+23) ={line1 +1,line1 +2,line1 +3,line1 +4};
Plane Surface(surf1+24) ={surf1+ 22,surf1+ 23};
Line Loop(surf1+25) ={line1+42,line1+43,line1+44,line1+45};
Plane Surface(surf1+26) ={surf1+ 25};
Line Loop(surf1+27) ={ line1+38, line1+47,-line1-42,-line1-46};
Ruled Surface(surf1+28) ={surf1+ 27};
Line Loop(surf1+29) ={ line1+39, line1+48,-line1-43,-line1-47};
Ruled Surface(surf1+30) ={surf1+ 29};
Line Loop(surf1+31) ={ line1+40, line1+49,-line1-44,-line1-48};
Ruled Surface(surf1+32) ={surf1+ 31};
Line Loop(surf1+33) ={ line1+41, line1+46,-line1-45,-line1-49};
Ruled Surface(surf1+34) ={surf1+ 33};

//**** Embed points in surfaces

Point {p1+ 35} In Surface {surf1+ 18};
Point {p1+ 36} In Surface {surf1+ 18};
Point {p1+ 37} In Surface {surf1+ 18};
Point {p1+ 38} In Surface {surf1+ 18};
Point {p1+ 39} In Surface {surf1+ 18};
Point {p1+ 40} In Surface {surf1+ 18};
Point {p1+ 41} In Surface {surf1+ 18};
Point {p1+ 42} In Surface {surf1+ 18};
Point {p1+ 43} In Surface {surf1+ 18};
Point {p1+ 44} In Surface {surf1+ 18};
Point {p1+ 45} In Surface {surf1+ 18};
Point {p1+ 46} In Surface {surf1+ 18};
Point {p1+ 47} In Surface {surf1+ 18};
Point {p1+ 48} In Surface {surf1+ 18};
Point {p1+ 49} In Surface {surf1+ 18};
Point {p1+ 50} In Surface {surf1+ 18};
Point {p1+ 51} In Surface {surf1+ 18};
Point {p1+ 52} In Surface {surf1+ 18};
Point {p1+ 57} In Surface {surf1+ 19};
Point {p1+ 58} In Surface {surf1+ 19};
Point {p1+ 59} In Surface {surf1+ 19};
Point {p1+ 60} In Surface {surf1+ 19};
Point {p1+ 61} In Surface {surf1+ 19};
Point {p1+ 62} In Surface {surf1+ 19};
Point {p1+ 63} In Surface {surf1+ 19};
Point {p1+ 64} In Surface {surf1+ 19};
Point {p1+ 65} In Surface {surf1+ 19};
Point {p1+ 66} In Surface {surf1+ 19};
Point {p1+ 67} In Surface {surf1+ 19};
Point {p1+ 68} In Surface {surf1+ 19};
Point {p1+ 69} In Surface {surf1+ 19};
Point {p1+ 70} In Surface {surf1+ 19};
Point {p1+ 71} In Surface {surf1+ 19};
Point {p1+ 72} In Surface {surf1+ 19};
Point {p1+ 73} In Surface {surf1+ 19};
Point {p1+ 74} In Surface {surf1+ 19};
Point {p1+ 79} In Surface {surf1+ 20};
Point {p1+ 80} In Surface {surf1+ 20};
Point {p1+ 81} In Surface {surf1+ 20};
Point {p1+ 82} In Surface {surf1+ 20};
Point {p1+ 83} In Surface {surf1+ 20};
Point {p1+ 84} In Surface {surf1+ 20};
Point {p1+ 85} In Surface {surf1+ 20};
Point {p1+ 86} In Surface {surf1+ 20};
Point {p1+ 87} In Surface {surf1+ 20};
Point {p1+ 88} In Surface {surf1+ 20};
Point {p1+ 89} In Surface {surf1+ 20};
Point {p1+ 90} In Surface {surf1+ 20};
Point {p1+ 91} In Surface {surf1+ 20};
Point {p1+ 92} In Surface {surf1+ 20};
Point {p1+ 93} In Surface {surf1+ 20};
Point {p1+ 94} In Surface {surf1+ 20};
Point {p1+ 95} In Surface {surf1+ 20};
Point {p1+ 96} In Surface {surf1+ 20};
Point {p1+101} In Surface {surf1+ 21};
Point {p1+102} In Surface {surf1+ 21};
Point {p1+103} In Surface {surf1+ 21};
Point {p1+104} In Surface {surf1+ 21};
Point {p1+105} In Surface {surf1+ 21};
Point {p1+106} In Surface {surf1+ 21};
Point {p1+107} In Surface {surf1+ 21};
Point {p1+108} In Surface {surf1+ 21};
Point {p1+109} In Surface {surf1+ 21};
Point {p1+110} In Surface {surf1+ 21};
Point {p1+111} In Surface {surf1+ 21};
Point {p1+112} In Surface {surf1+ 21};
Point {p1+113} In Surface {surf1+ 21};
Point {p1+114} In Surface {surf1+ 21};
Point {p1+115} In Surface {surf1+ 21};
Point {p1+116} In Surface {surf1+ 21};
Point {p1+117} In Surface {surf1+ 21};
Point {p1+118} In Surface {surf1+ 21};

//**** Embed lines in surfaces

Line {line1+ 17} In Surface {surf1+ 13};
Line {line1+ 19} In Surface {surf1+ 13};
Line {line1+ 21} In Surface {surf1+ 13};
Line {line1+ 23} In Surface {surf1+ 13};
Line {line1+ 18} In Surface {surf1+ 12};
Line {line1+ 20} In Surface {surf1+ 12};
Line {line1+ 22} In Surface {surf1+ 12};
Line {line1+ 24} In Surface {surf1+ 12};

//**** Transfinite thruster box lines

Transfinite Line {line1+  1,line1+  2,line1+  3,line1+  4} =   4 Using Progression 1.000000;
Transfinite Line {line1+  5,line1+  6,line1+  7,line1+  8} =   4 Using Progression 1.000000;

//**** Transfinite exit channel lines

Transfinite Line {line1+  9,line1+ 10,line1+ 11,line1+ 12} =   7 Using Progression 1.000000;
Transfinite Line {line1+ 13,line1+ 14,line1+ 15,line1+ 16} =   7 Using Progression 1.000000;

//**** Transfinite longitudinal lines

prog=1.097392;
prog2=1.097659;

Transfinite Line {line1+ 33} =  35 Using Progression 1.097659;
Transfinite Line {line1+ 34,line1+ 35,line1+ 36,line1+ 37} =  35 Using Progression 1.097392;

//**** Transfinite radial lines

Transfinite Line {line1+ 17,line1+ 19,line1+ 21,line1+ 23} =   9 Using Progression 1.000000;
Transfinite Line {line1+ 18,line1+ 20,line1+ 22,line1+ 24} =   3 Using Progression 1.000000;
Transfinite Line {line1+ 25,line1+ 26,line1+ 27,line1+ 28} =  11 Using Progression 1.000000;

//**** Transfinite far field lines

Transfinite Line {line1+ 46,line1+ 47,line1+ 48,line1+ 49} =  21 Using Progression 1.000000;
Transfinite Line {line1+ 38,line1+ 39,line1+ 40,line1+ 41} =  13 Using Progression 1.000000;
Transfinite Line {line1+ 42,line1+ 43,line1+ 44,line1+ 45} =  13 Using Progression 1.000000;

//**** Transfinite surfaces

//Transfinite Surface {surf1+ 18} = {p1+  1,p1+ 12,p1+ 21,p1+ 20};
//Transfinite Surface {surf1+ 19} = {p1+  1,p1+ 14,p1+ 22,p1+ 20};
//Transfinite Surface {surf1+ 20} = {p1+  1,p1+ 16,p1+ 23,p1+ 20};
//Transfinite Surface {surf1+ 21} = {p1+  1,p1+ 18,p1+ 24,p1+ 20};

vol1 = newreg;
Surface Loop (vol1)= {surf1+  5,surf1+  6,surf1+  7,surf1+  8,surf1+ 10,surf1+ 12,surf1+ 13,surf1+ 24,surf1+ 26,surf1+ 28,surf1+ 30,surf1+ 32,surf1+ 34 };
Volume(vol1+1) ={vol1};

//**** Embed surfaces in volume

Surface {surf1+ 18} In Volume {vol1+1};
Surface {surf1+ 19} In Volume {vol1+1};
Surface {surf1+ 20} In Volume {vol1+1};
Surface {surf1+ 21} In Volume {vol1+1};


