[Bf-blender-cvs] [754e366] hair_system: Improved damping calculation.
Lukas Tönne
noreply at git.blender.org
Sun Aug 17 20:32:46 CEST 2014
Commit: 754e3664daea82be7b1119a121052e582a9b0a77
Author: Lukas Tönne
Date: Sun Aug 17 20:30:52 2014 +0200
Branches: hair_system
https://developer.blender.org/rB754e3664daea82be7b1119a121052e582a9b0a77
Improved damping calculation.
Now actually uses the damping substeps setting. Damping happens as an
additional step _after_ the main force calculations. Note that it
updates the next.velocity based on the next.velocity itself, unlike the
actual semi-implicit Euler step, which works with the previous velocity.
This makes damping much more accurate by letting it work with the
changed velocities right away before they introduce too much energy into
the system.
===================================================================
M source/blender/hair/intern/HAIR_smoothing.h
M source/blender/hair/intern/HAIR_solver.cpp
===================================================================
diff --git a/source/blender/hair/intern/HAIR_smoothing.h b/source/blender/hair/intern/HAIR_smoothing.h
index 3959b46..30724eb 100644
--- a/source/blender/hair/intern/HAIR_smoothing.h
+++ b/source/blender/hair/intern/HAIR_smoothing.h
@@ -118,6 +118,11 @@ protected:
HAIR_CXX_CLASS_ALLOC(SmoothingIterator)
};
+/* XXX optimizing the frame iterator could be quite rewarding in terms of performance ...
+ * See "Parallel Transport Approach to Curve Framing" (Hanson et al. 1995)
+ * ftp://cgi.soic.indiana.edu/pub/techreports/TR425.pdf
+ */
+
template <typename WalkerT>
struct FrameIterator {
FrameIterator()
diff --git a/source/blender/hair/intern/HAIR_solver.cpp b/source/blender/hair/intern/HAIR_solver.cpp
index defaf17..b90a148 100644
--- a/source/blender/hair/intern/HAIR_solver.cpp
+++ b/source/blender/hair/intern/HAIR_solver.cpp
@@ -227,12 +227,13 @@ static float3 calc_stretch_force(const HairParams ¶ms, const Point *point0,
return params.stretch_stiffness * (length - rest_length) * dir;
}
-static float3 calc_stretch_damping(const HairParams ¶ms, const Point *point0, const Point *point1, float time)
+static float3 calc_stretch_damping(const HairParams ¶ms, const Point *point0, const Point *point1)
{
float3 dir;
float3 edge = point1->cur.co - point0->cur.co;
normalize_v3_v3(dir, edge);
- float3 dvel = point1->cur.vel - point0->cur.vel;
+ /* NB: for damping, use the updated next.vel velocities! */
+ float3 dvel = point1->next.vel - point0->next.vel;
return params.stretch_damping * dot_v3v3(dvel, dir) * dir;
}
@@ -257,12 +258,13 @@ static float3 calc_bend_force(const HairParams ¶ms, const Point *point0, con
return params.bend_stiffness * (edge - target);
}
-static float3 calc_bend_damping(const HairParams ¶ms, const Point *point0, const Point *point1, const Frame &frame, float time)
+static float3 calc_bend_damping(const HairParams ¶ms, const Point *point0, const Point *point1)
{
float3 dir;
float3 edge = point1->cur.co - point0->cur.co;
normalize_v3_v3(dir, edge);
- float3 dvel = point1->cur.vel - point0->cur.vel;
+ /* NB: for damping, use the updated next.vel velocities! */
+ float3 dvel = point1->next.vel - point0->next.vel;
return params.bend_damping * (dvel - dot_v3v3(dvel, dir) * dir);
}
@@ -274,7 +276,7 @@ struct SolverTaskData {
int startpoint, totpoints;
};
-static void do_internal_forces(const HairParams ¶ms, const SolverForces &forces, float time, float timestep, const Point *point0, const Point *point1, const Frame &frame,
+static void accum_internal_forces(const HairParams ¶ms, const SolverForces &forces, float time, float timestep, const Point *point0, const Point *point1, const Frame &frame,
float3 &force0, float3 &force1)
{
float3 stretch, bend;
@@ -292,7 +294,7 @@ static void do_internal_forces(const HairParams ¶ms, const SolverForces &for
force1 = - stretch - bend;
}
-static void do_external_forces(const HairParams ¶ms, const SolverForces &forces, float time, float timestep, const Point *point0, const Point *point1, const Frame &frame,
+static void accum_external_forces(const HairParams ¶ms, const SolverForces &forces, float time, float timestep, const Point *point0, const Point *point1, const Frame &frame,
float3 &force)
{
float3 acc = float3(0.0f, 0.0f, 0.0f);
@@ -304,29 +306,6 @@ static void do_external_forces(const HairParams ¶ms, const SolverForces &for
force = acc;
}
-static void do_damping(const HairParams ¶ms, const SolverForces &forces, float time, float timestep, const Point *point0, const Point *point1, const Frame &frame,
- float3 &force0, float3 &force1)
-{
- const int totsteps = params.substeps_damping;
- float dt = timestep / (float)totsteps;
- float3 stretch, bend;
-
- force0 = float3(0.0f, 0.0f, 0.0f);
- force1 = float3(0.0f, 0.0f, 0.0f);
-
- if (point1) {
- for (int step = 0; step < totsteps; ++step) {
- stretch = calc_stretch_damping(params, point0, point1, time);
- bend = calc_bend_damping(params, point0, point1, frame, time);
-
- force0 = force0 + stretch + bend;
- force1 = force1 - stretch - bend;
-
- time += dt;
- }
- }
-}
-
static void do_collision(const HairParams ¶ms, const SolverForces &forces, float time, float timestep, float restitution_scale, const SolverTaskData &data, const PointContactCache &contacts)
{
int start = data.startpoint;
@@ -389,7 +368,7 @@ static void calc_forces(const HairParams ¶ms, const SolverForces &forces, fl
Frame rest_frame(root_nor, root_tan, cross_v3_v3(root_nor, root_tan));
FrameIterator<SolverDataLocWalker> frame_iter(SolverDataLocWalker(curve), curve->avg_rest_length, params.bend_smoothing, rest_frame);
- float3 intern_force, intern_force_next, extern_force, damping, damping_next;
+ float3 intern_force, intern_force_next, extern_force;
float3 acc_prev; /* reactio from previous point */
/* Root point animation */
@@ -397,12 +376,11 @@ static void calc_forces(const HairParams ¶ms, const SolverForces &forces, fl
int k = 0;
if (numpoints > 0) {
Point *point_next = k < numpoints-1 ? point+1 : NULL;
- do_internal_forces(params, forces, time, timestep, point, point_next, frame_iter.frame(), intern_force, intern_force_next);
- do_external_forces(params, forces, time, timestep, point, point_next, frame_iter.frame(), extern_force);
- do_damping(params, forces, time, timestep, point, point_next, frame_iter.frame(), damping, damping_next);
+ accum_internal_forces(params, forces, time, timestep, point, point_next, frame_iter.frame(), intern_force, intern_force_next);
+ accum_external_forces(params, forces, time, timestep, point, point_next, frame_iter.frame(), extern_force);
frame_iter.next();
- acc_prev = intern_force_next + damping_next;
+ acc_prev = intern_force_next;
++k;
++point;
}
@@ -410,14 +388,13 @@ static void calc_forces(const HairParams ¶ms, const SolverForces &forces, fl
/* Integrate free points */
for (; k < numpoints; ++k, ++point) {
Point *point_next = k < numpoints-1 ? point+1 : NULL;
- do_internal_forces(params, forces, time, timestep, point, point_next, frame_iter.frame(), intern_force, intern_force_next);
- do_external_forces(params, forces, time, timestep, point, point_next, frame_iter.frame(), extern_force);
- do_damping(params, forces, time, timestep, point, point_next, frame_iter.frame(), damping, damping_next);
+ accum_internal_forces(params, forces, time, timestep, point, point_next, frame_iter.frame(), intern_force, intern_force_next);
+ accum_external_forces(params, forces, time, timestep, point, point_next, frame_iter.frame(), extern_force);
- point->force_accum = point->force_accum + intern_force + extern_force + damping + acc_prev;
+ point->force_accum = point->force_accum + intern_force + extern_force + acc_prev;
frame_iter.next();
- acc_prev = intern_force_next + damping_next;
+ acc_prev = intern_force_next;
}
}
@@ -425,40 +402,52 @@ static void calc_forces(const HairParams ¶ms, const SolverForces &forces, fl
do_collision(params, forces, time, timestep, restitution_scale, data, contacts);
}
-static void calc_velocity_step(const SolverTaskData &data, float time, float timestep, float t0, float t1)
+static void damping_forces(const HairParams ¶ms, const Point *point0, const Point *point1,
+ float3 &force0, float3 &force1)
+{
+ float3 stretch, bend;
+
+ stretch = calc_stretch_damping(params, point0, point1);
+ bend = calc_bend_damping(params, point0, point1);
+
+ force0 = stretch + bend;
+ force1 = -(stretch + bend);
+}
+
+static void apply_acceleration(const SolverTaskData &data, float time, float timestep, float t0, float t1)
{
Curve *curve = data.curves;
for (int i = 0; i < data.totcurves; ++i, ++curve) {
- int numpoints = curve->totpoints;
+ int totpoints = curve->totpoints;
/* note: roots are evaluated at the end of the timestep: time + timestep
- * so the hair points align perfectly with them
- */
+ * so the hair points align perfectly with them
+ */
float3 root_vel = get_root_velocity(t0, t1, time + timestep, curve);
/* Root point animation */
Point *point = curve->points;
int k = 0;
- if (numpoints > 0) {
- point->next.vel = root_vel;
+ if (totpoints > 0) {
+ curve->points->next.vel = root_vel;
++k;
++point;
}
/* Integrate free points */
- for (; k < numpoints; ++k, ++point) {
+ for (; k < totpoints; ++k, ++point) {
/* semi-implicit Euler step */
point->next.vel = point->cur.vel + point->force_accum * timestep;
}
}
}
-static void calc_position_step(const SolverTaskData &data, float time, float timestep, float t0, float t1)
+static void apply_velocity(const SolverTaskData &data, float time, float timestep, float t0, float t1)
{
Curve *curve = data.curves;
for (int i = 0; i < data.totcurves; ++i, ++curve) {
- int numpoints = curve->totpoints;
+ int totpoints = curve->totpoints;
/* note: roots are evaluated at the end of the timestep: time + timestep
* so the hair points align perfectly with them
@@ -468,7 +457,7 @@ static void calc_position_step(const SolverTaskData &data, float time, float tim
/* Root point animation */
Point *point = curve->points;
int k = 0;
- if (numpoints > 0) {
+ if (totpoints > 0) {
point->next.co = root_co;
++k;
@@ -476,13 +465,50 @@ static void calc_position_step(const SolverTaskData &data, float time, float tim
}
/* Integrate free points */
- for (; k < numpoints; ++k, ++point) {
+ for (; k < totpoints; ++k, ++point) {
/* semi-implicit Euler step */
point->next.co = point->cur.co + point->next.vel * timestep;
}
}
}
+static void do_damping(const HairParams ¶ms, float time, float timestep, const SolverTaskData &data)
+{
+ const int totsteps = params.substeps_damping;
+ float dt = timestep / (float)totsteps;
+
+ for (int step = 0; step < totsteps; ++step) {
+ Curve *curve = data.curves;
+ for (int i = 0; i < data.totcurves; ++i, ++curve) {
+ int totpoints = curve->totpoints;
@@ Diff output truncated at 10240 characters. @@
More information about the Bf-blender-cvs
mailing list