// Gmsh project created on Tue May 24 14:19:42 2014

//#include <math.h>

ld = 1; // diameter of cylinders
lx1 = 50;
lx2 = 50;
ly1  = 50;
ly2  = 50;

linnerRec = 4*ld;
ldistance = 2.0*ld; // spacing ratio
lcircular = 0.636396*ld; // diameter of the outer circle
thita = 45.0*Pi/180;
ratio = 0.001*ld;


//delta = (ldistance-0.5)/2*sin(45.0);
// cylinder 1
x11 = -ldistance/2; y11 = -ldistance/2; 
Point(11) = {x11*Cos(thita)-y11*Sin(thita), x11*Sin(thita)+y11*Cos(thita), 0, ratio};

x12 = -Sqrt( 2 )/4*ld-ldistance/2; y12 = Sqrt( 2 )/4*ld-ldistance/2; 
x13 =  Sqrt( 2 )/4*ld-ldistance/2; y13 = Sqrt( 2 )/4*ld-ldistance/2; 
x14 =  Sqrt( 2 )/4*ld-ldistance/2; y14 = -Sqrt( 2 )/4*ld-ldistance/2; 
x15 = -Sqrt( 2 )/4*ld-ldistance/2; y15 = -Sqrt( 2 )/4*ld-ldistance/2; 

Point(12) = {x12*Cos(thita)-y12*Sin(thita), x12*Sin(thita)+y12*Cos(thita), 0, ratio};
Point(13) = {x13*Cos(thita)-y13*Sin(thita), x13*Sin(thita)+y13*Cos(thita), 0, ratio};
Point(14) = {x14*Cos(thita)-y14*Sin(thita), x14*Sin(thita)+y14*Cos(thita), 0, ratio};
Point(15) = {x15*Cos(thita)-y15*Sin(thita), x15*Sin(thita)+y15*Cos(thita), 0, ratio};

x16 = -lcircular-ldistance/2; y16 = lcircular-ldistance/2; 
x17 = lcircular-ldistance/2; y17 = lcircular-ldistance/2; 
x18 = lcircular-ldistance/2; y18 = -lcircular-ldistance/2; 
x19 = -lcircular-ldistance/2; y19 = -lcircular-ldistance/2; 

Point(16) = {x16*Cos(thita)-y16*Sin(thita), x16*Sin(thita)+y16*Cos(thita), 0, ratio};
Point(17) = {x17*Cos(thita)-y17*Sin(thita), x17*Sin(thita)+y17*Cos(thita), 0, ratio};
Point(18) = {x18*Cos(thita)-y18*Sin(thita), x18*Sin(thita)+y18*Cos(thita), 0, ratio};
Point(19) = {x19*Cos(thita)-y19*Sin(thita), x19*Sin(thita)+y19*Cos(thita), 0, ratio};


// cylinder 3

x31 = +ldistance/2; y31 = +ldistance/2; 
Point(31) = {x31*Cos(thita)-y31*Sin(thita), x31*Sin(thita)+y31*Cos(thita), 0, ratio};

x32 = -Sqrt( 2 )/4*ld+ldistance/2; y32 = Sqrt( 2 )/4*ld+ldistance/2; 
x33 =  Sqrt( 2 )/4*ld+ldistance/2; y33 = Sqrt( 2 )/4*ld+ldistance/2; 
x34 =  Sqrt( 2 )/4*ld+ldistance/2; y34 = -Sqrt( 2 )/4*ld+ldistance/2; 
x35 = -Sqrt( 2 )/4*ld+ldistance/2; y35 = -Sqrt( 2 )/4*ld+ldistance/2; 

Point(32) = {x32*Cos(thita)-y32*Sin(thita), x32*Sin(thita)+y32*Cos(thita), 0, ratio};
Point(33) = {x33*Cos(thita)-y33*Sin(thita), x33*Sin(thita)+y33*Cos(thita), 0, ratio};
Point(34) = {x34*Cos(thita)-y34*Sin(thita), x34*Sin(thita)+y34*Cos(thita), 0, ratio};
Point(35) = {x35*Cos(thita)-y35*Sin(thita), x35*Sin(thita)+y35*Cos(thita), 0, ratio};

x36 = -lcircular+ldistance/2; y36 = lcircular+ldistance/2; 
x37 = lcircular+ldistance/2; y37 = lcircular+ldistance/2; 
x38 = lcircular+ldistance/2; y38 = -lcircular+ldistance/2; 
x39 = -lcircular+ldistance/2; y39 = -lcircular+ldistance/2; 

Point(36) = {x36*Cos(thita)-y36*Sin(thita), x36*Sin(thita)+y36*Cos(thita), 0, ratio};
Point(37) = {x37*Cos(thita)-y37*Sin(thita), x37*Sin(thita)+y37*Cos(thita), 0, ratio};
Point(38) = {x38*Cos(thita)-y38*Sin(thita), x38*Sin(thita)+y38*Cos(thita), 0, ratio};
Point(39) = {x39*Cos(thita)-y39*Sin(thita), x39*Sin(thita)+y39*Cos(thita), 0, ratio};


// outer rectangular

lcos=ldistance*Cos( thita );
lsin=ldistance*Sin( thita );

Point(1) = {-linnerRec-lcos/2, -linnerRec-lsin/2, 0, ratio};
Point(2) = {linnerRec+lcos/2, -linnerRec-lsin/2, 0, ratio};
Point(3) = {linnerRec+lcos/2, linnerRec+lsin/2, 0, ratio};
Point(4) = {-linnerRec-lcos/2, linnerRec+lsin/2, 0, ratio};

// cylinder 1
Circle(1) = {12, 11, 13};
Circle(2) = {13, 11, 14};
Circle(3) = {14, 11, 15};
Circle(4) = {15, 11, 12};
Circle(5) = {16, 11, 17};
Circle(6) = {17, 11, 18};
Circle(7) = {18, 11, 19};
Circle(8) = {19, 11, 16};


// cylinder 3
Circle(17) = {32, 31, 33};
Circle(18) = {33, 31, 34};
Circle(19) = {34, 31, 35};
Circle(20) = {35, 31, 32};
Circle(21) = {36, 31, 37};
Circle(22) = {37, 31, 38};
Circle(23) = {38, 31, 39};
Circle(24) = {39, 31, 36};



// inner rec & related
Line(49) = {4, 1};
Line(50) = {4, 3};
Line(51) = {3, 2};
Line(52) = {2, 1};


// circle
Transfinite Line{1:8, 17:24} = 5;

// rec
Transfinite Line{49:52} = 10;//40 Using Bump 2.5;

// inner to outer circles

// cylinder 3
Line(61) = {32, 36};
Line(62) = {33, 37};
Line(63) = {34, 38};
Line(64) = {35, 39};

// cylinder 1
Line(65) = {12, 16};
Line(66) = {13, 17};
Line(67) = {14, 18};
Line(68) = {15, 19};

Transfinite Line{61:68} = 5 Using Progression 1.086;

// cylinder 3, surfaces
Line Loop(87) = {21, -62, -17, 61};
Plane Surface(88) = {87};
Line Loop(89) = {62, 22, -63, -18};
Plane Surface(90) = {89};
Line Loop(91) = {19, 64, -23, -63};
Plane Surface(92) = {91};
Line Loop(93) = {20, 61, -24, -64};
Plane Surface(94) = {93};


// cylinder 1, surfaces
Line Loop(103) = {66, 6, -67, -2};
Plane Surface(104) = {103};
Line Loop(105) = {3, 68, -7, -67};
Plane Surface(106) = {105};
Line Loop(107) = {4, 65, -8, -68};
Plane Surface(108) = {107};
Line Loop(109) = {65, 5, -66, -1};
Plane Surface(110) = {109};

Transfinite Surface{88:94:2, 104:110:2};
Recombine Surface{88:94:2, 104:110:2};

// inner rect
Line Loop(81) = {50, 51, 52, -49};
Line Loop(82) = {21, 22, 23, 24};
Line Loop(84) = {5, 6, 7, 8};
Plane Surface(86) = {81, 82, 84};
Recombine Surface{86};
//Mesh.SubdivisionAlgorithm=2;
//Mesh.Smoothing = 0;
Smoother Surface{86} = 0;

