Liquid droplets

The 3D model I propose to study in this project is, in a certain sense, very different from those that can be found in the vast online repositories under the category of mathematical objects: liquid droplets on a flat surface or suspended, numerically modeled by solving an ordinary differential equation and printed on a 3D printer. 

It's a physical manifestation of a simulation of a real physical phenomenon. 

And that's a very cool thing!

We have all seen examples of complex mathematical structures printed on a 3D printer. Almost all of these structures were generated or built without conceptual restrictions or morphological limitations inspired by the natural or physical world, i.e., they are models that are morphologically limited only by the conceptual restrictions of the artist. 

These droplets models are a very different kind of objects because  they are simultaneously a physical manifestation and a numerical simulation of a concrete physical phenomenon that exists.

Characteristics of Liquid Droplets

When a droplet of a fluid, which we can always assume to be a water droplet without loss of generality since they are the ones we are most accustomed to, rests on a horizontal surface or is suspended from a horizontal plane, with the rest of the droplet being free, the surface generated by the fluid boundary is evidently a surface of revolution around a vertical axis and intersects this axis at right angles. Let us take the axis of revolution as the z-axis and the point where it intersects the free surface as the origin. 

Obviously this topic has quite a few challenges, but lets look at the mathematics and physics of droplets! I'm using Francis Bashforth's text “An attempt to test the theories of capillary action by comparing the theoretical and measured forms of drops of fluid. With an explanation of the method of integration employed in constucting the tables which give the theoretical forms of such drops”. With this title it really deserves a read!

Bashforth's 1883 paper can be found on Internet Archive: https://archive.org/details/attempttest00bashrich

P.S.
If you're only interested on printing the model you can safely skip the math part. 

Lets math!

Let us consider the liquid in equilibrium under the action of the gravitational force and a uniform surface tension T. Let R and R' be the principal radii of curvature at a point on the liquid surface and p the pressure at that point, thus:

          1/R + 1/R' = p/T          (1)

The pressure p is given by p=gσz +C where g is the acceleration due to gravity, σ is the fluid density and z is the vertical coordinate.

Let x be the horizontal coordinate and z the vertical coordinate at any point on the meridian of the fluid surface, ρ the radius of curvature of the meridional section and ϕ the angle that the normal to the surface makes with the axis of rotation. Then

          R = ρ,

          R' = x/sin(ϕ).          (2)

Equation (1) becomes:  

          1/ρ + sin(ϕ)/x =(C+g σ z)/T.          (3)

In the limit ϕ = 0 (x = z = 0) both radii of curvature are equal and we have R = R′ = b where b is a real number and therefore C/T = 2/b. The equation governing the shape of the droplet can then be written as:

          1/ρ + sin(ϕ)/x = 2/b +g σ z/T.          (4)

Defining the quantities (x, z, ρ) →(x/b, z/b, ρ/b) and the constant β = gσb²/T equation (4) can be written as:

          1/ρ + sin(ϕ)/x = 2/b +β.          (5)

This shows how the shape of the droplet depends only on the parameter β and that the constant b defines the scaling parameter of the problem, rescaling all relevant quantities. 

The most relevant form of these last equations for their numerical solution can be constructed in terms of the meridional arc length s. The problem is thus formulated as a system of ordinary differential equations whose numerical solution allows us to determine the shape of a droplet.

The problem is thus formulated as a system of ordinary differential equations whose numerical solution allows us to determine the shape of a drop. The equations take the form

          dx/ds = cos(ϕ),
          dz/ds = sin(ϕ),
          dϕ/ds = 2 + βz - sin(ϕ)/x.          (6)

These set of ordinary differential equations should be solve numerically. 

I've used GNU/Octave and a simple script to generate the STL files directly, although this could be also done using only OpenSCAD. The code for solving the ODEs is given at the bottom of this text. Probably I'll do a version 2 of this project :D, stay tuned!

The next image shows the profiles of droplets on a horizontal plane for positive values of β and hanging droplets for negative values of β.

3D printed droplets!

There's no special requisites for printing these models. Layer height of 0.2mm, and 3 tops/bottoms, and 20% infill works great.

Here's some pics. 

Printed suspended droplets (β = −1, −0.5, −0.35).

Here's the OpenSCAD screen shot with a red ping-pong ball.

Printed drops on a plane (β = 1, 10, 100).

Here's the OpenSCAD screen shot with a red ping-pong ball.

All the drops have a volume of 33510 mm^2, the same volume as the typical ping-pong ball ;)

Here's photo  with a ping-pong ball for comparing. 

Conclusion

There are many examples of incredible 3D prints that have been created by mathematical generation but this set of 3D models is quite different in that it is simulating a physical phenomenon and transforming that simulation into a physical reality.
 

Have fun and enjoy your BIG drops!

Citation
Charters, T. (2025). Liquid droplets: a physical manifestation of a simulation of a real physical phenomenon. Zenodo. https://doi.org/10.5281/zenodo.17604747

Appendix

Here's the GNU/Octave code for solving the ODEs.

clear all
x0=[eps 0 0]

n=1000;
t=linspace (0,30,n);
choice="u";

if (choice=="u")
 y=lsode ("fu",x0,t);
 cte=1;
else
 y=lsode ("fl",x0,t);
 cte=-1;
end

figure(1)
#clf
hold on
for i=1:n-1
 if y(i,2)-y(i+1,2)<=0
   u(i)=y(i,1);
   v(i)=y(i,2);
   rho(i)=2+y(i,2)-sin(y(i,3))./y(i,1);
 else
   break
 end
end
zmin=max(v);
z=cte*(v-zmin);
plot(u,z)
plot(-u,z)

axis("equal")
grid on
xlabel("x")
ylabel("z")

And two auxiliary functions:

function xdot=fu(x,t)
 beta=-.5;
 xdot(1)=cos(x(3)); #x
 xdot(2)=sin(x(3)); #z
 xdot(3)=2+beta*x(2)-sin(x(3))/x(1); #phi
endfunction
 

function xdot=fl(x,t)
 beta=.2;
 xdot(1)=cos(x(3)); #x
 xdot(2)=sin(x(3)); #z
 xdot(3)=2+beta*x(2)-sin(x(3))/x(1); #phi
endfunction