GitHub-Proxy/eigen
1
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2025-08-13 04:09:10 +08:00
Code Issues Packages Projects Releases Wiki Activity

beautify functors for lmdif, lmstr, hybrj, hybrd

Browse Source
This commit is contained in:
Thomas Capricelli 2009-08-23 04:57:48 +02:00
parent acd757737a
commit 134dea76d3
7 changed files with 76 additions and 81 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
unsupported/Eigen/src/NonLinear/fdjac1.h
View File

@ -34,7 +34,7 @@ int ei_fdjac1(
if (h == 0.) if (h == 0.)
h = eps; h = eps;
x[j] = temp + h; x[j] = temp + h;
iflag = Functor::f(x, wa1, 1); iflag = Functor::f(x, wa1);
if (iflag < 0) if (iflag < 0)
goto L30; goto L30;
x[j] = temp; x[j] = temp;
@ -59,7 +59,7 @@ L40:
x[j] = wa2[j] + h; x[j] = wa2[j] + h;
/* L60: */ /* L60: */
} }
iflag = Functor::f(x, wa1, 1); iflag = Functor::f(x, wa1);
if (iflag < 0) { if (iflag < 0) {
/* goto L100; */ /* goto L100; */
return iflag; return iflag;

2
unsupported/Eigen/src/NonLinear/fdjac2.h
View File

@ -26,7 +26,7 @@ int ei_fdjac2(
h = eps; h = eps;
} }
x[j] = temp + h; x[j] = temp + h;
iflag = Functor::f(x, wa, 1); iflag = Functor::f(x, wa);
if (iflag < 0) { if (iflag < 0) {
/* goto L30; */ /* goto L30; */
return iflag; return iflag;

8
unsupported/Eigen/src/NonLinear/hybrd.h
View File

@ -75,7 +75,7 @@ L20:
/* evaluate the function at the starting point */ /* evaluate the function at the starting point */
/* and calculate its norm. */ /* and calculate its norm. */
iflag = Functor::f(x, fvec, 1); iflag = Functor::f(x, fvec);
nfev = 1; nfev = 1;
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -215,7 +215,7 @@ L180:
} }
iflag = 0; iflag = 0;
if ((iter - 1) % nprint == 0) { if ((iter - 1) % nprint == 0) {
iflag = Functor::f(x, fvec, 0); iflag = Functor::debug(x, fvec);
} }
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -244,7 +244,7 @@ L190:
/* evaluate the function at x + p and calculate its norm. */ /* evaluate the function at x + p and calculate its norm. */
iflag = Functor::f(wa2, wa4, 0); iflag = Functor::f(wa2, wa4);
++nfev; ++nfev;
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -404,7 +404,7 @@ L300:
info = iflag; info = iflag;
} }
if (nprint > 0) { if (nprint > 0) {
iflag = Functor::f(x, fvec, 0); iflag = Functor::debug(x, fvec);
} }
return info; return info;

10
unsupported/Eigen/src/NonLinear/hybrj.h
View File

@ -70,7 +70,7 @@ L20:
/* evaluate the function at the starting point */ /* evaluate the function at the starting point */
/* and calculate its norm. */ /* and calculate its norm. */
iflag = Functor::f(x, fvec, fjac, 1); iflag = Functor::f(x, fvec);
nfev = 1; nfev = 1;
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -92,7 +92,7 @@ L30:
/* calculate the jacobian matrix. */ /* calculate the jacobian matrix. */
iflag = Functor::f(x, fvec, fjac, 2); iflag = Functor::df(x, fjac);
++njev; ++njev;
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -203,7 +203,7 @@ L180:
} }
iflag = 0; iflag = 0;
if ((iter - 1) % nprint == 0) if ((iter - 1) % nprint == 0)
iflag = Functor::f(x, fvec, fjac, 0); iflag = Functor::debug(x, fvec, fjac);
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
} }
@ -231,7 +231,7 @@ L190:
/* evaluate the function at x + p and calculate its norm. */ /* evaluate the function at x + p and calculate its norm. */
iflag = Functor::f(wa2, wa4, fjac, 1); iflag = Functor::f(wa2, wa4);
++nfev; ++nfev;
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -395,7 +395,7 @@ L300:
info = iflag; info = iflag;
} }
if (nprint > 0) if (nprint > 0)
iflag = Functor::f(x, fvec, fjac, 0); iflag = Functor::debug(x, fvec, fjac);
return info; return info;
/* last card of subroutine hybrj. */ /* last card of subroutine hybrj. */

8
unsupported/Eigen/src/NonLinear/lmdif.h
View File

@ -66,7 +66,7 @@ L20:
/* evaluate the function at the starting point */ /* evaluate the function at the starting point */
/* and calculate its norm. */ /* and calculate its norm. */
iflag = Functor::f(x, fvec, 1); iflag = Functor::f(x, fvec);
nfev = 1; nfev = 1;
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -97,7 +97,7 @@ L30:
} }
iflag = 0; iflag = 0;
if ((iter - 1) % nprint == 0) { if ((iter - 1) % nprint == 0) {
iflag = Functor::f(x, fvec, 0); iflag = Functor::debug(x, fvec);
} }
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -239,7 +239,7 @@ L200:
/* evaluate the function at x + p and calculate its norm. */ /* evaluate the function at x + p and calculate its norm. */
iflag = Functor::f(wa2, wa4, 1); iflag = Functor::f(wa2, wa4);
++nfev; ++nfev;
if (iflag < 0) { if (iflag < 0) {
goto L300; goto L300;
@ -380,7 +380,7 @@ L300:
} }
iflag = 0; iflag = 0;
if (nprint > 0) { if (nprint > 0) {
iflag = Functor::f(x, fvec, 0); iflag = Functor::debug(x, fvec);
} }
return info; return info;

10
unsupported/Eigen/src/NonLinear/lmstr.h
View File

@ -66,7 +66,7 @@ L20:
/* evaluate the function at the starting point */ /* evaluate the function at the starting point */
/* and calculate its norm. */ /* and calculate its norm. */
iflag = Functor::f(x, fvec, wa3, 1); iflag = Functor::f(x, fvec);
nfev = 1; nfev = 1;
if (iflag < 0) { if (iflag < 0) {
goto L340; goto L340;
@ -89,7 +89,7 @@ L30:
} }
iflag = 0; iflag = 0;
if ((iter - 1) % nprint == 0) { if ((iter - 1) % nprint == 0) {
iflag = Functor::f(x, fvec, wa3, 0); iflag = Functor::debug(x, fvec, wa3);
} }
if (iflag < 0) { if (iflag < 0) {
goto L340; goto L340;
@ -111,7 +111,7 @@ L40:
} }
iflag = 2; iflag = 2;
for (i = 0; i < m; ++i) { for (i = 0; i < m; ++i) {
if (Functor::f(x, fvec, wa3, iflag) < 0) { if (Functor::df(x, wa3, iflag) < 0) {
goto L340; goto L340;
} }
temp = fvec[i]; temp = fvec[i];
@ -264,7 +264,7 @@ L240:
/* evaluate the function at x + p and calculate its norm. */ /* evaluate the function at x + p and calculate its norm. */
iflag = Functor::f(wa2, wa4, wa3, 1); iflag = Functor::f(wa2, wa4);
++nfev; ++nfev;
if (iflag < 0) { if (iflag < 0) {
goto L340; goto L340;
@ -406,7 +406,7 @@ L340:
} }
iflag = 0; iflag = 0;
if (nprint > 0) { if (nprint > 0) {
iflag = Functor::f(x, fvec, wa3, 0); iflag = Functor::debug(x, fvec, wa3);
} }
return info; return info;

115
unsupported/test/NonLinear.cpp
View File

@ -214,40 +214,36 @@ void testLmder()
} }
struct hybrj_functor { struct hybrj_functor {
static int f(const VectorXd &x, VectorXd &fvec, MatrixXd &fjac, int iflag)
static int debug(const VectorXd & /* x */, const VectorXd & /* fvec */, const MatrixXd & /* fjac */) { return 0;}
static int f(const VectorXd &x, VectorXd &fvec)
{ {
/* subroutine fcn for hybrj1 example. */ double temp, temp1, temp2;
int j, k;
double one=1, temp, temp1, temp2, three=3, two=2, zero=0, four=4;
const int n = x.size(); const int n = x.size();
assert(fvec.size()==n); assert(fvec.size()==n);
for (int k = 0; k < n; k++)
{
temp = (3. - 2.*x[k])*x[k];
temp1 = 0.;
if (k) temp1 = x[k-1];
temp2 = 0.;
if (k != n-1) temp2 = x[k+1];
fvec[k] = temp - temp1 - 2.*temp2 + 1.;
}
return 0;
}
static int df(const VectorXd &x, MatrixXd &fjac)
{
const int n = x.size();
assert(fjac.rows()==n); assert(fjac.rows()==n);
assert(fjac.cols()==n); assert(fjac.cols()==n);
for (int k = 0; k < n; k++)
if (iflag != 2)
{ {
for (k = 0; k < n; k++) for (int j = 0; j < n; j++)
{ fjac(k,j) = 0.;
temp = (three - two*x[k])*x[k]; fjac(k,k) = 3.- 4.*x[k];
temp1 = zero; if (k) fjac(k,k-1) = -1.;
if (k) temp1 = x[k-1]; if (k != n-1) fjac(k,k+1) = -2.;
temp2 = zero;
if (k != n-1) temp2 = x[k+1];
fvec[k] = temp - temp1 - two*temp2 + one;
}
}
else
{
for (k = 0; k < n; k++)
{
for (j = 0; j < n; j++)
fjac(k,j) = zero;
fjac(k,k) = three - four*x[k];
if (k) fjac(k,k-1) = -one;
if (k != n-1) fjac(k,k+1) = -two;
}
} }
return 0; return 0;
} }
@ -319,24 +315,21 @@ void testHybrj()
} }
struct hybrd_functor { struct hybrd_functor {
static int f(const VectorXd &x, VectorXd &fvec, int /*iflag*/) static int debug(const VectorXd & /* x */, const VectorXd & /* fvec */) { return 0;}
static int f(const VectorXd &x, VectorXd &fvec)
{ {
/* subroutine fcn for hybrd1 example. */ double temp, temp1, temp2;
int k;
double one=1, temp, temp1, temp2, three=3, two=2, zero=0;
const int n = x.size(); const int n = x.size();
assert(fvec.size()==n); assert(fvec.size()==n);
for (int k=0; k < n; k++)
for (k=0; k < n; k++)
{ {
temp = (three - two*x[k])*x[k]; temp = (3. - 2.*x[k])*x[k];
temp1 = zero; temp1 = 0.;
if (k) temp1 = x[k-1]; if (k) temp1 = x[k-1];
temp2 = zero; temp2 = 0.;
if (k != n-1) temp2 = x[k+1]; if (k != n-1) temp2 = x[k+1];
fvec[k] = temp - temp1 - two*temp2 + one; fvec[k] = temp - temp1 - 2.*temp2 + 1.;
} }
return 0; return 0;
} }
@ -400,41 +393,42 @@ void testHybrd()
} }
struct lmstr_functor { struct lmstr_functor {
static int f(const VectorXd &x, VectorXd &fvec, VectorXd &fjrow, int iflag) static int debug(const VectorXd & /* x */, const VectorXd & /* fvec */, const VectorXd & /* fjac */) { return 0;}
static int f(const VectorXd &x, VectorXd &fvec)
{ {
/* subroutine fcn for lmstr1 example. */ /* subroutine fcn for lmstr1 example. */
int i; double tmp1, tmp2, tmp3;
double tmp1, tmp2, tmp3, tmp4;
double y[15]={1.4e-1, 1.8e-1, 2.2e-1, 2.5e-1, 2.9e-1, 3.2e-1, 3.5e-1, double y[15]={1.4e-1, 1.8e-1, 2.2e-1, 2.5e-1, 2.9e-1, 3.2e-1, 3.5e-1,
3.9e-1, 3.7e-1, 5.8e-1, 7.3e-1, 9.6e-1, 1.34, 2.1, 4.39}; 3.9e-1, 3.7e-1, 5.8e-1, 7.3e-1, 9.6e-1, 1.34, 2.1, 4.39};
assert(15==fvec.size()); assert(15==fvec.size());
assert(3==x.size()); assert(3==x.size());
assert(fjrow.size()==x.size());
if (iflag < 2) for (int i=0; i<15; i++)
{ {
for (i=0; i<15; i++)
{
tmp1 = i+1;
tmp2 = 16 - i - 1;
tmp3 = (i>=8)? tmp2 : tmp1;
fvec[i] = y[i] - (x[0] + tmp1/(x[1]*tmp2 + x[2]*tmp3));
}
}
else
{
i = iflag-2;
tmp1 = i+1; tmp1 = i+1;
tmp2 = 16 - i - 1; tmp2 = 16 - i - 1;
tmp3 = (i>=8)? tmp2 : tmp1; tmp3 = (i>=8)? tmp2 : tmp1;
tmp4 = (x[1]*tmp2 + x[2]*tmp3); tmp4 = tmp4*tmp4; fvec[i] = y[i] - (x[0] + tmp1/(x[1]*tmp2 + x[2]*tmp3));
fjrow[0] = -1;
fjrow[1] = tmp1*tmp2/tmp4;
fjrow[2] = tmp1*tmp3/tmp4;
} }
return 0; return 0;
} }
static int df(const VectorXd &x, VectorXd &jac_row, int rownb)
{
assert(x.size()==3);
assert(jac_row.size()==x.size());
double tmp1, tmp2, tmp3, tmp4;
int i = rownb-2;
tmp1 = i+1;
tmp2 = 16 - i - 1;
tmp3 = (i>=8)? tmp2 : tmp1;
tmp4 = (x[1]*tmp2 + x[2]*tmp3); tmp4 = tmp4*tmp4;
jac_row[0] = -1;
jac_row[1] = tmp1*tmp2/tmp4;
jac_row[2] = tmp1*tmp3/tmp4;
return 0;
}
}; };
void testLmstr1() void testLmstr1()
@ -495,7 +489,8 @@ void testLmstr()
} }
struct lmdif_functor { struct lmdif_functor {
static int f(const VectorXd &x, VectorXd &fvec, int /*iflag*/) static int debug(const VectorXd & /* x */, const VectorXd & /* fvec */) { return 0;}
static int f(const VectorXd &x, VectorXd &fvec)
{ {
/* function fcn for lmdif1 example */ /* function fcn for lmdif1 example */