[Bf-blender-cvs] [6a535b7] temp_hair_flow: Gradient calculation and normalization from the pressure field.

[Lukas Tönne]

Commit: 6a535b7b184357e33e7ad4ebfb5e0a524e7ebff6
Author: Lukas Tönne
Date:   Thu Jan 8 10:45:20 2015 +0100
Branches: temp_hair_flow
https://developer.blender.org/rB6a535b7b184357e33e7ad4ebfb5e0a524e7ebff6

Gradient calculation and normalization from the pressure field.

===================================================================

M	source/blender/physics/intern/grid.cpp
M	source/blender/physics/intern/grid.h
M	source/blender/physics/intern/hair_flow.cpp

===================================================================

diff --git a/source/blender/physics/intern/grid.cpp b/source/blender/physics/intern/grid.cpp
index ff188f0..0fc8917 100644
--- a/source/blender/physics/intern/grid.cpp
+++ b/source/blender/physics/intern/grid.cpp
@@ -66,19 +66,14 @@ void Grid::resize(float _cellsize, const float _offset[3], const int _res[3])
 	copy_v3_v3(offset, _offset);
 	copy_v3_v3_int(res, _res);
 	num_cells = _res[0] * _res[1] * _res[2];
-	
-	divergence.resize(num_cells);
-	pressure.resize(num_cells);
 }
 
 void Grid::init()
 {
-	divergence.setZero();
 }
 
 void Grid::clear()
 {
-	pressure.setZero();
 }
 
 int Grid::set_inner_cells(GridHash<bool> &bounds, GridHash<float3> &normal, Object *ob) const
@@ -224,6 +219,54 @@ void Grid::calc_divergence(GridHash<float> &divergence, const GridHash<bool> &so
 	}
 }
 
+/* Calculate velocity = grad(p) */
+void Grid::calc_gradient(GridHash<float3> &velocity, const GridHash<float> &pressure) const
+{
+	const float inv_flowfac = 1.0f / cellsize;
+	
+	velocity.clear();
+	
+	for (int z = 0; z < res[2]; ++z) {
+		for (int y = 0; y < res[1]; ++y) {
+			for (int x = 0; x < res[0]; ++x) {
+				bool is_margin = !(x > 0 && x < res[0]-1 && y > 0 && y < res[1]-1 && z > 0 && z < res[2]-1);
+				if (is_margin)
+					continue;
+				
+				/*const float *p  = pressure.get(x, y, z);*/
+				const float *pl = pressure.get(x-1, y, z);
+				const float *pr = pressure.get(x+1, y, z);
+				const float *pb = pressure.get(x, y-1, z);
+				const float *pt = pressure.get(x, y+1, z);
+				const float *pd = pressure.get(x, y, z-1);
+				const float *pu = pressure.get(x, y, z+1);
+				
+				/* finite difference estimate of pressure gradient */
+				float dvel[3];
+				dvel[0] = *pr - *pl;
+				dvel[1] = *pt - *pb;
+				dvel[2] = *pu - *pd;
+				mul_v3_fl(dvel, -0.5f * inv_flowfac);
+				
+				velocity.add(x, y, z) = float3(dvel);
+			}
+		}
+	}
+}
+
+void Grid::normalize(GridHash<float3> &velocity) const
+{
+	for (int z = 0; z < res[2]; ++z) {
+		for (int y = 0; y < res[1]; ++y) {
+			for (int x = 0; x < res[0]; ++x) {
+				float3 *v = velocity.get(x, y, z);
+				if (v)
+					normalize_v3(*v);
+			}
+		}
+	}
+}
+
 /* Main Poisson equation system:
  * This is derived from the discretezation of the Poisson equation
  *   div(grad(p)) = div(v)
@@ -234,7 +277,7 @@ void Grid::calc_divergence(GridHash<float> &divergence, const GridHash<bool> &so
  * For a good overview of eulerian fluid sim methods, see
  * http://www.proxyarch.com/util/techpapers/papers/Fluid%20flow%20for%20the%20rest%20of%20us.pdf
  */
-void Grid::solve_pressure(GridHash<float> &pressure, const GridHash<float> &divergence)
+void Grid::solve_pressure(GridHash<float> &pressure, const GridHash<float> &divergence) const
 {
 	int stride[3] = { 1, res[0], res[0] * res[1] };
 	
diff --git a/source/blender/physics/intern/grid.h b/source/blender/physics/intern/grid.h
index 206225f..0b39f02 100644
--- a/source/blender/physics/intern/grid.h
+++ b/source/blender/physics/intern/grid.h
@@ -140,18 +140,18 @@ typedef struct Grid {
 	void clear();
 	
 	int set_inner_cells(GridHash<bool> &bounds, GridHash<float3> &normal, struct Object *ob) const;
+	
 	void calc_divergence(GridHash<float> &divergence, const GridHash<bool> &source, const GridHash<float3> &source_normal) const;
+	void calc_gradient(GridHash<float3> &velocity, const GridHash<float> &pressure) const;
+	void normalize(GridHash<float3> &velocity) const;
 	
-	void solve_pressure(GridHash<float> &pressure, const GridHash<float> &divergence);
+	void solve_pressure(GridHash<float> &pressure, const GridHash<float> &divergence) const;
 	
 	float cellsize, inv_cellsize;
 	float offset[3];
 	int res[3];
 	int num_cells;
 	
-	lVector divergence;
-	lVector pressure;
-	
 	struct SimDebugData *debug_data;
 	float debug1, debug2;
 	int debug3, debug4;
diff --git a/source/blender/physics/intern/hair_flow.cpp b/source/blender/physics/intern/hair_flow.cpp
index e0c0696..dcf7ca7 100644
--- a/source/blender/physics/intern/hair_flow.cpp
+++ b/source/blender/physics/intern/hair_flow.cpp
@@ -148,6 +148,11 @@ HairFlowData *BPH_strands_solve_hair_flow(Scene *scene, Object *ob, float max_le
 	pressure.resize(data->grid.res);
 	data->grid.solve_pressure(pressure, divergence);
 	
+	GridHash<float3> velocity;
+	velocity.resize(data->grid.res);
+	data->grid.calc_gradient(velocity, pressure);
+	data->grid.normalize(velocity);
+	
 	{
 		float col0[3] = {0.0, 0.0, 0.0};
 		float colp[3] = {0.0, 1.0, 1.0};
@@ -185,6 +190,8 @@ HairFlowData *BPH_strands_solve_hair_flow(Scene *scene, Object *ob, float max_le
 					}
 					if (fac > 0.05f)
 						BKE_sim_debug_data_add_circle(debug_data, vec, 0.02f, col[0], col[1], col[2], "hair_flow", 5522, x, y, z);
+					
+					BKE_sim_debug_data_add_vector(debug_data, vec, *velocity.get(x, y, z), 1,1,0, "hair_flow", 957, x, y, z);
 				}
 			}
 		}