tadd simulation parameter value check - 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 abf7d095f654ff0108832d8055a3ce72bf333920
(DIR) parent d71a34ad2042f7d67856a76cc130c204a69a3550
(HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
Date: Wed, 16 Mar 2016 13:01:36 -0700
add simulation parameter value check
Diffstat:
M main.c | 12 +++++++++---
M simulation.c | 63 +++++++++++++++++++++++++++++--
M simulation.h | 9 +++++----
M slider.c | 6 +++---
M slider.h | 4 ----
5 files changed, 77 insertions(+), 17 deletions(-)
---
(DIR) diff --git a/main.c b/main.c
t@@ -23,6 +23,8 @@ simulation setup_simulation()
sim.sliders[i].mass = 1.;
sim.sliders[i].moment_of_inertia = 1.;
}
+
+ return sim;
}
int main(int argc, char** argv)
t@@ -40,16 +42,20 @@ int main(int argc, char** argv)
find_and_bond_to_neighbors_n2(sim.sliders, sim.N, sim.bond_length_limit);
// determine time step length
- sim.dt = find_time_step(sim.sliders, sim.N);
+ sim.dt = find_time_step(sim.sliders, sim.N, 0.07);
// check simulation parameters for missing or wrong parameter values
- if (check_simulation(sim)) {
+ if (check_simulation_values(sim)) {
fprintf(stderr, "Error: Simulation configuration invalid.\n");
return EXIT_FAILURE;
}
// temporal loop
- for (sim.t=0.0; sim.t<=sim.t_end; sim.t+=sim.dt) {
+ sim.iteration = 0;
+ for (sim.time = 0.0;
+ sim.time <= sim.time_end;
+ sim.time += sim.dt) {
+
for (i=0; i<sim.N; i++) {
update_kinematics(sim.sliders[i], sim.dt, sim.iteration);
}
(DIR) diff --git a/simulation.c b/simulation.c
t@@ -1,16 +1,73 @@
+#include <stdio.h>
+#include <math.h>
#include "slider.h"
#include "simulation.h"
int check_simulation_values(const simulation sim)
{
+ int return_status = 0;
+ // Slider-specific parameters
+ int i;
+ for (i=0; i<sim.N; i++) {
+ if (sim.sliders[i].mass <= 0.0) {
+ fprintf(stderr, "Error: Mass of slider %d is zero or negative "
+ "(%f kg)\n", i, sim.sliders[i].mass);
+ return_status = 1;
+ }
+
+ if (sim.sliders[i].moment_of_inertia <= 0.0) {
+ fprintf(stderr, "Error: Moment of inertia of slider %d is "
+ "zero or negative (%f kg*m*m)\n",
+ i, sim.sliders[i].moment_of_inertia);
+ return_status = 1;
+ }
+ }
+
+ // General parameters
+ if (sim.N <= 0) {
+ fprintf(stderr, "Error: The number of sliders (N = %d) must be a "
+ "positive number.\n", sim.N);
+ return_status = 1;
+ }
+
+ if (sim.dt <= 0.0) {
+ fprintf(stderr, "Error: The numerical time step (dt = %f) must be a "
+ "positive value.\n", sim.dt);
+ return_status = 1;
+ }
+
+ if (sim.time > sim.time_end) {
+ fprintf(stderr, "Error: Current time (time = %f s) exceeds "
+ "total time (time_end = %f s)\n",
+ sim.time, sim.time_end);
+ return_status = 1;
+ }
+
+ if (sim.bond_length_limit <= 0.0) {
+ fprintf(stderr, "Error: The inter-slider bond length limit "
+ "(bond_length_limit = %f) must be a positive value.\n",
+ sim.bond_length_limit);
+ return_status = 1;
+ }
+
+
+ return return_status;
}
-Float find_time_step(const slider* sliders, const int N)
+Float find_time_step(const slider* sliders, const int N, const Float safety)
{
- Float largest_dt = 0.;
+ Float smallest_mass = 1.0e20;
+ Float largest_stiffness = 0.0;
+ int i;
+ for (i=0; i<N; i++) {
+ if (sliders[i].mass < smallest_mass)
+ smallest_mass = sliders[i].mass;
+ if (sliders[i].bond_parallel_stiffness > largest_stiffness)
+ largest_stiffness = sliders[i].bond_parallel_stiffness;
+ }
- return largest_dt;
+ return safety/sqrt(largest_stiffness/smallest_mass);
}
(DIR) diff --git a/simulation.h b/simulation.h
t@@ -4,16 +4,17 @@
typedef struct {
slider* sliders;
int N;
- Float t;
- Float t_end;
+ Float time;
+ Float time_end;
Float dt;
long int iteration;
Float bond_length_limit;
} simulation;
-Float find_time_step(const slider* sim.sliders, const int sim.N);
-int check_simulation(const simulation sim);
+Float find_time_step(const slider* sliders, const int N, const Float safety);
+
+int check_simulation_values(const simulation sim);
// user-defined function which sets up the simulation
simulation setup_simulation();
(DIR) diff --git a/slider.c b/slider.c
t@@ -17,9 +17,9 @@ void initialize_slider_values(slider s)
int i;
for (i=0; i<MAX_NEIGHBORS; i++) {
s.neighbors[i] = -1;
- s.neighbor_distance[i] = zeros_float3();
- s.neighbor_relative_distance_displacement[i] = zeros_float3();
- s.neighbor_relative_distance_velocity[i] = zeros_float3();
+ s.neighbor_distance[i] = zeroes_float3();
+ s.neighbor_relative_distance_displacement[i] = zeroes_float3();
+ s.neighbor_relative_distance_velocity[i] = zeroes_float3();
}
}
(DIR) diff --git a/slider.h b/slider.h
t@@ -33,10 +33,6 @@ typedef struct {
Float bond_parallel_stiffness; // Hookean elastic stiffness [N/m]
Float bond_parallel_viscosity; // viscosity [N/(m*s)]
- // slider-surface contact model parameters
- // should this value belong to the lower surface instead?
- Float coulomb_friction; // Coulomb frictional value [-]
-
// The uniquely identifying indexes of the slider neighbors which are bonded
// to this slider. A value of -1 denotes an empty field. Preallocated for
// speed.