tAdd automatic timestep - cngf-pf - continuum model for granular flows with pore-pressure dynamics (renamed from 1d_fd_simple_shear)
(HTM) git clone git://src.adamsgaard.dk/cngf-pf
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) README
(DIR) LICENSE
---
(DIR) commit 7b1e7bc6cab88ecaeb000a5c84ff018460b22387
(DIR) parent 7dd82629ef1e91bacda6383fc6223e4c806e4049
(HTM) Author: Anders Damsgaard <anders@adamsgaard.dk>
Date: Fri, 28 Jun 2019 10:34:43 +0200
Add automatic timestep
Diffstat:
M arrays.c | 13 +++++++++++++
M arrays.h | 1 +
M fluid.c | 16 ++++++++++++----
3 files changed, 26 insertions(+), 4 deletions(-)
---
(DIR) diff --git a/arrays.c b/arrays.c
t@@ -65,6 +65,19 @@ linspace(const double lower, const double upper, const int n)
return x;
}
+/* Return an array of `n-1` values with the intervals between `x` values */
+double*
+spacing(const double* x, const int n)
+{
+ int i;
+ double *dx;
+
+ dx = malloc((n-1)*sizeof(double));
+ for (i=0; i<n-1; ++i)
+ dx[i] = x[i+i] - x[i];
+ return dx;
+}
+
/* Return an array of `n` values with the value 0.0 */
double*
zeros(const int n)
(DIR) diff --git a/arrays.h b/arrays.h
t@@ -19,6 +19,7 @@ unsigned int idx2g(
unsigned int idx1g(const unsigned int i);
+double* spacing(const double* x, const int n);
double* linspace(const double lower, const double upper, const int n);
double* zeros(const int n);
double* ones(const int n);
(DIR) diff --git a/fluid.c b/fluid.c
t@@ -20,13 +20,21 @@ hydrostatic_fluid_pressure_distribution(struct simulation* sim)
void
set_largest_fluid_timestep(struct simulation* sim, const double safety)
{
- double dx_min, diff_max;
+ int i;
+ double dx_min, diff, diff_max;
+ double *dx;
- /* determine smallest cell size */
- dx_min = 0.0;
+ dx = spacing(sim->z, sim->nz);
+ dx_min = INFINITY;
+ for (i=0; i<sim->nz-1; ++i)
+ if (dx[i] < dx_min) dx_min = dx[i];
/* determine largest diffusion size */
- diff_max = 0.0;
+ diff_max = -INFINITY;
+ for (i=0; i<sim->nz; ++i) {
+ diff = sim->k[i]/(sim->phi[i]*sim->beta_f*sim->mu_f);
+ if (diff > diff_max) diff_max = diff;
+ }
sim->dt = safety*0.5*dx_min*dx_min/diff_max;
}