tsimulation.h - 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
---
tsimulation.h (5250B)
---
1 #ifndef SIMULATION_
2 #define SIMULATION_
3
4 #include "arrays.h"
5
6 #define DEFAULT_SIMULATION_NAME "unnamed_simulation"
7 #define PI 3.14159265358979323846
8 #define DEG2RAD(x) (x*PI/180.0)
9
10 extern struct simulation sim;
11
12 /* Simulation settings */
13 struct simulation {
14
15 /* simulation name to use for output files */
16 char name[100];
17
18 /* gravitational acceleration magnitude [m/s^2] */
19 double G;
20
21 /* normal stress from the top wall [Pa] */
22 double P_wall;
23
24 /* optionally fix top shear velocity to this value [m/s] */
25 double v_x_fix;
26
27 /* optionally fix top shear velocity to this value [m/s] */
28 double v_x_limit;
29
30 /* bottom velocity along x [m/s] */
31 double v_x_bot;
32
33 /* stress ratio at top wall */
34 double mu_wall;
35
36 /* nonlocal amplitude [-] */
37 double A;
38
39 /* rate dependence beyond yield [-] */
40 double b;
41
42 /* bulk and critical state static yield friction coefficient [-] */
43 double mu_s;
44
45 /* material cohesion [Pa] */
46 double C;
47
48 /* representative grain size [m] */
49 double d; /* ohlala */
50
51 /* grain material density [kg/m^3] */
52 double rho_s;
53
54 /* nodes along z */
55 int nz;
56
57 /* origo of axis [m] */
58 double origo_z;
59
60 /* length of domain [m] */
61 double L_z;
62
63 /* array of cell coordinates */
64 double *z;
65
66 /* cell spacing [m] */
67 double dz;
68
69 /* current time [s] */
70 double t;
71
72 /* end time [s] */
73 double t_end;
74
75 /* time step length [s] */
76 double dt;
77
78 /* interval between output files [s] */
79 double file_dt;
80
81 /* output file number */
82 int n_file;
83
84 double transient;
85 double phi_min;
86 double phi_max;
87 double dilatancy_constant;
88
89 /* Fluid parameters */
90 int fluid; /* flag to switch fluid on (1) or off (0) */
91 double p_f_top; /* fluid pressure at the top [Pa] */
92 double p_f_mod_ampl; /* amplitude of fluid pressure variations [Pa] */
93 double p_f_mod_freq; /* frequency of fluid pressure variations [s^-1] */
94 double p_f_mod_phase; /* phase of fluid pressure variations [s^-1] */
95 double p_f_mod_pulse_time; /* single pressure pulse at this time [s] */
96 int p_f_mod_pulse_shape; /* waveform for fluid-pressure pulse */
97 double beta_f; /* adiabatic fluid compressibility [Pa^-1] */
98 double alpha; /* adiabatic grain compressibility [Pa^-1] */
99 double mu_f; /* fluid dynamic viscosity [Pa*s] */
100 double rho_f; /* fluid density [kg/m^3] */
101 double D; /* diffusivity [m^2/s], overrides k, beta_f, alpha, mu_f */
102
103 /* arrays */
104 double *mu; /* static yield friction [-] */
105 double *mu_c; /* critical-state static yield friction [-] */
106 double *sigma_n_eff; /* effective normal pressure [Pa] */
107 double *sigma_n; /* normal stress [Pa] */
108 double *p_f_ghost; /* fluid pressure [Pa] */
109 double *p_f_next_ghost; /* fluid pressure for next iteration [Pa] */
110 double *p_f_dot; /* fluid pressure change [Pa/s] */
111 double *p_f_dot_expl; /* fluid pressure change (explicit solution) [Pa/s] */
112 double *p_f_dot_impl; /* fluid pressure change (implicit solution) [Pa/s] */
113 double *p_f_dot_impl_r_norm; /* normalized residual fluid pressure change [-] */
114 double *k; /* hydraulic permeability [m^2] */
115 double *phi; /* porosity [-] */
116 double *phi_c; /* critical-state porosity [-] */
117 double *phi_dot; /* porosity change [s^-1] */
118 double *xi; /* cooperativity length */
119 double *gamma_dot_p; /* plastic shear strain rate [s^-1] */
120 double *v_x; /* shear velocity [m/s] */
121 double *d_x; /* cumulative shear displacement [m] */
122 double *g_local; /* local fluidity */
123 double *g_ghost; /* fluidity with ghost nodes */
124 double *g_r_norm; /* normalized residual of fluidity field */
125 double *I; /* inertia number [-] */
126 double *tan_psi; /* tan(dilatancy_angle) [-] */
127 double *old_val; /* temporary storage for iterative solvers */
128 double *fluid_old_val;/* temporary storage for fluid iterative solver */
129 double *tmp_ghost; /* temporary storage for iterative solvers */
130 double *p_f_dot_old; /* temporary storage for old p_f_dot */
131 };
132
133 void init_sim(struct simulation *sim);
134 void prepare_arrays(struct simulation *sim);
135 void free_arrays(struct simulation *sim);
136 void check_simulation_parameters(struct simulation *sim);
137 void lithostatic_pressure_distribution(struct simulation *sim);
138 void compute_effective_stress(struct simulation *sim);
139
140 void set_bc_neumann(double *a,
141 const int nz,
142 const int boundary,
143 const double df,
144 const double dx);
145
146 void set_bc_dirichlet(double *a,
147 const int nz,
148 const int boundary,
149 const double value);
150
151 double residual(double new_val, double old_val);
152 double kozeny_carman(const double diameter, const double porosity);
153
154 void write_output_file(struct simulation *sim, const int normalize);
155 void print_output(struct simulation *sim, FILE *fp, const int normalize);
156
157 int coupled_shear_solver(struct simulation *sim,
158 const int max_iter,
159 const double rel_tol);
160
161 void set_coupled_fluid_transient_timestep(struct simulation *sim, const double safety);
162
163 double find_flux(const struct simulation *sim);
164
165 #endif