// capillary.fe

// Demonstration of curved constraint - capillary meniscus in 
// vertical tube.

parameter radius = 1  // of tube
parameter height = 4  // of tube
parameter angle = 45  // interior contact angle, degrees

#define WALLT (-cos(angle*pi/180))

constraint 1  // the tube wall, for contact line
formula: x^2 + y^2 = radius^2
energy:
e1:  -WALLT*z*y/radius
e2:   WALLT*z*x/radius
e3:   0

constraint 2  // the tube wall, for display 
formula: x^2 + y^2 = radius^2

vertices
// the meniscus
1  radius*cos(0*pi/3) radius*sin(0*pi/3) height/2 constraint 1
2  radius*cos(1*pi/3) radius*sin(1*pi/3) height/2 constraint 1
3  radius*cos(2*pi/3) radius*sin(2*pi/3) height/2 constraint 1
4  radius*cos(3*pi/3) radius*sin(3*pi/3) height/2 constraint 1
5  radius*cos(4*pi/3) radius*sin(4*pi/3) height/2 constraint 1
6  radius*cos(5*pi/3) radius*sin(5*pi/3) height/2 constraint 1
// bottom of tube
11  radius*cos(0*pi/3) radius*sin(0*pi/3) 0 constraint 2
12  radius*cos(1*pi/3) radius*sin(1*pi/3) 0 constraint 2
13  radius*cos(2*pi/3) radius*sin(2*pi/3) 0 constraint 2
14  radius*cos(3*pi/3) radius*sin(3*pi/3) 0 constraint 2
15  radius*cos(4*pi/3) radius*sin(4*pi/3) 0 constraint 2
16  radius*cos(5*pi/3) radius*sin(5*pi/3) 0 constraint 2
// top of tube
21  radius*cos(0*pi/3) radius*sin(0*pi/3) height constraint 2
22  radius*cos(1*pi/3) radius*sin(1*pi/3) height constraint 2
23  radius*cos(2*pi/3) radius*sin(2*pi/3) height constraint 2
24  radius*cos(3*pi/3) radius*sin(3*pi/3) height constraint 2
25  radius*cos(4*pi/3) radius*sin(4*pi/3) height constraint 2
26  radius*cos(5*pi/3) radius*sin(5*pi/3) height constraint 2

edges
// the meniscus
1   1 2 constraint 1
2   2 3 constraint 1
3   3 4 constraint 1
4   4 5 constraint 1
5   5 6 constraint 1
6   6 1 constraint 1
// the bottom of the tube
7  11 12  constraint 2 fixed bare no_refine
8  12 13  constraint 2 fixed bare no_refine
9  13 14  constraint 2  fixed bare no_refine
10  14 15  constraint 2 fixed bare no_refine
11  15 16  constraint 2 fixed bare no_refine
12  16 11  constraint 2 fixed bare no_refine
// top of the tube
13  21 22  constraint 2 fixed bare no_refine
14  22 23  constraint 2 fixed bare no_refine
15  23 24  constraint 2 fixed bare no_refine
16  24 25  constraint 2 fixed bare no_refine
17  25 26  constraint 2 fixed bare no_refine
18  26 21  constraint 2 fixed bare no_refine
// Verticals sides of tube
19  11 21  constraint 2 fixed bare no_refine
20  12 22  constraint 2 fixed bare no_refine
21  13 23  constraint 2 fixed bare no_refine
22  14 24  constraint 2 fixed bare no_refine
23  15 25  constraint 2 fixed bare no_refine
24  16 26  constraint 2 fixed bare no_refine

faces
1   1 2 3 4 5 6

bodies 
1    1  volume pi*radius^2*height/2

read

// get top and bottom round-looking
refine edge where original >= 7 and original <= 18
refine edge where original >= 7 and original <= 18
refine edge where original >= 7 and original <= 18