%Units:
mm = 0.001;

%Shrinkages (293 K --> 4.2 K):
spphr_para_shrnkg = 1.0-721e-6;
spphr_perp_shrnkg = 1.0-599e-6;
cppr_shrnkg = 1-41.4e-4;
stnlss_stl_shrnkg = 1.0 %1-29.7e-4;

%Dimensions:
belly_diameter = spphr_perp_shrnkg * 50.0 * mm
neck_diameter = spphr_perp_shrnkg * 20.0 * mm
radial_lwr_cutoff = 0.1 * mm
neck_height =  spphr_para_shrnkg * 8.0 * mm
slope_height = spphr_para_shrnkg * 5.0 * mm
belly_height = spphr_para_shrnkg * 20.00 * mm 
can_diameter = cppr_shrnkg * 90.0 * mm
can_height = cppr_shrnkg * 50.0 * mm

%Coordinate macros
%radial:
r_neck = neck_diameter /2;
r_belly = belly_diameter /2;
r_can = can_diameter /2;

%axial:
h_lwr_edge = neck_height + slope_height;
h_uppr_edge = h_lwr_edge + belly_height;
h_uppr_neck = h_uppr_edge + slope_height;

100 * belly_diameter
100 * neck_diameter
100 * belly_height
100 * (can_height - h_uppr_neck)

%Point coordinates:
p1 = [radial_lwr_cutoff, 0.0];
p2 = [radial_lwr_cutoff, neck_height];
p3 = [radial_lwr_cutoff, h_uppr_neck];
p4 = [radial_lwr_cutoff, can_height];

p5 = [r_neck, 0.0];
p6 = [r_neck, neck_height];
p7 = [r_neck, h_uppr_neck];
p8 = [r_neck, can_height];

p9 = [r_belly, 0.0];
p10= [r_belly, h_lwr_edge];
p11 = [r_belly, h_uppr_edge];
p12 = [r_belly, can_height];

p13 = [r_can, 0.0];
p14 = [r_can, h_lwr_edge];
p15 = [r_can, h_uppr_edge];
p16 = [r_can, can_height];

%Subdomains:
s1=poly2vert([p1;p5;p6;p2]);
s2=poly2vert([p2;p6;p7;p3]);
s3=poly2vert([p3;p7;p8;p4]);
v1=poly2vert([p5;p9;p10;p6]);
s4=poly2vert([p6;p10;p11;p7]);
v2=poly2vert([p7;p11;p12;p8]);
v3=poly2vert([p9;p13;p14;p10]);
v4=poly2vert([p10;p14;p15;p11]);
v5=poly2vert([p11;p15;p16;p12]);

geom=geomcsg({s1},{s2},{s3},{s4},{s4},{v1},{v2},{v3},{v4},{v5})

geomplot(geom)

strng = 's1 s2 s3 s4 v1 v2 v3 v4 v5'
