Cloth: implement support for a hydrostatic pressure gradient.

When a fluid is put under influence of gravity or acceleration, it
forms an internal pressure gradient, which causes observable effects
like buoyancy. Since now cloth has support for simulating pressure
changes caused by fluid compression or expansion, it makes sense to
also support the effects of gravity.

This is intended for better simulation of objects filled or
surrounded by fluids, especially when constrained by collisions
or pinned vertices, and should result in more realistic shapes.

Obviously, this doesn't actually simulate fluid dynamics; instead
it is assumed that the fluid immediately adapts to changes in the
shape or acceleration of the object without friction or turbulence,
and instantly reaches a new static equilibrium.

Differential Revision: https://developer.blender.org/D6442

1 files changed, 45 insertions, 7 deletions

diff --git a/source/blender/physics/intern/implicit_blender.c b/source/blender/physics/intern/implicit_blender.c
index b1afaa52ef2..063c224f158 100644
--- a/source/blender/physics/intern/implicit_blender.c
+++ b/source/blender/physics/intern/implicit_blender.c
@@ -1246,6 +1246,11 @@ void BPH_mass_spring_get_position(struct Implicit_Data *data, int index, float x
   root_to_world_v3(data, index, x, data->X[index]);
 }
 
+void BPH_mass_spring_get_velocity(struct Implicit_Data *data, int index, float v[3])
+{
+  root_to_world_v3(data, index, v, data->V[index]);
+}
+
 void BPH_mass_spring_get_new_position(struct Implicit_Data *data, int index, float x[3])
 {
   root_to_world_v3(data, index, x, data->Xnew[index]);
@@ -1488,21 +1493,54 @@ float BPH_tri_tetra_volume_signed_6x(Implicit_Data *data, int v1, int v2, int v3
   return volume_tri_tetrahedron_signed_v3_6x(data->X[v1], data->X[v2], data->X[v3]);
 }
 
-void BPH_mass_spring_force_pressure(
-    Implicit_Data *data, int v1, int v2, int v3, float pressure_difference, float weights[3])
+float BPH_tri_area(struct Implicit_Data *data, int v1, int v2, int v3)
+{
+  float nor[3];
+
+  return calc_nor_area_tri(nor, data->X[v1], data->X[v2], data->X[v3]);
+}
+
+void BPH_mass_spring_force_pressure(Implicit_Data *data,
+                                    int v1,
+                                    int v2,
+                                    int v3,
+                                    float common_pressure,
+                                    const float *vertex_pressure,
+                                    const float weights[3])
 {
   float nor[3], area;
-  float factor;
+  float factor, base_force;
+  float force[3];
 
   /* calculate face normal and area */
   area = calc_nor_area_tri(nor, data->X[v1], data->X[v2], data->X[v3]);
   /* The force is calculated and split up evenly for each of the three face verts */
-  factor = pressure_difference * area / 3.0f;
+  factor = area / 3.0f;
+  base_force = common_pressure * factor;
+
+  /* Compute per-vertex force values from local pressures.
+   * From integrating the pressure over the triangle and deriving
+   * equivalent vertex forces, it follows that:
+   *
+   * force[idx] = (sum(pressure) + pressure[idx]) * area / 12
+   *
+   * Effectively, 1/4 of the pressure acts just on its vertex,
+   * while 3/4 is split evenly over all three.
+   */
+  if (vertex_pressure) {
+    copy_v3_fl3(force, vertex_pressure[v1], vertex_pressure[v2], vertex_pressure[v3]);
+    mul_v3_fl(force, factor / 4.0f);
+
+    base_force += force[0] + force[1] + force[2];
+  }
+  else {
+    zero_v3(force);
+  }
 
   /* add pressure to each of the face verts */
-  madd_v3_v3fl(data->F[v1], nor, factor * weights[0]);
-  madd_v3_v3fl(data->F[v2], nor, factor * weights[1]);
-  madd_v3_v3fl(data->F[v3], nor, factor * weights[2]);
+  madd_v3_v3fl(data->F[v1], nor, (base_force + force[0]) * weights[0]);
+  madd_v3_v3fl(data->F[v2], nor, (base_force + force[1]) * weights[1]);
+  madd_v3_v3fl(data->F[v3], nor, (base_force + force[2]) * weights[2]);
 }
 
 static void edge_wind_vertex(const float dir[3],