// Benchmark 1 Wang et al., Int. J. Solids. Struc., 51, p.1132-1148, 2014

lx=50.*1.e-3;
ly=lx;
lz=100.*1.e-3;
lyh=ly/2.;
lxh=lyh;
mesh_size = lz/60.;

// Flaw parameters:
a_flaw=24.*1.e-3/2.;
b_flaw=33.*1.e-3/2.;
alpha_flaw = (30./180.)*Pi;
beta_flaw = (38./180.)*Pi;
width_flaw = 0.5*1.e-3/2.; // 5.*1.e-3/2.;  // 0.5*1.e-3/2.; 

He = (lz - 4*a_flaw*Cos(alpha_flaw)-2.*b_flaw*Sin(beta_flaw))/2.;
ycut = Sqrt(lyh^2 - (lyh - 2.*a_flaw*Sin(alpha_flaw))^2);


P1=newp; Point(P1) = {0,0,0,mesh_size};
P2=newp; Point(P2) = {lyh-2.*a_flaw*Sin(alpha_flaw),ycut,0.,mesh_size};
// P3=newp; Point(P3) = {-lyh,0.,0.,mesh_size};

P4=newp; Point(P4) = {0,0,lz,mesh_size};
P5=newp; Point(P5) = {lyh-2.*a_flaw*Sin(alpha_flaw),ycut,lz,mesh_size};
// P6=newp; Point(P6) = {-lyh,0.,lz,mesh_size};

P7=newp; Point(P7) = {lyh-2.*a_flaw*Sin(alpha_flaw),-ycut,0.,mesh_size};
P8=newp; Point(P8) = {lyh-2.*a_flaw*Sin(alpha_flaw),-ycut,lz,mesh_size};


P9=newp; Point(P9) = {-lyh+2.*a_flaw*Sin(alpha_flaw),-ycut,0,mesh_size};
P10=newp; Point(P10) = {-lyh+2.*a_flaw*Sin(alpha_flaw),-ycut,lz,mesh_size};
P11=newp; Point(P11) = {-lyh+2.*a_flaw*Sin(alpha_flaw),ycut,0.,mesh_size};
P12=newp; Point(P12) = {-lyh+2.*a_flaw*Sin(alpha_flaw),ycut,lz,mesh_size};

Pmid1=newp; Point(Pmid1) = {0,0,He + 2.*a_flaw*Cos(alpha_flaw) - width_flaw,mesh_size};
P13=newp; Point(P13) = {lyh-2.*a_flaw*Sin(alpha_flaw),ycut,He + 2.*a_flaw*Cos(alpha_flaw) - width_flaw,mesh_size};

Pmid2=newp; Point(Pmid2) = {0,0,He + 2.*a_flaw*Cos(alpha_flaw) + width_flaw,mesh_size};
P14=newp; Point(P14) = {lyh-2.*a_flaw*Sin(alpha_flaw),ycut,He + 2.*a_flaw*Cos(alpha_flaw) + width_flaw,mesh_size};

P15=newp; Point(P15) = {lyh-2.*a_flaw*Sin(alpha_flaw),-ycut,He + 2.*a_flaw*Cos(alpha_flaw) - width_flaw,mesh_size};
P16=newp; Point(P16) = {lyh-2.*a_flaw*Sin(alpha_flaw),-ycut,He + 2.*a_flaw*Cos(alpha_flaw) + width_flaw,mesh_size};


Pmid3=newp; Point(Pmid3) = {0,0,lz-He-2.*a_flaw*Cos(alpha_flaw) - width_flaw,mesh_size};
Pmid4=newp; Point(Pmid4) = {0,0,lz-He-2.*a_flaw*Cos(alpha_flaw) + width_flaw,mesh_size};


P18=newp; Point(P18) = {-lyh+2.*a_flaw*Sin(alpha_flaw),ycut,lz-He-2.*a_flaw*Cos(alpha_flaw) - width_flaw,mesh_size};
P19=newp; Point(P19) = {-lyh+2.*a_flaw*Sin(alpha_flaw),ycut,lz-He-2.*a_flaw*Cos(alpha_flaw) + width_flaw,mesh_size};

P20=newp; Point(P20) = {-lyh+2.*a_flaw*Sin(alpha_flaw),-ycut,lz-He-2.*a_flaw*Cos(alpha_flaw) - width_flaw,mesh_size};
P21=newp; Point(P21) = {-lyh+2.*a_flaw*Sin(alpha_flaw),-ycut,lz-He-2.*a_flaw*Cos(alpha_flaw) + width_flaw,mesh_size};

L1314Join=newl; Line(L1314Join) = {P13,P14};
L1516Join=newl; Line(L1516Join) = {P15,P16};
L1819Join=newl; Line(L1819Join) = {P18,P19};
L2021Join=newl; Line(L2021Join) = {P20,P21};

C1 = newl; Circle(C1) = {P2,P1,P11};
C2 = newl; Circle(C2) = {P5,P4,P12};

C3 = newl; Circle(C3) = {P9,P1,P7};
C4 = newl; Circle(C4) = {P10,P4,P8};

L4 = newl; Line(L4) = {P2,P5};
L10 = newl; Line(L10) = {P7,P8};



// defining flaw 1 at front

// mesh for flaw 1 
nmesh=100;
lmesh=2*a_flaw/(nmesh*1.);

Pcut1Up = newp; Point(Pcut1Up) = {lyh,0,lz-He + width_flaw,lmesh};
Pcut1Down = newp; Point(Pcut1Down) = {lyh,0,lz-He - width_flaw,lmesh};
Pcut2Up = newp; Point(Pcut2Up) = {lyh-2.*a_flaw*Sin(alpha_flaw),-ycut,lz-He-2.*a_flaw*Cos(alpha_flaw) + width_flaw,lmesh};
Pcut2Down = newp; Point(Pcut2Down) = {lyh-2.*a_flaw*Sin(alpha_flaw),-ycut,lz-He-2.*a_flaw*Cos(alpha_flaw) - width_flaw,lmesh};
Pcut3Up = newp; Point(Pcut3Up) = {lyh-2.*a_flaw*Sin(alpha_flaw),ycut,lz-He-2.*a_flaw*Cos(alpha_flaw) + width_flaw,lmesh};
Pcut3Down = newp; Point(Pcut3Down) = {lyh-2.*a_flaw*Sin(alpha_flaw),ycut,lz-He-2.*a_flaw*Cos(alpha_flaw) - width_flaw,lmesh};

L4a1 = newl; Line(L4a1) ={P2,P13};
L4a2 = newl; Line(L4a2) ={P14,Pcut3Down};
L4b = newl; Line(L4b) ={Pcut3Up,P5};

L10a1 = newl; Line(L10a1) ={P7,P15};
L10a2 = newl; Line(L10a2) = {P16,Pcut2Down};
L10b = newl; Line(L10b) ={Pcut2Up,P8};


// finding spline points for curve joining Pcut1Up/Down with Pcut3UpDown 


For i In {1:nmesh-1}
    If (i==1)
       PSpline13Up[i-1] = newp;
       Point(PSpline13Up[i-1]) = { lyh-lmesh*(i*1./2.)*Sin(alpha_flaw),Sqrt(lyh^2 - (lyh - lmesh*(i*1./2.)*Sin(alpha_flaw))^2), lz - He + width_flaw - lmesh*(i*1./2.)*Cos(alpha_flaw),lmesh};  
    EndIf

    PSpline13Up[i] = newp;
    Point(PSpline13Up[i]) = { lyh-lmesh*(i*1.)*Sin(alpha_flaw),Sqrt(lyh^2 - (lyh - lmesh*(i*1.)*Sin(alpha_flaw))^2), lz - He + width_flaw - lmesh*(i*1.)*Cos(alpha_flaw),lmesh}; 

EndFor


For i In {1:nmesh-1}
    If (i==1)
       PSpline13Down[i-1] = newp;
       Point(PSpline13Down[i-1]) = { lyh-lmesh*(i*1./2.)*Sin(alpha_flaw),Sqrt(lyh^2 - (lyh - lmesh*(i*1./2.)*Sin(alpha_flaw))^2), lz - He - width_flaw - lmesh*(i*1./2.)*Cos(alpha_flaw),lmesh};
    EndIf

    PSpline13Down[i] = newp;
    Point(PSpline13Down[i]) = { lyh-lmesh*(i*1.)*Sin(alpha_flaw),Sqrt(lyh^2 - (lyh - lmesh*(i*1.)*Sin(alpha_flaw))^2), lz - He - width_flaw - lmesh*(i*1.)*Cos(alpha_flaw),lmesh};

EndFor


// finding spline points for curve joining Pcut2Up/Down with Pcut1UpDown 

j=0;
For i In {nmesh-1:1:-1}

    PSpline21Up[j] = newp;
    Point(PSpline21Up[j]) = { lyh-lmesh*(i*1.)*Sin(alpha_flaw),-Sqrt(lyh^2 - (lyh - lmesh*(i*1.)*Sin(alpha_flaw))^2), lz - He + width_flaw - lmesh*(i*1.)*Cos(alpha_flaw),lmesh};

    If (i==1)
       PSpline21Up[j+1] = newp;
       Point(PSpline21Up[j+1]) = { lyh-lmesh*(i*1./2.)*Sin(alpha_flaw),-Sqrt(lyh^2 - (lyh - lmesh*(i*1./2.)*Sin(alpha_flaw))^2), lz - He + width_flaw - lmesh*(i*1./2.)*Cos(alpha_flaw),lmesh};
    EndIf
    j += 1; 


EndFor

j=0;
For i In {nmesh-1:1:-1}

    PSpline21Down[j] = newp;
    Point(PSpline21Down[j]) = { lyh-lmesh*(i*1.)*Sin(alpha_flaw),-Sqrt(lyh^2 - (lyh - lmesh*(i*1.)*Sin(alpha_flaw))^2), lz - He - width_flaw - lmesh*(i*1.)*Cos(alpha_flaw),lmesh};

    If (i==1)
       PSpline21Down[j+1] = newp;
       Point(PSpline21Down[j+1]) = { lyh-lmesh*(i*1./2.)*Sin(alpha_flaw),-Sqrt(lyh^2 - (lyh - lmesh*(i*1./2.)*Sin(alpha_flaw))^2), lz - He - width_flaw - lmesh*(i*1./2.)*Cos(alpha_flaw),lmesh};
    EndIf
    j += 1;

EndFor


// Complete splines joining Pcut2Up/Down with Pcut3Up/Down

Sp23Up=newl; Spline(Sp23Up) = {Pcut2Up,PSpline21Up[],Pcut1Up,PSpline13Up[],Pcut3Up};
Sp23Down=newl; Spline(Sp23Down) = {Pcut2Down,PSpline21Down[],Pcut1Down,PSpline13Down[],Pcut3Down};

L32Up = newl; Line(L32Up) = {Pcut3Up,Pcut2Up};
L32Down = newl; Line(L32Down) = {Pcut3Down,Pcut2Down};

L2Join = newl; Line(L2Join) = {Pcut2Down,Pcut2Up};
L3Join = newl; Line(L3Join) = {Pcut3Down,Pcut3Up};

llcut1 = newll;
Line Loop(llcut1) = {Sp23Up,L32Up};
Scut1 = news; 
Plane Surface(Scut1) = {llcut1};

llcut2 = newll;
Line Loop(llcut2) = {Sp23Down,L32Down};
Scut2 = news;
Plane Surface(Scut2) = {-llcut2};


llcut3 = newll;
Line Loop(llcut3) = {Sp23Down,L3Join,-Sp23Up,-L2Join};
Scut3 = news;
Ruled Surface(Scut3) = {llcut3};

llcut4 = newll;
Line Loop(llcut4) = {-L32Down,L3Join,L32Up,-L2Join};
Scut4 = news;
Plane Surface(Scut4) = {-llcut4};

sFlaw1 = newreg;
Surface Loop(sFlaw1) = {Scut1,Scut2,Scut3,Scut4};
vFlaw1 = newreg;

Volume(vFlaw1) = {sFlaw1};
Physical Volume("Flaw1") = {vFlaw1};


//  mesh for flaw 2


Pcut4Up = newp; Point(Pcut4Up) = {-lyh,0,He+width_flaw,lmesh};
Pcut4Down = newp; Point(Pcut4Down) = {-lyh,0,He-width_flaw,lmesh};

Pcut5Up = newp; Point(Pcut5Up) = {-lyh + 2.*a_flaw*Sin(alpha_flaw),-ycut, He+2.*a_flaw*Cos(alpha_flaw)+width_flaw,lmesh};
Pcut5Down = newp; Point(Pcut5Down) = {-lyh + 2.*a_flaw*Sin(alpha_flaw),-ycut, He+2.*a_flaw*Cos(alpha_flaw)-width_flaw,lmesh};

Pcut6Up = newp; Point(Pcut6Up) = {-lyh + 2.*a_flaw*Sin(alpha_flaw),ycut, He+2.*a_flaw*Cos(alpha_flaw)+width_flaw,lmesh};
Pcut6Down = newp; Point(Pcut6Down) = {-lyh + 2.*a_flaw*Sin(alpha_flaw),ycut, He+2.*a_flaw*Cos(alpha_flaw)-width_flaw,lmesh};

L13a = newl; Line(L13a) = {P9,Pcut5Down};
L5Join = newl; Line(L5Join) = {Pcut5Down,Pcut5Up};
L13b1 = newl; Line(L13b1) = {Pcut5Up,P20};
L13b2 = newl; Line(L13b2) = {P21,P10};

L14a = newl; Line(L14a) = {P11,Pcut6Down};
L6Join = newl; Line(L6Join) = {Pcut6Down,Pcut6Up};
L14b1 = newl; Line(L14b1) = {Pcut6Up,P18};

L14b2 = newl; Line(L14b2) = {P19,P12};


L65Up = newl; Line(L65Up) = {Pcut6Up,Pcut5Up};
L65Down = newl; Line(L65Down) = {Pcut6Down,Pcut5Down};

// finding spline points for curve joining Pcut5Up/Down with Pcut4UpDown


For i In {1:nmesh-1}

    PSpline54Up[i-1] = newp;
    Point(PSpline54Up[i-1]) = { -lyh + 2.*a_flaw*Sin(alpha_flaw) -lmesh*(i*1.)*Sin(alpha_flaw),-Sqrt(lyh^2 - (lyh - 2.*a_flaw*Sin(alpha_flaw) +  lmesh*(i*1.)*Sin(alpha_flaw))^2), He + 2.*a_flaw*Cos(alpha_flaw) + width_flaw - lmesh*(i*1.)*Cos(alpha_flaw),lmesh};

    If (i==nmesh-1)
       PSpline54Up[i] = newp;
       Point(PSpline54Up[i]) = { -lyh + 2.*a_flaw*Sin(alpha_flaw) - lmesh*(i+1./2.)*Sin(alpha_flaw), -Sqrt(lyh^2 - (lyh - 2.*a_flaw*Sin(alpha_flaw) + lmesh*(i+1./2.)*Sin(alpha_flaw))^2), He + 2.*a_flaw*Cos(alpha_flaw) + width_flaw - lmesh*(i+1./2.)*Cos(alpha_flaw),lmesh};
    EndIf

EndFor

For i In {1:nmesh-1}

    PSpline54Down[i-1] = newp;
    Point(PSpline54Down[i-1]) = { -lyh + 2.*a_flaw*Sin(alpha_flaw) -lmesh*(i*1.)*Sin(alpha_flaw),-Sqrt(lyh^2 - (lyh - 2.*a_flaw*Sin(alpha_flaw) +  lmesh*(i*1.)*Sin(alpha_flaw))^2), He + 2.*a_flaw*Cos(alpha_flaw) - width_flaw - lmesh*(i*1.)*Cos(alpha_flaw),lmesh};

    If (i==nmesh-1)
       PSpline54Down[i] = newp;
       Point(PSpline54Down[i]) = { -lyh + 2.*a_flaw*Sin(alpha_flaw) - lmesh*(i+1./2.)*Sin(alpha_flaw), -Sqrt(lyh^2 - (lyh - 2.*a_flaw*Sin(alpha_flaw) + lmesh*(i+1./2.)*Sin(alpha_flaw))^2), He + 2.*a_flaw*Cos(alpha_flaw) - width_flaw - lmesh*(i+1./2.)*Cos(alpha_flaw),lmesh};
    EndIf

EndFor

// finding spline points for curve joining Pcut4Up/Down with Pcut6UpDown


For i In {1:nmesh-1}

    If (i==1)
       PSpline46Up[i-1] = newp;
       Point(PSpline46Up[i-1]) = { -lyh + lmesh*(i*1./2.)*Sin(alpha_flaw), Sqrt(lyh^2 - (lyh - lmesh*(i*1./2.)*Sin(alpha_flaw))^2), He + width_flaw + lmesh*(i*1./2.)*Cos(alpha_flaw),lmesh};
    EndIf


    PSpline46Up[i] = newp;
    Point(PSpline46Up[i]) = { -lyh + lmesh*(i*1.)*Sin(alpha_flaw),Sqrt(lyh^2 - (lyh - lmesh*(i*1.)*Sin(alpha_flaw))^2), He + width_flaw + lmesh*(i*1.)*Cos(alpha_flaw),lmesh};


EndFor

For i In {1:nmesh-1}

    If (i==1)
       PSpline46Down[i-1] = newp;
       Point(PSpline46Down[i-1]) = { -lyh + lmesh*(i*1./2.)*Sin(alpha_flaw), Sqrt(lyh^2 - (lyh - lmesh*(i*1./2.)*Sin(alpha_flaw))^2), He - width_flaw + lmesh*(i*1./2.)*Cos(alpha_flaw),lmesh};
    EndIf


    PSpline46Down[i] = newp;
    Point(PSpline46Down[i]) = { -lyh + lmesh*(i*1.)*Sin(alpha_flaw),Sqrt(lyh^2 - (lyh - lmesh*(i*1.)*Sin(alpha_flaw))^2), He - width_flaw + lmesh*(i*1.)*Cos(alpha_flaw),lmesh};


EndFor

// Complete splines joining Pcut2Up/Down with Pcut3Up/Down

Sp56Up=newl; Spline(Sp56Up) = {Pcut5Up,PSpline54Up[],Pcut4Up,PSpline46Up[],Pcut6Up};
Sp56Down=newl; Spline(Sp56Down) = {Pcut5Down,PSpline54Down[],Pcut4Down,PSpline46Down[],Pcut6Down};

llcut5 = newll;
Line Loop(llcut5) = {Sp56Up,L65Up};
Scut5 = news;
Plane Surface(Scut5) = {-llcut5};

llcut6 = newll;
Line Loop(llcut6) = {Sp56Down,L65Down};
Scut6 = news;
Plane Surface(Scut6) = {llcut6};

llcut7 = newll;
Line Loop(llcut7) = {Sp56Down,L6Join,-Sp56Up,-L5Join};
Scut7 = news;
Ruled Surface(Scut7) = {-llcut7};

llcut8 = newll;
Line Loop(llcut8) = {-L65Down,L6Join,L65Up,-L5Join};
Scut8 = news;
Plane Surface(Scut8) = {llcut8};

sFlaw2 = newreg;
Surface Loop(sFlaw2) = {Scut5,Scut6,Scut7,Scut8};
vFlaw2 = newreg;

Volume(vFlaw2) = {sFlaw2};
Physical Volume("Flaw2") = {vFlaw2};






// adding mesh around flaw
// vInner


// need to break sInner into multiple regions as Ruled Surface can have max 4 lines in line loop

C5 = newl; Circle(C5) = {P13,Pmid1,Pcut6Down};

llInner5 = newll;
Line Loop(llInner5) = {C1,L14a,-C5,-L4a1};
sIntact1 = newreg;
Ruled Surface(sIntact1) = {llInner5};

C6 = newl; Circle(C6) = {P14,Pmid2,Pcut6Up};

llInner6 = newll;
Line Loop(llInner6) = {C5,L6Join,-C6,-L1314Join};
sIntact2 = newreg;
Ruled Surface(sIntact2) = {llInner6};

C7 = newl; Circle(C7) = {P11,P1,P9};

llInner7 = newll;
Line Loop(llInner7) = {C7,L13a,Sp56Down,-L14a};
sIntact3 = newreg;
Ruled Surface(sIntact3) = {llInner7};

C8 = newl; Circle(C8) = {Pcut5Down,Pmid1,P15};

llInner8 = newll;
Line Loop(llInner8) = {C3,L10a1,-C8,-L13a};
sIntact4 = newreg;
Ruled Surface(sIntact4) = {llInner8};

C9 = newl; Circle(C9) = {P7,P1,2};
C10 = newl; Circle(C10) = {P15,Pmid1,P13};
llInner9 = newll;
Line Loop(llInner9) = {C9,L4a1,-C10,-L10a1};
sIntact5 = newreg;
Ruled Surface(sIntact5) = {llInner9};

C11 = newl; Circle(C11) = {Pcut3Down,Pmid3,P18};

llInner10 = newll;
Line Loop(llInner10) = {C6,L14b1,-C11,-L4a2};
sIntact6 = newreg;
Ruled Surface(sIntact6) = {llInner10};


C12 = newl; Circle(C12) = {Pcut3Up,Pmid4,P19};

llInner11 = newll;
Line Loop(llInner11) = {C11,L1819Join,-C12,-L3Join};
sIntact7 = newreg;
Ruled Surface(sIntact7) = {llInner11};

llInner12 = newll;
Line Loop(llInner12) = {C12,L14b2,-C2,-L4b};
sIntact8 = newreg;
Ruled Surface(sIntact8) = {llInner12};

C13 = newl; Circle(C13) = {P18,Pmid3,P20};
C14 = newl; Circle(C14) = {P19,Pmid4,P21};
C15 = newl; Circle(C15) = {P12,P4,P10};


llInner13 = newll;
Line Loop(llInner13) = {-Sp56Up,L13b1,-C13,-L14b1};
sIntact9 = newreg;
Ruled Surface(sIntact9) = {llInner13};

llInner14 = newll;
Line Loop(llInner14) = {C13,L2021Join,-C14,-L1819Join};
sIntact10 = newreg;
Ruled Surface(sIntact10) = {llInner14};

llInner15 = newll;
Line Loop(llInner15) = {C14,L13b2,-C15,-L14b2};
sIntact11 = newreg;
Ruled Surface(sIntact11) = {llInner15};

C16 = newl; Circle(C16) = {Pcut5Up,Pmid2,P16};
C17 = newl; Circle(C17) = {P20,Pmid3,Pcut2Down};
C18 = newl; Circle(C18) = {P21,Pmid4,Pcut2Up};

llInner16 = newll;
Line Loop(llInner16) = {C8,L1516Join,-C16,-L5Join};
sIntact12 = newreg;
Ruled Surface(sIntact12) = {llInner16};

llInner17 = newll;
Line Loop(llInner17) = {C16,L10a2,-C17,-L13b1};
sIntact13 = newreg;
Ruled Surface(sIntact13) = {llInner17};

llInner18 = newll;
Line Loop(llInner18) = {C17,L2Join,-C18,-L2021Join};
sIntact14 = newreg;
Ruled Surface(sIntact14) = {llInner18};

llInner19 = newll;
Line Loop(llInner19) = {C18,L10b,-C4,-L13b2};
sIntact15 = newreg;
Ruled Surface(sIntact15) = {llInner19};

C19 = newl; Circle(C19) = {P16,Pmid2,P14};
C20 = newl; Circle(C20) = {P8,P4,P5};

llInner20 = newll;
Line Loop(llInner20) = {C10,L1314Join,-C19,-L1516Join};
sIntact16 = newreg;
Ruled Surface(sIntact16) = {llInner20};

llInner21 = newll;
Line Loop(llInner21) = {C19,L4a2,-Sp23Down,-L10a2};
sIntact17 = newreg;
Ruled Surface(sIntact17) = {llInner21};

llInner22 = newll;
Line Loop(llInner22) = {Sp23Up,L4b,-C20,-L10b};
sIntact18 = newreg;
Ruled Surface(sIntact18) = {llInner22};

Physical Surface("Sides") = {Scut3,Scut7,sIntact5, sIntact16,sIntact17,sIntact18,sIntact4,sIntact12,sIntact13,sIntact14,sIntact15,sIntact3,sIntact9,sIntact10,sIntact11,sIntact1,sIntact2,sIntact6,sIntact7,sIntact8};


llTop = newll;
Line Loop(llTop) = {C2,C15,C4,C20};
sTop = newreg;
Plane Surface(sTop) = {llTop};
Physical Surface("Top") = {sTop};


llBottom = newll;
Line Loop(llBottom) = {-C1,-C9,-C3,-C7};
sBottom = newreg;
Plane Surface(sBottom) = {llBottom};
Physical Surface("Bottom") = {sBottom};

sInner = newreg;
Surface Loop(sInner) = {sIntact5, sIntact16,sIntact17,-Scut1,-Scut2,-Scut4,sIntact18,sIntact4,sIntact12,sIntact13,sIntact14,sIntact15,sIntact3,-Scut5,-Scut6,-Scut8,sIntact9,sIntact10,sIntact11,sIntact1,sIntact2,sIntact6,sIntact7,sIntact8,sTop,sBottom};

vInner = newreg;

Volume(vInner) = {sInner};
Physical Volume("Inner") = {vInner};

