From f4a4130dfe8e48c89e464e0724a28ec2cd041438 Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Fri, 27 Apr 2018 16:20:53 +0200 Subject: Simplify step size choice. --- init.c | 38 ++++++++++++-------------------------- 1 file changed, 12 insertions(+), 26 deletions(-) diff --git a/init.c b/init.c index ff995b6..c512cca 100644 --- a/init.c +++ b/init.c @@ -343,8 +343,7 @@ int td_solve(TDContext *td, double *coeffs_init[3]) const double delta = 1e-5; const double a1 = a0 - delta; - double cur_amplitude, new_amplitude, inverse_step; - int dir; + double cur_amplitude, new_amplitude, step; ret = constraint_eval_alloc(td, a0, &ce); if (ret < 0) @@ -382,16 +381,11 @@ int td_solve(TDContext *td, double *coeffs_init[3]) } cur_amplitude = a1; - dir = SGN(td->amplitude - cur_amplitude); - inverse_step = 1.0 / td->amplitude - 1.0 / cur_amplitude; - while (fabs(cur_amplitude - td->amplitude) > DBL_EPSILON) { - //double scale_factor; - //scale_factor = cur_amplitude; - //cur_amplitude = MIN(td->amplitude, 2 * cur_amplitude); - //scale_factor = cur_amplitude / scale_factor; + step = (td->amplitude - cur_amplitude) / 128.0; - new_amplitude = 1.0 / ((1.0 / cur_amplitude) + inverse_step); - if (dir * (td->amplitude - new_amplitude) < 0) + while (fabs(cur_amplitude - td->amplitude) > DBL_EPSILON) { + new_amplitude = cur_amplitude + step; + if (SGN(step) != SGN(td->amplitude - new_amplitude)) new_amplitude = td->amplitude; tdi_log(&s->logger, 2, "Trying amplitude %g\n", new_amplitude); @@ -404,27 +398,19 @@ int td_solve(TDContext *td, double *coeffs_init[3]) ce, s->coeffs_tmp, 1); tdi_constraint_eval_free(&ce); if (ret == -EDOM) { - inverse_step = 0.5 * inverse_step; - //if (fabs(inverse_step) < 1e-2) - // return ret; + step *= 0.5; + if (fabs(step) < 1e-5) + return ret; continue; } else if (ret < 0) return ret; + + if (ret <= s->solver->maxiter / 2) + step *= 1.75; + cur_amplitude = new_amplitude; memcpy(s->coeffs, s->coeffs_tmp, sizeof(*s->coeffs) * N); } - //while (1) { - - - // if (fabs(cur_amplitude - td->amplitude) < DBL_EPSILON) - // goto finish; - - // cblas_dscal(td->nb_coeffs[0] * td->nb_coeffs[1], SQR(scale_factor), td->coeffs[0], 1); - // cblas_dscal(td->nb_coeffs[0] * td->nb_coeffs[1], scale_factor, td->coeffs[1], 1); - // cblas_dscal(td->nb_coeffs[0] * td->nb_coeffs[1], scale_factor, td->coeffs[2], 1); - - //} - } finish: tdi_nlsolve_print_stats(s->solver); -- cgit v1.2.3