Page 320 - The Combined Finite-Discrete Element Method
P. 320
SORTING CONTACT DETECTION ALGORITHM 303
.........continued from previous listing.............
forco=R0; uforc=R0; vforc=R0; /* force and center of force */
for(i=1;i<(nbpoin-1);i++)
{ penetr=penetb[0]+penetb[i]+penetb[i+1];
if(penetr>EPSILON)
{ force=((ub[i]-ub[0])*(vb[i+1]-vb[0])-
(vb[i]-vb[0])*(ub[i+1]-ub[0]))*penetr*penalty;
fact0=(RP5*penetb[0]+
RP25*(penetb[i]+penetb[i+1])) / penetr;
facti=(RP5*penetb[i]+
RP25*(penetb[0]+penetb[i+1])) / penetr;
fact1=R1-fact0-facti;
if(ABS(force+forco)>EPSILON)
{ uforc=(forco*uforc+force*
(fact0*ub[0]+facti*ub[i]+fact1*ub[i+1])) / (forco+force);
vforc=(forco*vforc+force*
(fact0*vb[0]+facti*vb[i]+fact1*vb[i+1])) / (forco+force);
forco=forco+force;
}}}
Listing 10.17 Integration of contact force over the B-polygon.
.........continued from previous listing.............
/* resultant at C-points */
for(i=0;i<4;i++)
{ fc[i]=R0; ft[i]=R0;
}
tmp=((uc[1]-uc[0])*(vc[2]-vc[0])-(vc[1]-vc[0])*(uc[2]-uc[0]));
for(i=0;i<3;i++)
{ j=i+1; if(j>2)j=0; k=j+1; if(k>2)k=0;
fc[k]=forco*
(((uc[j]-uc[i])*(vforc-vc[i])-(vc[j]-vc[i])*(uforc-uc[i])) / tmp);
}
Listing 10.18 Calculating resultant force at C-points.
dimensional arrays forming a part of a single two-dimensional array and using dynamic
memory allocation (see Section 10.2).
10.5 SORTING CONTACT DETECTION ALGORITHM
As explained in Chapter 3, the sorting contact detection algorithm is based on sorting
one-dimensional arrays. A general function for sorting such arrays is given in Listings
10.21 and 10.22. In Listing 10.21, local variables are given, while in Listing 10.22, the
code performing actual sorting is given. The sorting follows procedure, which is described