timplement bond-parallel interaction - slidergrid - grid of elastic sliders on a frictional surface
(HTM) git clone git://src.adamsgaard.dk/slidergrid
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) README
(DIR) LICENSE
---
(DIR) commit de268460b46629fdcfcfd155a1699fd3bd69b5c0
(DIR) parent 73fcb298a3755f169f631813de24a4e2904510d3
(HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
Date: Wed, 16 Mar 2016 11:41:29 -0700
implement bond-parallel interaction
Diffstat:
M slider.c | 78 +++++++++++++++++++++++++------
M slider.h | 5 +++--
M vector_math.h | 2 +-
3 files changed, 67 insertions(+), 18 deletions(-)
---
(DIR) diff --git a/slider.c b/slider.c
t@@ -62,7 +62,8 @@ void update_kinematics(slider s, Float dt, long int iteration)
}
// Find the relative displacement and velocity between two sliders
-void slider_displacement(slider s1, const slider s2, const Float dt, const int i)
+void slider_displacement(slider s1, const slider s2,
+ const int idx_neighbor, const int iteration)
{
// vector pointing from neighbor (s2) position to this slider position (s1)
t@@ -78,52 +79,99 @@ void slider_displacement(slider s1, const slider s2, const Float dt, const int i
const Float3 vel_linear = subtract_float3(s1.vel, s2.vel);
// relative contact interface velocity with rolling
- const Float3 vel = add_float3(vel_linear,
+ /*const Float3 vel = add_float3(vel_linear,
add_float3(
multiply_scalar_float3(dist_norm/2.,
cross_float3(dist, s1.angvel)),
multiply_scalar_float3(dist_norm/2.,
- cross_float3(dist, s2.angvel))));
+ cross_float3(dist, s2.angvel))));*/
// Normal component of the relative contact interface velocity
- const Float vel_n = multiply_scalar_float3(-1.0, dot_float3(vel_linear, n));
+ const Float vel_n = -1.0*dot_float3(vel_linear, n);
// Tangential component of the relative contact interface velocity
// Hinrichsen and Wolf 2004, eq. 13.9
- const Float3 = subtract_float3(vel, multiply_float3(n, dot_float3(n, vel)));
+ //const Float3 vel_t =
+ //subtract_float3(vel, multiply_float3_scalar(n, dot_float3(n, vel)));
// read previous inter-slider distance vector
- const Float3 dist0 = s1.neighbor_distance[i];
+ Float3 dist0;
+ if (iteration == 0)
+ dist0 = make_float3(0., 0., 0.);
+ else
+ dist0 = s1.neighbor_distance[idx_neighbor];
// increment in inter-slider distance
const Float3 ddist = subtract_float3(dist, dist0);
// save current inter-slider distance
- s1.neighbor_distance[i] = dist;
+ s1.neighbor_distance[idx_neighbor] = dist;
// total relative displacement in inter-slider distance
- s1.neighbor_relative_distance_displacement[i] =
- add_float3(s1.neighbor_relative_distance_displacement[i], ddist);
+ if (iteration == 0)
+ s1.neighbor_relative_distance_displacement[idx_neighbor] = ddist;
+ else
+ s1.neighbor_relative_distance_displacement[idx_neighbor] =
+ add_float3(s1.neighbor_relative_distance_displacement[idx_neighbor],
+ ddist);
// total relative displacement in inter-slider distance
- s1.neighbor_relative_distance_velocity[i] = vel_n;
-
+ s1.neighbor_relative_distance_velocity[idx_neighbor] =
+ multiply_float3_scalar(n, vel_n);
}
// Resolve bond mechanics between a slider and one of its neighbors based on the
// incremental linear-elastic model by Potyondy and Cundall 2004
+void slider_interaction(slider s1, const slider s2, const int idx_neighbor)
+{
+
+ // determine the bond stiffness from the harmonic mean.
+ // differs from Potyondy & Stack 2004, eq. 6.
+ const Float spring_stiffness =
+ 2.*s1.spring_stiffness*s2.spring_stiffness/
+ (s1.spring_stiffness + s2.spring_stiffness);
+
+
+ // bond-parallel elasticity
+ const Float3 f_n_elastic =
+ multiply_scalar_float3(bond_normal_stiffness,
+ s1.neighbor_relative_distance_displacement[idx_neighbor]);
+
+ // bond-parallel viscosity
+ const Float3 f_n_viscous =
+ multiply_scalar_float3(bond_normal_viscosity,
+ s1.neighbor_relative_distance_displacement[idx_neighbor]);
+
+ // bond-parallel force, counteracts tension and compression
+ const Float3 f_n = multiply_float3(f_n_elastic, f_n_viscous);
+
+ // add bond-parallel force to sum of forces on slider
+ s1.force = add_float3(s1.force, f_n);
+}
// Resolve interaction between a slider and all of its neighbors
void slider_neighbor_interaction(
slider s,
const slider* sliders,
- const int N)
+ const int N,
+ const int iteration)
{
- int i;
- for (i=0; i<MAX_NEIGHBORS; i++) {
- if (s.neighbors[i] != -1) {
+ int idx_neighbor;
+ for (idx_neighbor=0; idx_neighbor<MAX_NEIGHBORS; idx_neighbor++) {
+
+ // reset sum of forces on slider
+ s.force = make_float3(0., 0., 0.);
+
+ if (s.neighbors[idx_neighbor] != -1) {
+
+ slider_displacement(
+ s, sliders[s.neighbors[idx_neighbor]],
+ idx_neighbor, iteration);
+
+ slider_interaction(
+ s, sliders[s.neighbors[idx_neighbor]], idx_neighbor);
}
(DIR) diff --git a/slider.h b/slider.h
t@@ -30,7 +30,8 @@ typedef struct {
Float moment_of_inertia;
// inter-slider contact model parameters
- Float spring_stiffness; // elastic compressibility [N/m]
+ Float bond_normal_stiffness; // Hookean elastic stiffness [N/m]
+ Float bond_normal_viscosity; // viscosity [N/(m*s)]
// slider-surface contact model parameters
Float coulomb_friction; // Coulomb frictional value [-]
t@@ -43,7 +44,7 @@ typedef struct {
// relative spatial movement between this slider and its neighbors
Float3 neighbor_distance[MAX_NEIGHBORS];
Float3 neighbor_relative_distance_displacement[MAX_NEIGHBORS];
- Float neighbor_relative_distance_velocity[MAX_NEIGHBORS];
+ Float3 neighbor_relative_distance_velocity[MAX_NEIGHBORS];
} slider;
(DIR) diff --git a/vector_math.h b/vector_math.h
t@@ -13,7 +13,7 @@ Float3 subtract_float3(Float3 v1, Float3 v2);
Float3 multiply_float3(Float3 v1, Float3 v2);
Float3 divide_float3(Float3 v1, Float3 v2);
Float3 cross_float3(Float3 v1, Float3 v2);
-Float dot_float3(Float3 v1, Float3 v2)
+Float dot_float3(Float3 v1, Float3 v2);
// vector-scalar operations
Float3 multiply_float3_scalar(Float3 v, Float s);