This example demonstrates the compression modelling of a volunteer's left
breast. The contact between the breast and the compression plates is defined as
frictionless contact. The contact between the breast and the rib is defined as
tied contact. The tied contact is a modifed version of frictional contact, and
it does not involve computation of the slip function (i.e. only stick
condition). The elastic slip is controlled by the tied stiffness value (higher
the value, more the stick condition is enforced). No fixed displacement B.C. is
necessary as the breast is 'attached' to the rib through tied contact. The
compression plates and the rib are all defined in region 2. Like example 5J2,
this example does not solve for the entire load steps, but only solves for the
last load step, by reading in previous load step solution (ipinit) and its contact
pressure data file (ipdata).
if (!defined $example)
{
$example ='.';
}
set echo
$RUN_LOAD=10; # Simulation starts from load step 10 (=1 if solving from the beginning)
$CONT_PTS_1=4;
$CONT_PTS_2=4;
$load_step=10; # total number of load step
$load_step_size=2;
$BREAST=1;
$PLATE=2;
fem def para;r;contact;example
fem def coord 3,1
fem def region;r;region;example
fem def base;r;Base;example # need to define all basis functions in one file, otherwise gives error
fem def node;r;Initial_xi3;example reg $BREAST
fem def elem;r;Initial_xi3;example reg $BREAST
fem def fibr;r;Initial_xi3;example reg $BREAST
fem def elem;r;Initial_xi3;example fib reg $BREAST
#fem export node;Initial as Initial reg $BREAST
#fem export elem;Initial as Initial reg $BREAST
fem def node;r;Pad;example reg $PLATE
fem def elem;r;Pad;example reg $PLATE
fem def fibr;r;Pad;example reg $PLATE
fem def elem;r;Pad;example fib reg $PLATE
#fem export node;Pad as Pad reg $PLATE
#fem export elem;Pad as Pad reg $PLATE
# Breast node, elem & face grouping #
fem group elem all external s1=0 as XXX_1
fem group elem all external s3=0 as XXX_3
fem group node xi1=0 elem XXX_1 as FIX # skin & rib
fem group face all xi3 low external elem XXX_3 as RIB_FACE
fem group face 145,149,152,155,158,162,179,196,200,203,206,209,212,214,216,218 as BREAST_FACE
fem group node 1001..1016 as MASTER_MOVE reg $PLATE
fem group node 1113..1120 as MASTER_IMMOBILE reg $PLATE
fem group node 1005..1006,1011..1016,1113..1116,1218,1230,1298,1306 as MASTER_FIX reg $PLATE
fem group elem 1001..1004 as PLATE_ELEM reg $PLATE
fem group face all xi3 low external elem PLATE_ELEM as PLATE_FACE reg $PLATE
fem group elem 1005 as RIB_ELEM reg $PLATE
fem group face all xi3 high external elem RIB_ELEM as RIB_MASTER_FACE reg $PLATE
# Compression plate alignment
fem change nodes translate by -3,0,0 node MASTER_IMMOBILE reg $PLATE
fem change nodes translate by 0,10,0 node MASTER_IMMOBILE reg $PLATE # align
fem change nodes translate by -9,6,0 node MASTER_MOVE reg $PLATE # align
fem change nodes translate by 0,3,0 node MASTER_MOVE reg $PLATE # align to match with comp2 data (don't need this for comp1)
fem def equa;r;slave;example reg $BREAST
fem def mate;r;slave;example reg $BREAST
if($RUN_LOAD==1){
# if starting from the very beginning
fem def init;r;slave;example reg $BREAST
} else {
# if starting from an intermediate solution
$temp=$RUN_LOAD-1;
$name=sprintf("Compression.Load.%02d", $temp);
fem def init;r;$name;example
# Update YP(ny,3) with XP
fem update solution geom YP_index 3 reg $BREAST
# Read in contact pressure data from the previous load step
fem def data;r;$name;example contact
}
fem def equa;r;master;example reg $PLATE
fem def mate;r;master;example reg $PLATE
fem def init;r;master;example reg $PLATE
fem def contact;r;contact;example
fem def solv;r;LU;example reg $BREAST
# load steps
for ($i=1;$i<=$load_step;$i++){
fem change nodes translate by $load_step_size,0,0 node MASTER_MOVE reg $PLATE
$name=sprintf("plate_LoadStep_%02d",$i);
# fem export node;$name as Pad_LoadStep reg $PLATE
# fem export elem;$name as Pad_LoadStep reg $PLATE
if($i>=$RUN_LOAD){
$j=1;
$CONVERGED=0;
while ($CONVERGED==0){
######################################Projection
# Update XP from YP
fem up geom from sol to 1..3 reg $BREAST
# Set contact Xi points on specified faces
fem def xi;c contact_point faces BREAST_FACE points $CONT_PTS_1 by $CONT_PTS_2 boundary # define frictionless contact between breast and compression plates
fem def xi;c contact_tied faces RIB_FACE points $CONT_PTS_1 by $CONT_PTS_2 add # define tied contact between breast and rib
# Define data at xi positions
fem def data;c from_xi
# Update slave info
fem update data field to slave
# Project onto target face
fem def xi;c closest_face search_start 1000 faces PLATE_FACE,RIB_MASTER_FACE reg $PLATE
# Store projection gap in data fields
fem update data field from gap reg $PLATE
# Update master info
fem update data field to master
# Place initial geometry YP(ny,3) into XP
FEM up geom from sol YP_index 3 to 1..3 reg $BREAST
#####################################Mechanics problem
# Solve finite elasticity/contact problem
fem solve error 1.0D-08 iterate 150 reg 1
print " \n"
print " ==================== \n"
print " Load step = $i \n"
print " ==================== \n"
print " \n"
$j++;
}#CONVERGED
# SAVE EACH LOAD STEP SOLUTIONS
if($RUN_LOAD==1){
$name=sprintf("Compression.Load.%02d", $i);
} else {
$name=sprintf("Compression.Test.Load.%02d", $i);
}
fem def init;w;$name reg $BREAST
fem def data;w;$name contact
# Export final deformed meshes for current load step
# fem export node;$name as $name field reg $BREAST
# fem export elem;$name as $name field reg $BREAST
}
}