diff options
| -rw-r--r-- | unsupported/Eigen/src/NonLinear/fdjac1.h | 91 | ||||
| -rw-r--r-- | unsupported/Eigen/src/NonLinear/fdjac2.h | 25 |
2 files changed, 41 insertions, 75 deletions
diff --git a/unsupported/Eigen/src/NonLinear/fdjac1.h b/unsupported/Eigen/src/NonLinear/fdjac1.h index b9ec75bd4..e697e7ab6 100644 --- a/unsupported/Eigen/src/NonLinear/fdjac1.h +++ b/unsupported/Eigen/src/NonLinear/fdjac1.h @@ -22,68 +22,47 @@ int ei_fdjac1( eps = ei_sqrt((std::max(epsfcn,epsmch))); msum = ml + mu + 1; - if (msum < n) { - goto L40; - } - - /* computation of dense approximate jacobian. */ - - for (j = 0; j < n; ++j) { - temp = x[j]; - h = eps * ei_abs(temp); - if (h == 0.) - h = eps; - x[j] = temp + h; - iflag = Functor::f(x, wa1); - if (iflag < 0) - goto L30; - x[j] = temp; - for (i = 0; i < n; ++i) { - fjac(i,j) = (wa1[i] - fvec[i]) / h; - /* L10: */ + if (msum >= n) { + /* computation of dense approximate jacobian. */ + for (j = 0; j < n; ++j) { + temp = x[j]; + h = eps * ei_abs(temp); + if (h == 0.) + h = eps; + x[j] = temp + h; + iflag = Functor::f(x, wa1); + if (iflag < 0) + return iflag; + x[j] = temp; + fjac.col(j) = (wa1-fvec)/h; } - /* L20: */ - } -L30: - /* goto L110; */ - return iflag; -L40: - - /* computation of banded approximate jacobian. */ - for (k = 0; k < msum; ++k) { - for (j = k; msum< 0 ? j > n: j < n; j += msum) { - wa2[j] = x[j]; - h = eps * ei_abs(wa2[j]); - if (h == 0.) h = eps; - x[j] = wa2[j] + h; - /* L60: */ - } - iflag = Functor::f(x, wa1); - if (iflag < 0) { - /* goto L100; */ - return iflag; - } - for (j = k; msum< 0 ? j > n: j < n; j += msum) { - x[j] = wa2[j]; - h = eps * ei_abs(wa2[j]); - if (h == 0.) h = eps; - for (i = 0; i < n; ++i) { - fjac(i,j) = 0.; - if (i >= j - mu && i <= j + ml) { - fjac(i,j) = (wa1[i] - fvec[i]) / h; + }else { + /* computation of banded approximate jacobian. */ + for (k = 0; k < msum; ++k) { + for (j = k; msum< 0 ? j > n: j < n; j += msum) { + wa2[j] = x[j]; + h = eps * ei_abs(wa2[j]); + if (h == 0.) h = eps; + x[j] = wa2[j] + h; + } + iflag = Functor::f(x, wa1); + if (iflag < 0) { + return iflag; + } + for (j = k; msum< 0 ? j > n: j < n; j += msum) { + x[j] = wa2[j]; + h = eps * ei_abs(wa2[j]); + if (h == 0.) h = eps; + for (i = 0; i < n; ++i) { + fjac(i,j) = 0.; + if (i >= j - mu && i <= j + ml) { + fjac(i,j) = (wa1[i] - fvec[i]) / h; + } } - /* L70: */ } - /* L80: */ } - /* L90: */ } - /* L100: */ - /* L110: */ return iflag; - - /* last card of subroutine fdjac1. */ - } /* fdjac1_ */ diff --git a/unsupported/Eigen/src/NonLinear/fdjac2.h b/unsupported/Eigen/src/NonLinear/fdjac2.h index fd6cf1d8b..f2e627d4f 100644 --- a/unsupported/Eigen/src/NonLinear/fdjac2.h +++ b/unsupported/Eigen/src/NonLinear/fdjac2.h @@ -8,18 +8,15 @@ int ei_fdjac2( Matrix< Scalar, Dynamic, 1 > &wa) { /* Local variables */ - Scalar h; - int i, j; - Scalar eps, temp; + Scalar h, temp; int iflag; /* Function Body */ const Scalar epsmch = epsilon<Scalar>(); const int n = x.size(); - const int m = fvec.size(); + const Scalar eps = ei_sqrt((std::max(epsfcn,epsmch))); - eps = ei_sqrt((std::max(epsfcn,epsmch))); - for (j = 0; j < n; ++j) { + for (int j = 0; j < n; ++j) { temp = x[j]; h = eps * ei_abs(temp); if (h == 0.) { @@ -27,21 +24,11 @@ int ei_fdjac2( } x[j] = temp + h; iflag = Functor::f(x, wa); - if (iflag < 0) { - /* goto L30; */ + if (iflag < 0) return iflag; - } x[j] = temp; - for (i = 0; i < m; ++i) { - fjac(i,j) = (wa[i] - fvec[i]) / h; - /* L10: */ - } - /* L20: */ + fjac.col(j) = (wa-fvec)/h; } - /* L30: */ return iflag; - - /* last card of subroutine fdjac2. */ - -} /* fdjac2_ */ +} |