[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: Undefined behavior error in Nonlinear Least-Squares Fitting
|
From: |
Patrick Alken |
|
Subject: |
Re: Undefined behavior error in Nonlinear Least-Squares Fitting |
|
Date: |
Sun, 13 Dec 2020 15:03:14 -0700 |
|
User-agent: |
Mozilla/5.0 (X11; Linux x86_64; rv:68.0) Gecko/20100101 Thunderbird/68.10.0 |
Hello,
Thanks for your report. The problem is that on the first iteration,
the lambda value jumps to -12.8, and the value of exp(12.8 * 60) is too
large to represent in double precision, so infs and nans appear in the
residual vector and Jacobian matrix.
For this particular problem, the value of lambda should be > 0, however
there is no constraint placed to assure that, so when the first
iteration jumps to -12.8, bad numerical things happen.
I think the best solution for now is to adjust the model parameters a
little to prevent this type of behavior, and bring TMAX closer to 1
(instead of 40). A better solution would be to implement a constraint so
that lambda > 0, but this would require a lot of development effort.
I will modify the test program to make it less likely to run into this
problem in the future.
Thanks,
Patrick
On 12/7/20 3:18 PM, Zhoulai Fu@Gmail wrote:
> Running the attached program nlfit_bug.c triggers undefined behavior,
> calculating unreasonable
> results. The program "nlfit_debug.c" is adapted from the example program
> provided
> by GSL's doc/examples/nlfit.c. Their difference lies in only one parameter
> 40.0
> versus 60.0:
>
> 1d0
> <
> 12c11,12
> < #define TMAX (60.0) /* time variable in [0,TMAX] */ //40.0 instead of
> 60.0 in nlfit.c
> ---
>> //#define TMAX (40.0) /* time variable in [0,TMAX] */
>> #define TMAX (60.0) /* time variable in [0,TMAX] */
>
> The undefined behavior can be better visualized if GCC's Undefined
> Behavior Sanitizer is enabled (-fsanitize=undefined). Below is the
> execution results:
>
> data: 59.3939 1.03153 0.101317
> data: 60 1.14181 0.101239
> iter 0: A = 1.0000, lambda = 1.0000, b = 0.0000, cond(J) = inf,
> |f(x)| = 101.0200
> iter 1: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = 92.8216,
> |f(x)| = -nan
> iter 2: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 3: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 4: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
>
> *nielsen.c:98:7: runtime error: left shift of 4611686018427387904 by 1
> places cannot be represented in type 'long int'*iter 5: A = 3.5110, lambda
> = -12.8820, b = 1.2364, cond(J) = nan, |f(x)| = -nan
> iter 6: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 7: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 8: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 9: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 10: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 11: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 12: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 13: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 14: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 15: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 16: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 17: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 18: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 19: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 20: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 21: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 22: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 23: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 24: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 25: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 26: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 27: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 28: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 29: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 30: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 31: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 32: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 33: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 34: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 35: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 36: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 37: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 38: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 39: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 40: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 41: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 42: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 43: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 44: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 45: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 46: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 47: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 48: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 49: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 50: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 51: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 52: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 53: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 54: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 55: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 56: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 57: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 58: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 59: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 60: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 61: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 62: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 63: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 64: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 65: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 66: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 67: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 68: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 69: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 70: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 71: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 72: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 73: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 74: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 75: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 76: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 77: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 78: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 79: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 80: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 81: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 82: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 83: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 84: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 85: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 86: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 87: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 88: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 89: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 90: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 91: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 92: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 93: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 94: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 95: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 96: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 97: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 98: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 99: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> iter 100: A = 3.5110, lambda = -12.8820, b = 1.2364, cond(J) = nan,
> |f(x)| = -nan
> summary from method 'trust-region/levenberg-marquardt'
> number of iterations: 100
> function evaluations: 1586
> Jacobian evaluations: 2
> reason for stopping: small gradient
> initial |f(x)| = 101.020000
> final |f(x)| = inf
> chisq/dof = inf
> A = 3.51103 +/- nan
> lambda = -12.88201 +/- nan
> b = 1.23642 +/- nan
> status = exceeded max number of iterations
>
> Preliminary debugging: The undefined behavior occurs in
> multifit_nlinear/nielsen.c, line 98, triggered by overflow of a left shift
> of a long int.
>
> *nu <<= 1;
>
> BR,
> Zhoulai