cngf-pf

commit 41e7ad7df881868f8421486b4bfb8533a277b279
parent 031b626cdbd7f4ae6a27bfd6e621b751d96b9ed6
Author: Anders Damsgaard <anders@adamsgaard.dk>
Date:   Mon, 13 Jun 2022 21:07:38 +0200

add options to set solver iteration limit and tolerance criteria

Diffstat:
M
cngf-pf.1
 | 
9
+++++++++
M
cngf-pf.c
 | 
20
+++++++++++---------
M
simulation.c
 | 
3
+--

3 files changed, 21 insertions(+), 11 deletions(-)
diff --git a/cngf-pf.1 b/cngf-pf.1
@@ -49,6 +49,8 @@
 .Op Fl v
 .Op Fl Y Ar max-porosity
 .Op Fl y Ar min-porosity
+.Op Fl X Ar relative-tolerance
+.Op Fl x Ar max-iter
 .Op name
 .Sh DESCRIPTION
 The
@@ -201,6 +203,13 @@ Minimum granular material porosity [-] in transient simulations
 .Fl ( T )
 (default 0.20).
 .sp
+.It Fl X Ar relative-tolerance
+Sets the relative tolerance criteria for the granular solver (default
+1e-5).
+.It Fl x Ar max-iter
+Set the maximum number of iterations in the granular solver (default
+100000).
+Simulations with many cells require higher values.
 .El
 The final simulation state is written to stdout, see
 .Sx OUTPUT FORMAT
diff --git a/cngf-pf.c b/cngf-pf.c
@@ -10,10 +10,6 @@
 
 #include "arg.h"
 
-/* relative tolerance criteria for the solvers */
-#define RTOL 1e-5
-#define MAX_ITER_1D_FD_SIMPLE_SHEAR 100000
-
 /* uncomment to print time spent per time step to stdout */
 /* #define BENCHMARK_PERFORMANCE */
 
@@ -62,6 +58,8 @@ usage(void)
 	        "[-v] "
 	        "[-Y max-porosity] "
 	        "[-y min-porosity] "
+	        "[-X relative-tolerance] "
+	        "[-x max-iter] "
 	        "[name]\n", argv0);
 	exit(1);
 }
@@ -69,10 +67,10 @@ usage(void)
 int
 main(int argc, char *argv[])
 {
-	int i, normalize, dt_override, ret;
-	unsigned long iter;
+	int i, normalize = 0, dt_override = 0, ret, iter, max_iter = 100000;
 	double new_phi, new_k, filetimeclock;
 	struct simulation sim;
+	double rtol = 1e-5;
 
 #ifdef BENCHMARK_PERFORMANCE
 	clock_t t_begin, t_end;
@@ -85,10 +83,8 @@ main(int argc, char *argv[])
 #endif
 
 	init_sim(&sim);
-	normalize = 0;
 	new_phi = sim.phi[0];
 	new_k = sim.k[0];
-	dt_override = 0;
 
 	ARGBEGIN {
 	case 'A':
@@ -224,6 +220,12 @@ main(int argc, char *argv[])
 	case 'y':
 		sim.phi_min = atof(EARGF(usage()));
 		break;
+	case 'X':
+		rtol = atoi(EARGF(usage()));
+		break;
+	case 'x':
+		max_iter = atoi(EARGF(usage()));
+		break;
 	default:
 		usage();
 	} ARGEND;
@@ -280,7 +282,7 @@ main(int argc, char *argv[])
 		t_begin = clock();
 #endif
 
-		if (coupled_shear_solver(&sim, MAX_ITER_1D_FD_SIMPLE_SHEAR, RTOL)) {
+		if (coupled_shear_solver(&sim, max_iter, rtol)) {
 			free_arrays(&sim);
 			exit(10);
 		}
diff --git a/simulation.c b/simulation.c
@@ -10,7 +10,6 @@
 /* iteration limits for solvers */
 #define MAX_ITER_GRANULAR 100000
 #define MAX_ITER_DARCY 1000000
-#define MAX_ITER_STRESS 20000
 
 /* tolerance criteria when in velocity driven or velocity limited mode */
 #define RTOL_VELOCITY 1e-3
@@ -881,7 +880,7 @@ coupled_shear_solver(struct simulation *sim,
 			}
 			sim->mu_wall *= 1.0 + (vel_res_norm * 1e-3);
 		}
-		if (++stress_iter > MAX_ITER_STRESS) {
+		if (++stress_iter > max_iter) {
 			fprintf(stderr, "error: stress solution did not converge:\n");
 			fprintf(stderr, "v_x=%g, v_x_fix=%g, v_x_limit=%g, "
 			        "vel_res_norm=%g, mu_wall=%g\n",