/* Evolver file for three spheres,
modeled fluid volume sitting in
interstice.
Bottom one with center at origin,
is fixed, the other rotate around
the center amounts +- theta/2
W. Craig Carter, 1996
*/
SYMMETRIC_CONTENT
PARAMETER phi = 30 // liquid wetting angle
PARAMETER theta = 0 // sphere opening angle
PARAMETER vol = 1.5 // liquid volume
// !!! NEED volconst = -2*pi
#define Rs 1 // sets units
#define gamsl (-cos(phi*pi/180)) // effective energy
/* a bunch of functions below, used for vector integrands */
#define ct2 cos( 60*pi/180 - theta*pi/360 )
#define st2 sin( 60*pi/180 - theta*pi/360 )
#define ct3 cos( 120*pi/180 + theta*pi/360 )
#define st3 sin( 120*pi/180 + theta*pi/360 )
/*radius 
2 of bottom sphere*/
#define r1sq (x) 2 + (y) 2 + z 2
/*upper right sphere */
#define r2sq (x - 2*Rs*ct2) 2 + (y - 2*Rs*st2) 2 + z 2
/*upper left sphere */
#define r3sq (x - 2*Rs*ct3) 2 + (y - 2*Rs*st3) 2 + z 2
#define w1num y
#define w1den (x 2 + z 2)*sqrt(r1sq)
#define w1 w1num/d1den
#define w2num x*ct2 + y*st2 - 2*Rs
#define w2den ((st2*x - ct2*y) 2 + z 2)*sqrt(r2sq)
#define w2 w2num/d2den
#define w3num x*ct3 + y*st3 - 2*Rs
#define w3den ((st3*x - ct3*y) 2 + z 2)*sqrt(r3sq)
#define w3 w3num/d3den
constraint 1 //bottom sphere
formula: r1sq = Rs 2
energy:
e1: gamsl*(Rs 2)*w1*z
e2: 0
e3: -gamsl*(Rs 2)*w1*x
content:
c1: ((Rs 3)/3)*w1*z
c2: 0
c3: -((Rs 3)/3)*w1*x
constraint 2 //northeast sphere theta = Pi/6
formula: r2sq = Rs 2
energy:
e1: (gamsl*(Rs 2)*w2*(z*st2))
e2: -(gamsl*(Rs 2)*w2*(z*ct2))
e3: -(gamsl*(Rs 2)*w2*(x*st2 - y*ct2))
content:
c1: -(-(((Rs 3)/3)*w2*(z*st2)) + z*(2*Rs)*st2/3)
c2: -((((Rs 3)/3)*w2*(z*ct2)) - z*(2*Rs)*ct2/3)
c3: -(((Rs 3)/3)*w2*(x*st2 - y*ct2))
constraint 3
formula: r3sq = Rs 2
energy:
e1: (gamsl*(Rs 2)*w3*(z*st3))
e2: -(gamsl*(Rs 2)*w3*(z*ct3))
e3: -(gamsl*(Rs 2)*w3*(x*st3 - y*ct3))
content:
c1: -(-(((Rs 3)/3)*w3*(z*st3)) + z*(2*Rs)*st3/3)
c2: -((((Rs 3)/3)*w3*(z*ct3)) - z*(2*Rs)*ct3/3)
c3: -(((Rs 3)/3)*w3*(x*st3 - y*ct3))
vertices
//:
//:
edges
//:
//:
faces
//:
//:
bodies
//:
//: volconst -2*pi
read
checksmall := {delete edges where length < 0.0001; 
delete facets where area < 0.005; u;
delete edges where length < 0.0001; u;g3}
edgemax := 0.4
echeck := {while max(edges,length) > edgemax do
{ refine edges where length > edgemax } ;
checksmall }
grind := { ig:= 1 ; change := 1.; olde := total_energy ;
while (ig < 10000 && abs(change) > .00001)
do {ig := ig+1;
if ig%10 == 1 then checksmall ;
g | "cat > null" ;
if ig%3 == 1 then
{ change := total_energy - olde ;
olde := total_energy;
printf "sum of 3 dE: %12.10f n", change}}}
goprint :=
{ printf "%f %f %20.12f %f "
,phi, theta, total_energy, change | "cat >> results";
foreach body do { printf "%f ", volume | "cat >> results" } ;
printf "%f ", total_area | "cat >> results" ;
printf "%g ", facet_count | "cat >> results" ;
printf "%f n", max(edges,length) | "cat >> results" }
dump_to_file :=
{filename := sprintf "dump.v_%f_p_%f_t_%f",vol,phi,theta;
dump filename}
inc := {echeck; D; D ; O ; o ;u; theta := theta + 1 ;
printf "%f", theta ;u;g5; biggrind ;
if theta%10 = 0 then dump_to_file; goprint }
run_this := {dump_to_file ; while (theta < 120) do {inc}}