Include "two_dim_data.pro"; 

// Enclosure points
p0 = newp; Point(p0) = {0, 0, 0, lc_air};
p1 = newp; Point(p1) = {r_enc, 0, 0, lc_air};
p2 = newp; Point(p2) = {0, r_enc, 0, lc_air};
p3 = newp; Point(p3) = {-r_enc, 0, 0, lc_air};
p4 = newp; Point(p4) = {0, -r_enc, 0, lc_air};

L1 = newreg; Circle(L1) = {p1, p0, p2};
L2 = newreg; Circle(L2) = {p2, p0, p3};
L3 = newreg; Circle(L3) = {p3, p0, p4};
L4 = newreg; Circle(L4) = {p4, p0, p1};

p5 = newp; Point(p5) = {-width/2, -thick/2, 0, lc_coil};
pp1[] = Extrude{width,0,0}{Point{p5};Layers{80};};
ppp1[] = Extrude{0,thick,0}{Line{pp1[1]};Layers{5};Recombine;};
LLTape1 = newreg; Line Loop(LLTape1) = {Boundary{Surface{ppp1[1]};}};
    
LLAir = newreg; Line Loop(LLAir) = {L1, L2, L3, L4};

air_dom = newreg; Plane Surface(air_dom) = {LLAir, LLTape1};

Physical Surface(AirDom) = {air_dom};
Physical Surface(TapeDom) = {ppp1[1]};
//For Index In {0:#ppp[]-1}
//  Printf("%f\n",ppp[Index]);
//EndFor
Physical Line(IntBnd) = {pp1[1],ppp1[0],ppp1[2],-ppp1[3]};
Physical Line(OutBnd) = {L1, L2, L3, L4};
Mesh.Algorithm = 1;
Mesh.Binary = 1;

// Cohomology computation for the H-Phi formulation
Cohomology(1) {AirDom, {}};

General.ExpertMode = 1;
Mesh.Optimize = 1;

Coherence;
Coherence;
