commit 5419eee7ed46ffaae98e60ea4391d713d1731668
parent 40b8d91e4d2d99c2574ba1d3033f5f1bbc71f3ca
Author: Anders Damsgaard <anders@adamsgaard.dk>
Date:   Thu,  5 Sep 2019 18:02:29 +0200
Separate tolerance criteria, add triangular pressure pulse forcing
Diffstat:
6 files changed, 295 insertions(+), 50 deletions(-)
diff --git a/fluid.c b/fluid.c
@@ -54,10 +54,23 @@ static double
 sine_wave(const double time,
           const double amplitude,
           const double frequency,
-          const double phase,
-          const double base_value)
+          const double phase)
 {
-	return amplitude*sin(2.0*PI*frequency*time + phase) + base_value;
+	return amplitude*sin(2.0*PI*frequency*time + phase);
+}
+
+static double
+pulse(const double time,
+      const double peak_amplitude,
+      const double frequency,
+      const double peak_time)
+{
+	if (peak_time - 1.0/(2.0*frequency) < time && time <= peak_time)
+		return peak_amplitude*2.0*frequency*(time - peak_time) + peak_amplitude;
+	else if (peak_time < time && time < peak_time + 1.0/(2.0*frequency))
+		return peak_amplitude*2.0*frequency*(peak_time - time) + peak_amplitude;
+	else
+		return 0.0;
 }
 
 static void
@@ -167,11 +180,16 @@ darcy_solver_1d(struct simulation* sim,
 	/* TODO: values other than 1.0 do not work! */
 	theta = 1.0;
 
-	p_f_top = sine_wave(sim->t,
-	                    sim->p_f_mod_ampl,
-	                    sim->p_f_mod_freq,
-	                    sim->p_f_mod_phase,
-	                    sim->p_f_top);
+	if (isnan(sim->p_f_mod_pulse_time))
+		p_f_top = sim->p_f_top + sine_wave(sim->t,
+		                                   sim->p_f_mod_ampl,
+		                                   sim->p_f_mod_freq,
+		                                   sim->p_f_mod_phase);
+	else
+		p_f_top = sim->p_f_top + pulse(sim->t,
+		                               sim->p_f_mod_ampl,
+		                               sim->p_f_mod_freq,
+		                               sim->p_f_mod_pulse_time);
 
 	/* set fluid BCs (1 of 2) */
 	set_fluid_bcs(sim, p_f_top);
diff --git a/main.c b/main.c
@@ -12,8 +12,17 @@
 
 #include "parameter_defaults.h"
 
+/* relative tolerance criteria for granular fluidity solver */
+#define RTOL_GRANULAR 1e-5
+#define MAX_ITER_GRANULAR 10000
+
+/* relative tolerance criteria for fluid-pressure solver */
+#define RTOL_DARCY 1e-5
+#define MAX_ITER_DARCY 10000
+
 /* relative tolerance criteria when shear velocity is restricted */
-#define RTOL 1e-3
+#define RTOL_STRESS 1e-3
+#define MAX_ITER_STRESS 20000
 
 static void
 usage(void)
@@ -54,8 +63,11 @@ usage(void)
 	        " -k, --fluid-permeability VAL    fluid intrinsic permeability [m^2] (default %g)\n"
 	        " -O, --fluid-pressure-top VAL    fluid pressure at +z edge [Pa] (default %g)\n"
 	        " -a, --fluid-pressure-ampl VAL   amplitude of pressure variations [Pa] (default %g)\n"
-	        " -q, --fluid-pressure-freq VAL   frequency of pressure variations [s^-1] (default %g)\n"
+	        " -q, --fluid-pressure-freq VAL   frequency of sinusoidal pressure variations [s^-1]\n"
+	        "                                 (default %g)\n"
 	        " -H, --fluid-pressure-phase VAL  fluid pressure at +z edge [Pa] (default %g)\n"
+	        " -u, --fluid-pressure-pulse-time VAL fluid pressure pulse peak time [s]\n"
+	        "                                 (default %g)\n"
 	        " -t, --time VAL                  simulation start time [s] (default %g)\n"
 	        " -T, --time-end VAL              simulation end time [s] (default %g)\n"
 	        " -I, --file-interval VAL         interval between output files [s] (default %g)\n"
@@ -99,6 +111,7 @@ usage(void)
 	        sim.p_f_mod_ampl,
 	        sim.p_f_mod_freq,
 	        sim.p_f_mod_phase,
+	        sim.p_f_mod_pulse_time,
 	        sim.t,
 	        sim.t_end,
 	        sim.file_dt);
@@ -130,40 +143,41 @@ main(int argc, char* argv[])
 
 	norm = 0;
 
-	optstring = "hvNn:G:P:m:s:l:V:A:b:f:C:Fp:d:r:o:L:c:i:R:k:O:a:q:H:t:T:D:I:";
+	optstring = "hvNn:G:P:m:s:l:V:A:b:f:C:Fp:d:r:o:L:c:i:R:k:O:a:q:H:T:t:T:D:I:";
 	const struct option longopts[] = {
-		{"help",                 no_argument,       NULL, 'h'},
-		{"version",              no_argument,       NULL, 'v'},
-		{"normalize",            no_argument,       NULL, 'N'},
-		{"gravity",              required_argument, NULL, 'G'},
-		{"normal-stress",        required_argument, NULL, 'P'},
-		{"stress-ratio",         required_argument, NULL, 'm'},
-		{"set-shear-velocity",   required_argument, NULL, 's'},
-		{"limit-shear-velocity", required_argument, NULL, 'l'},
-		{"velocity-bottom",      required_argument, NULL, 'V'},
-		{"nonlocal-amplitude",   required_argument, NULL, 'A'},
-		{"rate-dependence",      required_argument, NULL, 'b'},
-		{"friction-coefficient", required_argument, NULL, 'f'},
-		{"cohesion",             required_argument, NULL, 'C'},
-		{"porosity",             required_argument, NULL, 'p'},
-		{"grain-size",           required_argument, NULL, 'd'},
-		{"density",              required_argument, NULL, 'r'},
-		{"resolution",           required_argument, NULL, 'n'},
-		{"origo",                required_argument, NULL, 'o'},
-		{"length",               required_argument, NULL, 'L'},
-		{"fluid",                no_argument,       NULL, 'F'},
-		{"fluid-compressiblity", required_argument, NULL, 'c'},
-		{"fluid-viscosity",      required_argument, NULL, 'i'},
-		{"fluid-density",        required_argument, NULL, 'R'},
-		{"fluid-permeability",   required_argument, NULL, 'k'},
-		{"fluid-pressure-top",   required_argument, NULL, 'O'},
-		{"fluid-pressure-ampl",  required_argument, NULL, 'a'},
-		{"fluid-pressure-freq",  required_argument, NULL, 'q'},
-		{"fluid-pressure-phase", required_argument, NULL, 'H'},
-		{"time",                 required_argument, NULL, 't'},
-		{"time-end",             required_argument, NULL, 'T'},
-		{"file-interval",        required_argument, NULL, 'I'},
-		{NULL,                   0,                 NULL, 0}
+		{"help",                      no_argument,       NULL, 'h'},
+		{"version",                   no_argument,       NULL, 'v'},
+		{"normalize",                 no_argument,       NULL, 'N'},
+		{"gravity",                   required_argument, NULL, 'G'},
+		{"normal-stress",             required_argument, NULL, 'P'},
+		{"stress-ratio",              required_argument, NULL, 'm'},
+		{"set-shear-velocity",        required_argument, NULL, 's'},
+		{"limit-shear-velocity",      required_argument, NULL, 'l'},
+		{"velocity-bottom",           required_argument, NULL, 'V'},
+		{"nonlocal-amplitude",        required_argument, NULL, 'A'},
+		{"rate-dependence",           required_argument, NULL, 'b'},
+		{"friction-coefficient",      required_argument, NULL, 'f'},
+		{"cohesion",                  required_argument, NULL, 'C'},
+		{"porosity",                  required_argument, NULL, 'p'},
+		{"grain-size",                required_argument, NULL, 'd'},
+		{"density",                   required_argument, NULL, 'r'},
+		{"resolution",                required_argument, NULL, 'n'},
+		{"origo",                     required_argument, NULL, 'o'},
+		{"length",                    required_argument, NULL, 'L'},
+		{"fluid",                     no_argument,       NULL, 'F'},
+		{"fluid-compressiblity",      required_argument, NULL, 'c'},
+		{"fluid-viscosity",           required_argument, NULL, 'i'},
+		{"fluid-density",             required_argument, NULL, 'R'},
+		{"fluid-permeability",        required_argument, NULL, 'k'},
+		{"fluid-pressure-top",        required_argument, NULL, 'O'},
+		{"fluid-pressure-ampl",       required_argument, NULL, 'a'},
+		{"fluid-pressure-freq",       required_argument, NULL, 'q'},
+		{"fluid-pressure-phase",      required_argument, NULL, 'H'},
+		{"fluid-pressure-pulse-time", required_argument, NULL, 'u'},
+		{"time",                      required_argument, NULL, 't'},
+		{"time-end",                  required_argument, NULL, 'T'},
+		{"file-interval",             required_argument, NULL, 'I'},
+		{NULL,                        0,                 NULL, 0}
 	};
 
 	new_phi = sim.phi[0];
@@ -262,6 +276,9 @@ main(int argc, char* argv[])
 			case 'H':
 				sim.p_f_mod_phase = atof(optarg);
 				break;
+			case 'u':
+				sim.p_f_mod_pulse_time = atof(optarg);
+				break;
 			case 't':
 				sim.t = atof(optarg);
 				break;
@@ -317,7 +334,7 @@ main(int argc, char* argv[])
 		stressiter = 0;
 		do {
 			if (sim.fluid) {
-				if (darcy_solver_1d(&sim, 10000, 1e-5))
+				if (darcy_solver_1d(&sim, MAX_ITER_DARCY, RTOL_DARCY))
 					exit(1);
 			}
 
@@ -325,7 +342,8 @@ main(int argc, char* argv[])
 			compute_friction(&sim);
 			compute_cooperativity_length(&sim);
 
-			if (implicit_1d_jacobian_poisson_solver(&sim, 10000, 1e-5))
+			if (implicit_1d_jacobian_poisson_solver(&sim, MAX_ITER_GRANULAR,
+			                                        RTOL_GRANULAR))
 				exit(1);
 
 			compute_shear_strain_rate_plastic(&sim);
@@ -344,7 +362,7 @@ main(int argc, char* argv[])
 				sim.mu_wall *= 1.0 + (res_norm*1e-2);
 			}
 
-			if (++stressiter > 10000) {
+			if (++stressiter > MAX_ITER_STRESS) {
 				fprintf(stderr, "error: stress solution did not converge:\n");
 				fprintf(stderr,
 				        "v_x=%g, v_x_fix=%g, v_x_limit=%g, "
@@ -356,7 +374,7 @@ main(int argc, char* argv[])
 			}
 
 		} while ((!isnan(sim.v_x_fix) || !isnan(sim.v_x_limit))
-		         && fabs(res_norm) > RTOL);
+		         && fabs(res_norm) > RTOL_STRESS);
 
 		if (!isnan(sim.v_x_limit))
 			sim.mu_wall = mu_wall_orig;
diff --git a/parameter_defaults.h b/parameter_defaults.h
@@ -17,7 +17,7 @@ struct simulation init_sim(void)
 
 	sim.G = 9.81;
 
-	sim.P_wall = 120e3; /* larger normal stress deepens the shear depth */
+	sim.P_wall = 120e3;
 	sim.mu_wall = 0.45;
 	sim.v_x_bot = 0.0;
 	sim.v_x_fix = NAN;
@@ -66,7 +66,7 @@ struct simulation init_sim(void)
 
 	/* spatial settings */
 	sim.origo_z = 0.0;
-	sim.L_z = 0.7;       /* Damsgaard et al 2013 */
+	sim.L_z = 1.0;
 
 	/* temporal settings */
 	sim.t = 0.0;
@@ -105,6 +105,7 @@ struct simulation init_sim(void)
 	sim.p_f_mod_ampl = 0.0;
 	sim.p_f_mod_freq = 1.0;
 	sim.p_f_mod_phase = 0.0;
+	sim.p_f_mod_pulse_time = NAN;
 
 	return sim;
 }
diff --git a/simulation.h b/simulation.h
@@ -84,6 +84,7 @@ struct simulation {
 	double p_f_mod_ampl;  /* amplitude of fluid pressure variations [Pa] */
 	double p_f_mod_freq;  /* frequency of fluid pressure variations [s^-1] */
 	double p_f_mod_phase; /* phase of fluid pressure variations [s^-1] */
+	double p_f_mod_pulse_time; /* single pressure pulse at this time [s] */
 	double beta_f;        /* adiabatic fluid compressibility [Pa^-1] */
 	double mu_f;          /* fluid dynamic viscosity [Pa*s] */
 	double rho_f;         /* fluid density [kg/m^3] */
diff --git a/test/Makefile b/test/Makefile
@@ -5,7 +5,8 @@ TESTS := simple_shear_dry \
          simple_shear_dry_cohesive \
          simple_shear_wet \
          simple_shear_wet_norm \
-         simple_shear_wet_vari
+         simple_shear_wet_vari \
+         simple_shear_wet_vari_pulse
 
 STANDARDS := $(addsuffix .std,$(TESTS))
 
@@ -22,6 +23,8 @@ SIMPLE_SHEAR_WET_VARI_OPTS =     --length 8.0 \
                                  --fluid-pressure-freq $$( echo "1.0/(3600*24)" | bc -l ) \
                                  --file-interval $$( echo "60*1000" | bc -l ) \
                                  --time-end $$( echo "60*20" | bc -l )
+SIMPLE_SHEAR_WET_VARI_PULSE_OPTS = $(SIMPLE_SHEAR_WET_VARI_OPTS) \
+                                 --fluid-pressure-pulse-time 500
 
 test: $(TESTS)
 
@@ -65,6 +68,10 @@ simple_shear_wet_norm.std: $(BIN)
 simple_shear_wet_vari.std: $(BIN)
 	./$< $(SIMPLE_SHEAR_WET_VARI_OPTS) > $@
 
+simple_shear_wet_vari_pulse.std: $(BIN)
+	./$< $(SIMPLE_SHEAR_WET_VARI_PULSE_OPTS) > $@
+
+
 
 $(BIN):
 	make -C ..
diff --git a/test/simple_shear_wet_vari_pulse.std b/test/simple_shear_wet_vari_pulse.std
@@ -0,0 +1,200 @@
+0	0	233416	69620	0.19278884052507114	2.1720851657206607e-70
+0.040201005025125629	2.9283675330214125e-71	232745.56783919601	69521.407035175886	0.19334417586456693	7.2843142383907636e-70
+0.080402010050251257	1.1858964030866926e-70	232075.1356783919	69422.814070351757	0.19390271977417128	2.2214858788390719e-69
+0.12060301507537688	3.8800577912168823e-70	231404.70351758791	69324.221105527642	0.19446450014176042	6.701726452973846e-69
+0.16080402010050251	1.1987094751429542e-69	230734.27135678395	69225.628140703513	0.19502954517934004	2.0166254438528989e-68
+0.20100502512562815	3.6350445182896457e-69	230063.8391959799	69127.035175879399	0.19559788342776782	6.0603834198273949e-68
+0.24120603015075376	1.0945670413958273e-68	229393.40703517586	69028.442211055284	0.19616954376155879	1.8185181915475712e-67
+0.28140703517587939	3.2856096691234426e-68	228722.9748743719	68929.849246231155	0.19674455539377558	5.4502185364724428e-67
+0.32160804020100503	9.8423296022403995e-68	228052.54271356785	68831.256281407041	0.19732294788100496	1.6309840833628432e-66
+0.36180904522613067	2.9439476878115838e-67	227382.1105527638	68732.663316582912	0.19790475112842174	4.8747903848740155e-66
+0.4020100502512563	8.7922574846911803e-67	226711.67839195981	68634.070351758797	0.19848999539494341	1.4547670619737995e-65
+0.44221105527638194	2.6223511876590728e-66	226041.24623115576	68535.477386934668	0.19907871129847607	4.3360245299850127e-65
+0.48241206030150752	7.809757917913618e-66	225370.81407035171	68436.884422110554	0.19967092982125362	1.2903674241508181e-64
+0.52261306532663321	2.3227555061629377e-65	224700.38190954775	68338.291457286425	0.20026668231527339	3.8351770394992951e-64
+0.56281407035175879	6.8979566422595218e-65	224029.94974874376	68239.69849246231	0.20086600050782868	1.1380812826040254e-63
+0.60301507537688448	2.0457258398493498e-64	223359.51758793966	68141.105527638196	0.20146891650714141	3.3728763118632009e-63
+0.64321608040201006	6.0578375081819407e-64	222689.08542713569	68042.512562814067	0.20207546280809568	9.9801277749773199e-63
+0.68341708542713564	1.7913782996902266e-63	222018.6532663317	67943.919597989952	0.2026856722980766	2.9491664403191805e-62
+0.72361809045226133	5.2892252859868059e-63	221348.22110552766	67845.326633165823	0.20329957826291392	8.7008943784127398e-62
+0.76381909547738691	1.5594996338350084e-62	220677.78894472361	67746.733668341709	0.2039172143929347	2.5635605492753655e-61
+0.8040201005025126	4.5909785308031583e-62	220007.35678391962	67648.140703517594	0.20453861478912627	7.5408037562082724e-61
+0.84422110552763818	1.349591106513562e-61	219336.92462311557	67549.547738693465	0.2051638139694128	2.2151019679151997e-60
+0.88442211055276387	3.9611150098743099e-61	218666.49246231152	67450.954773869351	0.20579284687504656	6.4961657096098614e-60
+0.92462311557788945	1.1609078010863288e-60	217996.06030150756	67352.361809045222	0.20642574887711768	1.9024307964960084e-59
+0.96482412060301503	3.3969378083284191e-60	217325.62814070351	67253.768844221107	0.20706255578318433	5.5621246430146995e-59
+1.0050251256281406	9.9249915656232056e-60	216655.19597989946	67155.175879396993	0.20770330384402574	1.6238533721270781e-58
+1.0452261306532664	2.8951648062444419e-59	215984.7638190955	67056.582914572864	0.20834802976052097	4.7328808035842766e-58
+1.085427135678392	8.4325004595633886e-59	215314.33165829146	66957.989949748735	0.20899677069065697	1.377412243763088e-57
+1.1256281407035176	2.4520579338657535e-58	214643.89949748741	66859.39698492462	0.20964956425666675	4.0019096211746685e-57
+1.1658291457286432	7.1192116970981506e-58	213973.46733668342	66760.804020100506	0.21030644855230254	1.1609544986040587e-56
+1.206030150753769	2.0635484333470994e-57	213303.03517587937	66662.211055276377	0.21096746215024637	3.3621728182977443e-56
+1.2462311557788945	5.9718764970784075e-57	212632.60301507538	66563.618090452263	0.21163264410965971	9.7219660585316269e-56
+1.2864321608040201	1.7253549127826681e-56	211962.17085427136	66465.025125628134	0.21230203398387765	2.8063160668986332e-55
+1.3266331658291457	4.9767491410850598e-56	211291.73869346737	66366.432160804019	0.21297567182824878	8.0878431429021999e-55
+1.3668341708542713	1.4330916616435051e-55	210621.30653266332	66267.839195979905	0.21365359820812507	2.3268491594933103e-54
+1.4070351758793971	4.119911830902165e-55	209950.8743718593	66169.246231155776	0.21433585420700471	6.6834651710309167e-54
+1.4472361809045227	1.1823653844616179e-54	209280.44221105531	66070.653266331661	0.21502248143483163	1.9163058259113612e-53
+1.4874371859296482	3.3875607696170614e-54	208610.01005025127	65972.060301507532	0.21571352203645514	5.4854235205741663e-53
+1.5276381909547738	9.6885914970914972e-54	207939.57788944722	65873.467336683418	0.2164090187002525	1.5673813934592659e-52
+1.5678391959798996	2.7662505959306222e-53	207269.14572864323	65774.874371859303	0.21710901466692012	4.4710112224759123e-52
+1.6080402010050252	7.884029202710329e-53	206598.71356783918	65676.281407035174	0.21781355373843461	1.273047428436452e-51
+1.6482412060301508	2.2430959805609797e-52	205928.28140703513	65577.68844221106	0.21852268028718988	3.6185489874712429e-51
+1.6884422110552764	6.370305269598026e-52	205257.84924623117	65479.095477386931	0.21923643926531236	1.0266433106479653e-50
+1.7286432160804019	1.8059306458715307e-51	204587.41708542712	65380.502512562816	0.21995487621416077	2.9076390457929238e-50
+1.7688442211055277	5.1102337336364442e-51	203916.98492462307	65281.909547738695	0.22067803727401142	8.2194539308152227e-50
+1.8090452261306533	1.4434253200897425e-50	203246.55276381911	65183.316582914573	0.22140596919393687	2.3193498424811686e-49
+1.8492462311557789	4.0694152630664609e-50	202576.12060301506	65084.723618090451	0.22213871934187998	6.5321499831545872e-49
+1.8894472361809045	1.1451672034675639e-49	201905.68844221107	64986.130653266329	0.22287633571492854	1.8363363719377829e-48
+1.9296482412060301	3.2164315426731421e-49	201235.25628140703	64887.537688442215	0.22361886694979571	5.1522700437738756e-48
+1.9698492462311559	9.017044050715112e-49	200564.82412060298	64788.944723618086	0.22436636233351046	1.4429023613754398e-47
+2.0100502512562812	2.5229352342863439e-48	199894.39195979899	64690.351758793971	0.22511887181432289	4.0328116876718963e-47
+2.050251256281407	7.0455948917885156e-48	199223.95979899497	64591.75879396985	0.2258764460128305	1.1250115898036651e-46
+2.0904522613065328	1.9636494126404889e-47	198553.52763819098	64493.165829145728	0.2266391362333289	3.1319861846108232e-46
+2.1306532663316582	5.462167278701703e-47	197883.09547738693	64394.572864321613	0.22740699447539403	8.7025631918272689e-46
+2.170854271356784	1.5162914475903345e-46	197212.66331658291	64295.979899497484	0.22818007344570002	2.4130608653039085e-45
+2.2110552763819098	4.2008787747314292e-46	196542.23115577892	64197.38693467337	0.22895842657007953	6.6779109762634732e-45
+2.2512562814070352	1.161436627038477e-45	195871.79899497487	64098.793969849248	0.22974210800583134	1.8441050145437681e-44
+2.291457286432161	3.204612844439332e-45	195201.36683417085	64000.201005025127	0.23053117265428161	5.0824008407846278e-44
+2.3316582914572863	8.8233643309602772e-45	194530.93467336684	63901.608040201005	0.23132567617360517	1.397664432272085e-43
+2.3718592964824121	2.4243890492242722e-44	193860.50251256279	63803.015075376883	0.23212567499191256	3.8358558826190083e-43
+2.4120603015075379	6.6470540149805636e-44	193190.0703517588	63704.422110552769	0.2329312263206095	1.0503879102318776e-42
+2.4522613065326633	1.8186543265812015e-43	192519.63819095478	63605.82914572864	0.23374238816803602	2.8704479511443236e-42
+2.4924623115577891	4.9648905928617218e-43	191849.20603015073	63507.236180904525	0.23455921935339086	7.8262627123727927e-42
+2.5326633165829144	1.3525322328128644e-42	191178.77386934674	63408.643216080403	0.23538177952094932	2.1294073941476466e-41
+2.5728643216080402	3.6762371414751171e-42	190508.34170854272	63310.050251256282	0.23621012915458139	5.7802159602973535e-41
+2.613065326633166	9.9706131397302968e-42	189837.9095477387	63211.45728643216	0.23704432959257704	1.5657260295659757e-40
+2.6532663316582914	2.6979701880847041e-41	189167.47738693468	63112.864321608038	0.23788444304278719	4.2310108243527905e-40
+2.6934673366834172	7.2844787325365078e-41	188497.04522613063	63014.271356783924	0.23873053259808774	1.1408940004323862e-39
+2.7336683417085426	1.9621638752009283e-40	187826.61306532664	62915.678391959802	0.23958266225217439	3.0688685548438529e-39
+2.7738693467336684	5.2735268865334752e-40	187156.1809045226	62817.08542713568	0.24044089691569778	8.2370154906897113e-39
+2.8140703517587942	1.4139059557134978e-39	186485.74874371861	62718.492462311558	0.24130530243274545	2.2053012518121235e-38
+2.8542713567839195	3.7822595023747532e-39	185815.31658291459	62619.899497487437	0.2421759455976821	5.8912794473198728e-38
+2.8944723618090453	1.0092830807021785e-38	185144.88442211054	62521.306532663315	0.24305289417235429	1.5697546120309491e-37
+2.9346733668341707	2.6869997451207669e-38	184474.45226130658	62422.713567839193	0.24393621690367109	4.1733202027412378e-37
+2.9748743718592965	7.135589285356945e-38	183804.02010050253	62324.120603015079	0.24482598354156981	1.1065866481337494e-36
+3.0150753768844223	1.8904517863870576e-37	183133.58793969848	62225.527638190957	0.24572226485737519	2.9275209839052659e-36
+3.0552763819095476	4.9955686090643399e-37	182463.15577889449	62126.934673366835	0.24662513266256436	7.7239780964097383e-36
+3.0954773869346734	1.3169137443143347e-36	181792.72361809044	62028.341708542714	0.24753465982794698	2.0331752474771531e-35
+3.1356783919597992	3.4624898107609384e-36	181122.2914572864	61929.748743718592	0.24845092030327054	5.3371204652859273e-35
+3.1758793969849246	9.0813660563330111e-36	180451.85929648243	61831.15577889447	0.24937398913726344	1.3976954660860532e-34
+3.2160804020100504	2.3754565586509047e-35	179781.42713567839	61732.562814070348	0.25030394249812671	3.6499583831312893e-34
+3.2562814070351758	6.1980027681783656e-35	179110.99497487437	61633.969849246234	0.25124085769448484	9.5085836961995578e-34
+3.2964824120603016	1.6127365674997321e-34	178440.56281407035	61535.376884422112	0.25218481319681019	2.4699290230712152e-33
+3.3366834170854274	4.1856212098133067e-34	177770.13065326633	61436.78391959799	0.25313588865933129	6.4000505477550174e-33
+3.3768844221105527	1.0832695412563893e-33	177099.69849246231	61338.190954773869	0.25409416494243936	1.6534597079342083e-32
+3.4170854271356785	2.7962294218093856e-33	176429.26633165829	61239.597989949747	0.25505972413560529	4.2609877028755785e-32
+3.4572864321608039	7.1971696049260488e-33	175758.83417085427	61141.005025125633	0.25603264958082123	1.0947338705502698e-31
+3.4974874371859297	1.8474977824553566e-32	175088.40201005025	61042.412060301504	0.25701302589658082	2.8053547946323445e-31
+3.5376884422110555	4.7285804057402255e-32	174417.9698492462	60943.819095477389	0.25800093900241255	7.1666930254211113e-31
+3.5778894472361809	1.2069331058118852e-31	173747.53768844221	60845.226130653267	0.25899647614398064	1.8260117247791829e-30
+3.6180904522613067	3.0713669068932076e-31	173077.1055276382	60746.633165829146	0.25999972591877024	4.6377790801897888e-30
+3.658291457286432	7.7939347682014283e-31	172406.67336683415	60648.040201005024	0.26101077830237079	1.1747387555004199e-29
+3.6984924623115578	1.9717354330360353e-30	171736.24120603019	60549.447236180902	0.26202972467537566	2.9659506160870327e-29
+3.7386934673366836	4.9737814003759958e-30	171065.80904522608	60450.854271356817	0.26305665785091503	7.4675893437581518e-29
+3.778894472361809	1.2507238178729955e-29	170395.376884422	60352.261306532746	0.26409167210283635	1.8739473736155474e-28
+3.8190954773869348	3.1358153010919852e-29	169724.94472361793	60253.668341708719	0.26513486319455593	4.6891650645072368e-28
+3.8592964824120601	7.8369118821104833e-29	169054.5125628137	60155.075376884764	0.26618632840859457	1.1693977745283512e-27
+3.8994974874371859	1.9526318739066792e-28	168384.08040200942	60056.482412060912	0.26724616657681965	2.9077399556678817e-27
+3.9396984924623117	4.8491833497626294e-28	167713.64824120494	59957.889447237278	0.26831447811141301	7.2051717961916749e-27
+3.9798994974874371	1.2005014874613362e-27	167043.21608040005	59859.296482414014	0.26939136503658384	1.7800130918066199e-26
+4.0201005025125625	2.9620703881614054e-27	166372.78391959454	59760.703517591392	0.27047693102104864	4.3819026404914217e-26
+4.0603015075376883	7.2851270296058486e-27	165702.3517587879	59662.110552769875	0.27157128141129994	1.0753603395593054e-25
+4.1005025125628141	1.7855885067641982e-26	165031.91959797937	59563.517587950271	0.2726745232656857	2.6294760619614877e-25
+4.1407035175879399	4.3621054493494706e-26	164361.48743716758	59464.924623133928	0.27378676538932323	6.4090858946808642e-25
+4.1809045226130657	1.0618799678598143e-25	163691.05527635023	59366.331658323114	0.27490811836987106	1.5563526895256073e-24
+4.2211055276381906	2.5762315443901408e-25	163020.6231155236	59267.738693521664	0.27603869461418395	3.7669495466191871e-24
+4.2613065326633164	6.2275743923010977e-25	162350.19095468108	59169.145728735995	0.27717860838587766	9.0827153341785057e-24
+4.3015075376884422	1.5001721948951618e-24	161679.75879381219	59070.552763976724	0.2783279758438274	2.1825692047168171e-23
+4.341708542713568	3.6003618684720744e-24	161009.32663289941	58971.959799261429	0.27948691508163259	5.224221813022569e-23
+4.3819095477386938	8.6098523269333264e-24	160338.89447191372	58873.366834618944	0.2806555461680732	1.2461107515422364e-22
+4.4221105527638196	2.0510904972136419e-23	159668.46231080804	58774.773870096447	0.28183399118859009	2.9603868454942687e-22
+4.4623115577889445	4.8682549935803207e-23	158998.03014950582	58676.180905770569	0.28302237428782301	7.0076966847121138e-22
+4.5025125628140703	1.1509508255003373e-22	158327.59798788303	58577.587941765196	0.28422082171324231	1.6520117487789843e-21
+4.5427135678391961	2.7107726154323207e-22	157657.16582574064	58478.994978279465	0.28542946185991164	3.8800567024558085e-21
+4.5829145728643219	6.3588505250055793e-22	156986.73366276029	58380.402015631684	0.28664842531642981	9.0745938000634807e-21
+4.6231155778894477	1.4858301768714112e-21	156316.3014984359	58281.809054327889	0.28787784491210133	2.1142384968724977e-20
+4.6633165829145726	3.4574995730334582e-21	155645.86933196767	58183.216095168042	0.28911785576539922	4.9045276229530914e-20
+4.7035175879396984	8.0132902521512822e-21	154975.43716209781	58084.623139409734	0.29036859533380049	1.1332529314305586e-19
+4.7437185929648242	1.849328484519627e-20	154305.0049868605	57986.030189018908	0.29163020346509094	2.6068986550199406e-19
+4.78391959798995	4.2503246707952761e-20	153634.57280319967	57887.437247051595	0.29290282245027865	5.9724780133606767e-19
+4.8241206030150758	9.7260204078509937e-20	152964.14060639151	57788.844318231604	0.29418659707829398	1.3620793145926098e-18
+4.8643216080402008	2.2161548185225409e-19	152293.70838917559	57690.251409819379	0.29548167469273084	3.0933375346218854e-18
+4.9045226130653266	5.0271127014272074e-19	151623.27614045597	57591.658532910864	0.29678820525098221	6.9922577337253585e-18
+4.9447236180904524	1.1353592081731908e-18	150952.84384337155	57493.065704367167	0.29810634138626718	1.5737117458507943e-17
+4.9849246231155782	2.5523924243805962e-18	150282.41147244704	57394.472949663497	0.2994362384732584	3.5248701253159208e-17
+5.025125628140704	5.7121372027724211e-18	149611.97898941307	57295.880307069339	0.30077805469830937	7.8598651362496637e-17
+5.0653266331658289	1.2723064714376315e-17	148941.5463371113	57197.287833742972	0.30213195113570196	1.7439682185114688e-16
+5.1055276381909547	2.8207111769939575e-17	148271.11343066505	57098.695614561082	0.30349809183191317	3.8516567050713614e-16
+5.1457286432160805	6.2230839079620746e-17	147600.6801447687	57000.10377482928	0.30487664390071512	8.4634021682831918e-16
+5.1859296482412063	1.3663390282586345e-16	146930.24629550803	56901.512498461809	0.30626777763303697	1.8507762106877877e-15
+5.2261306532663321	2.9848664679932295e-16	146259.81161452626	56802.92205381544	0.30767166662705236	4.0260870063398042e-15
+5.266331658291457	6.4881957035895423e-16	145589.37571255281	56704.332830160769	0.30908848794603405	8.714531473545827e-15
+5.3065326633165828	1.4030234878666505e-15	144918.93802825327	56605.745388832191	0.31051842231432059	1.8760822448003947e-14
+5.3467336683417086	3.0182832749672654e-15	144248.49775696898	56507.160534488292	0.31196165436548506	4.0179587204127797e-14
+5.3869346733668344	6.4583287272549577e-15	143578.05375210484	56408.579413724307	0.31341837296175429	8.5571130625656332e-14
+5.4271356783919602	1.3745181663065431e-14	142907.60438958817	56310.003650612845	0.31488877161024303	1.8126046677828551e-13
+5.4673366834170851	2.9091348513181177e-14	142237.14738284383	56211.435531729046	0.31637304901003482	3.8173590039662922e-13
+5.5075376884422109	6.1229175785712276e-14	141566.6795319615	56112.878256983211	0.31787140977506895	7.9942845340421109e-13
+5.5477386934673367	1.2812939641046459e-13	140896.19638602331	56014.336277293274	0.31938406539173214	1.6641429880407137e-12
+5.5879396984924625	2.6657929749548306e-13	140225.69179173704	55915.815745951419	0.32091123548767314	3.4439412894898346e-12
+5.6281407035175883	5.5132454704350576e-13	139555.15729440836	55817.325117651919	0.32245314951039111	7.0830380825070653e-12
+5.6683417085427132	1.1333830928177438e-12	138884.58134870522	55718.875937726923	0.32401004894140123	1.447870632613417e-11
+5.708542713567839	2.3155337045977558e-12	138213.94828645734	55620.483874346697	0.3255821902051057	2.9405996468027806e-11
+5.7487437185929648	4.7012146338303628e-12	137543.23697676291	55522.170058412972	0.32716984847173891	5.934381311466111e-11
+5.7889447236180906	9.4835113120226213e-12	136872.41909986138	55423.962809686353	0.32877332260174519	1.189596298700324e-10
+5.8291457286432165	1.9006543222167855e-11	136201.45694057975	55325.899843339823	0.33039294153536136	2.3688541876486271e-10
+5.8693467336683414	3.7837691346393624e-11	135530.3005897892	55228.031068502256	0.33202907249650326	4.6842480959011603e-10
+5.9095477386934672	7.4817673498295838e-11	134858.88442348325	55130.422109180065	0.33368213145447062	9.198770560285676e-10
+5.949748743718593	1.4690935411705535e-10	134187.12270927755	55033.158697757623	0.33535259637017872	1.7932805553916426e-09
+5.9899497487437188	2.8643552745694476e-10	133514.90417003262	54936.352111374406	0.33704102384473894	3.4707135618297492e-09
+6.0301507537688446	5.5441752561765853e-10	132842.08531490254	54840.14584087635	0.33874806988558914	6.6660522042477553e-09
+6.0703517587939704	1.0652250401701472e-09	132168.48233071604	54744.723699434711	0.34047451560652425	1.2706336924493155e-08
+6.1105527638190953	2.0311009883348641e-09	131493.86131285713	54650.31959166551	0.34222129877936752	2.4026164210597341e-08
+6.1507537688442211	3.8429293846420915e-09	130817.9266067447	54557.22917214976	0.34398955225208327	4.5069231358142275e-08
+6.1909547738693469	7.2129790843121985e-09	130140.30703098714	54465.823622279197	0.34578065033522049	8.38299862792939e-08
+6.2311557788944727	1.3428683784404587e-08	129460.53976399067	54376.565763647486	0.3475962643291613	1.5461565441479818e-07
+6.2713567839195985	2.4790599494445707e-08	128778.0517001682	54290.028701841824	0.34943842841148703	2.8262765328727286e-07
+6.3115577889447234	4.537392599702378e-08	128092.13812297652	54206.917153405382	0.35130961711949227	5.1201024675162954e-07
+6.3517587939698492	8.2308657655036147e-08	127401.93860284236	54128.091547911412	0.35321283564044637	9.1875144999305768e-07
+6.391959798994975	1.4795025638593794e-07	126706.41011139398	54054.594913731649	0.35515172405593559	1.6328347684311823e-06
+6.4321608040201008	2.6342321683271918e-07	126004.29745158507	53987.682447912397	0.35713067657307851	2.8723898911136833e-06
+6.4723618090452266	4.6445619990930458e-07	125294.10123776048	53928.853536108865	0.35915497661463841	5.0006954540300621e-06
+6.5125628140703515	8.1059333904633736e-07	124574.04382082688	53879.885827414328	0.36123094843676162	8.6101613360336887e-06
+6.5527638190954773	1.3998338733327587e-06	123842.03374034319	53842.870782269878	0.36336612570777455	1.4657358290374727e-05
+6.5929648241206031	2.3908342088539787e-06	123095.62949461931	53820.249902365627	0.36556943723145763	2.4651133346090346e-05
+6.6331658291457289	4.0366362418005487e-06	122332.00364680415	53814.850624552622	0.36785140975801878	4.0939325569545934e-05
+6.6733668341708547	6.7332549412025123e-06	121547.90852212843	53829.920623600206	0.37022438762743104	6.7078390147623854e-05
+6.7135678391959797	1.1089056354883107e-05	120739.64498919026	53869.159030910247	0.37270276887122555	0.00010835056016530477
+6.7537688442211055	1.801755639017922e-05	119903.03604427798	53936.742850194387	0.37530325740361009	0.00017234643837799079
+6.7939698492462313	2.8857676853403598e-05	119033.40711785445	54037.346650989748	0.37804513110714955	0.00026964799652270645
+6.8341708542713571	4.5514755213871262e-05	118125.5751803381	54176.153462877977	0.3809505260084457	0.00041434482421663318
+6.8743718592964829	7.060963532400162e-05	117173.8488228944	54358.854694693538	0.38404473738860001	0.00062423514273949275
+6.9145728643216078	0.00010759670893942804	116172.0415104753	54591.636881484512	0.38735653961923638	0.00092005345618373217
+6.9547738693467336	0.00016078874097573756	115113.5001319235	54881.153134408145	0.39091852778717234	0.0013231517969031991
+6.9949748743718594	0.000235180297088849	113991.15079068157	55234.477350021938	0.39476748579047261	0.0018504899583136473
+7.0351758793969852	0.00033591149664039758	112797.56347789417	55659.039537181212	0.39894478756909496	0.002505688588844771
+7.075376884422111	0.00046721650664821086	111525.0368405421	56162.541048905106	0.40349684048381668	0.0032662121239443553
+7.1155778894472359	0.0006306609666981586	110165.7036997595	56752.849064059621	0.40847558258821559	0.0040656809437424491
+7.1557788944723617	0.00083078058556243257	108711.65729486763	57437.870343323317	0.41393904867021553	0.0049779755192488149
+7.1959798994974875	0.001074838683675616	107155.09744110543	58225.405071457382	0.4199520235118343	0.0060709451905654396
+7.2361809045226133	0.0013718541434935276	105488.49491617015	59122.982470764531	0.42658680490001022	0.0073882595629705542
+7.2763819095477391	0.0017317699546152028	103704.77146351588	60137.68079779065	0.43392410363520578	0.008952905801651671
+7.316582914572864	0.0021649060158569684	101797.49185760434	61275.935278074045	0.44205411330709982	0.010774259523388917
+7.3567839195979898	0.0026816707993767602	99761.063563127973	62543.338446922287	0.45107778919702851	0.012854523990054819
+7.3969849246231156	0.0032924898132355492	97590.93868629518	63944.43819812694	0.46110838368561985	0.015194122969737369
+7.4371859296482414	0.0040078737645461836	95283.81221325876	65482.53954553519	0.47227329548154084	0.017795175788852022
+7.4773869346733672	0.0048385570627486017	92837.810009375244	67159.51662379058	0.48471630250062631	0.02066324704278516
+7.5175879396984921	0.0057956561919645534	90252.659759627903	68975.641747909787	0.49860026418999276	0.023807840839246801
+7.557788944723618	0.0068908186183262955	87529.838003918048	70929.438377991508	0.51411039967862948	0.027242165355748343
+7.5979899497487438	0.0081363506546363994	84672.686689005757	73017.564567275636	0.53145827491314246	0.030982609403213835
+7.6381909547738696	0.0095453237481795614	81686.493234315232	75234.732896338028	0.55088666703954181	0.035048205701886145
+7.6783919597989954	0.011131663635046888	78578.52899694482	77573.672008080306	0.57267552058335969	0.039460204685824764
+7.7185929648241203	0.012910223992097831	75358.042159144621	80025.133720252372	0.59714927180521626	0.044241688881642198
+7.7587939698492461	0.014896828563412109	72036.202665341072	82577.948088427787	0.62468589868703539	0.049416788711442679
+7.7989949748743719	0.017108196121665177	68625.997376484855	85219.128251655842	0.65572817474882394	0.055007768011545034
+7.8391959798994977	0.019561377685956578	65142.0816182251	87934.018884287463	0.69079769761932419	0.061022891411748575
+7.8793969849246235	0.022271083468090636	61600.558577421354	90706.516799463076	0.73051285636383811	0.067403931330584713
+7.9195979899497493	0.025237340249817379	58018.872581856078	93519.177669400189	0.77560969383732226	0.073785637445452773
+7.9597989949748742	0.028385837618310192	54414.763499660039	96354.261625968094	0.82698144962591169	0.07831887204125873
+8	0.031257537688129805	50806.38688042853	99193.61311957147	0.8857154142038538	0.071433539236762961