URI: 
       tallow negative values of sigma_n_eff, but trunkate for granular solver - cngf-pf - continuum model for granular flows with pore-pressure dynamics (renamed from 1d_fd_simple_shear)
  HTML git clone git://src.adamsgaard.dk/cngf-pf
   DIR Log
   DIR Files
   DIR Refs
   DIR README
   DIR LICENSE
       ---
   DIR commit 30733e55ce95327404e4d55f8cd679f194c8447a
   DIR parent 32c8ce3cfe0cb095da8786a5eb2150b0777fc50a
  HTML Author: Anders Damsgaard <anders@adamsgaard.dk>
       Date:   Mon, 30 Aug 2021 09:29:10 +0200
       
       allow negative values of sigma_n_eff, but trunkate for granular solver
       
       Diffstat:
         M simulation.c                        |      19 ++++++++++++-------
       
       1 file changed, 12 insertions(+), 7 deletions(-)
       ---
   DIR diff --git a/simulation.c b/simulation.c
       t@@ -15,6 +15,9 @@
        /* tolerance criteria when in velocity driven or velocity limited mode */
        #define RTOL_VELOCITY 1e-3
        
       +/* lower limit for effective normal stress sigma_n_eff for granular solver */
       +#define SIGMA_N_EFF_MIN 1e-1
       +
        
        /* Simulation settings */
        void
       t@@ -430,7 +433,7 @@ compute_inertia_number(struct simulation *sim)
                for (i = 0; i < sim->nz; ++i)
                        sim->I[i] = inertia_number(sim->gamma_dot_p[i],
                                                   sim->d,
       -                                           sim->sigma_n_eff[i],
       +                                           fmax(sim->sigma_n_eff[i], SIGMA_N_EFF_MIN),
                                                   sim->rho_s);
        }
        
       t@@ -451,8 +454,9 @@ compute_critical_state_friction(struct simulation *sim)
        
                if (sim->fluid)
                        for (i = 0; i < sim->nz; ++i)
       -                        sim->mu_c[i] = sim->mu_wall / (sim->sigma_n_eff[i]
       -                                                       / (sim->P_wall - sim->p_f_top));
       +                        sim->mu_c[i] = sim->mu_wall
       +                                       / (fmax(sim->sigma_n_eff[i], SIGMA_N_EFF_MIN)
       +                                          / (sim->P_wall - sim->p_f_top));
                else
                        for (i = 0; i < sim->nz; ++i)
                                sim->mu_c[i] = sim->mu_wall
       t@@ -551,7 +555,7 @@ compute_effective_stress(struct simulation *sim)
                                                         / 2.0);
                                /* sim->sigma_n_eff[i] = sim->sigma_n[i] - sim->p_f_ghost[i + 1]; */
                                if (sim->sigma_n_eff[i] < 0)
       -                                errx(1, "%s: sigma_n_eff[%d] is negative with value %g\n",
       +                                warnx("%s: sigma_n_eff[%d] is negative with value %g\n",
                                             __func__, i, sim->sigma_n_eff[i]);
                        }
                else
       t@@ -579,7 +583,7 @@ compute_cooperativity_length(struct simulation *sim)
                        sim->xi[i] = cooperativity_length(sim->A,
                                                          sim->d,
                                                          sim->mu[i],
       -                                                  sim->sigma_n_eff[i],
       +                                                  fmax(sim->sigma_n_eff[i], SIGMA_N_EFF_MIN),
                                                          sim->mu_s,
                                                          sim->C);
        }
       t@@ -605,7 +609,8 @@ compute_local_fluidity(struct simulation *sim)
                int i;
        
                for (i = 0; i < sim->nz; ++i)
       -                sim->g_local[i] = local_fluidity(sim->sigma_n_eff[i],
       +                sim->g_local[i] = local_fluidity(fmax(sim->sigma_n_eff[i],
       +                                                      SIGMA_N_EFF_MIN),
                                                         sim->mu[i],
                                                         sim->mu_s,
                                                         sim->C,
       t@@ -770,7 +775,7 @@ print_output(struct simulation *sim, FILE *fp, const int norm)
                                sim->gamma_dot_p[i],
                                sim->phi[i],
                                sim->I[i],
       -                        sim->mu[i] * sim->sigma_n_eff[i],
       +                        sim->mu[i] * fmax(sim->sigma_n_eff[i], SIGMA_N_EFF_MIN),
                                sim->d_x[i]);
        
                free(v_x_out);