1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137
|
#define brown1_N 20
#define brown1_P 4
static double brown1_x0[brown1_P] = { 25, 5, -5, -1 };
static double brown1_epsrel = 1.0e-5;
static double brown1_J[brown1_N * brown1_P];
static void
brown1_checksol(const double x[], const double sumsq,
const double epsrel, const char *sname,
const char *pname)
{
size_t i;
const double sumsq_exact = 8.582220162635628e+04;
const double brown1_x[brown1_P] = {
-1.159443990239263e+01, 1.320363005221244e+01,
-4.034395456782477e-01, 2.367789088597534e-01 };
gsl_test_rel(sumsq, sumsq_exact, epsrel, "%s/%s sumsq",
sname, pname);
for (i = 0; i < brown1_P; ++i)
{
gsl_test_rel(x[i], brown1_x[i], epsrel, "%s/%s i=%zu",
sname, pname, i);
}
}
static int
brown1_f (const gsl_vector * x, void *params, gsl_vector * f)
{
double x0 = gsl_vector_get (x, 0);
double x1 = gsl_vector_get (x, 1);
double x2 = gsl_vector_get (x, 2);
double x3 = gsl_vector_get (x, 3);
size_t i;
for (i = 0; i < brown1_N; i++)
{
double ti = 0.2 * (i + 1);
double ui = x0 + x1 * ti - exp (ti);
double vi = x2 + x3 * sin (ti) - cos (ti);
gsl_vector_set (f, i, ui * ui + vi * vi);
}
(void)params; /* avoid unused parameter warning */
return GSL_SUCCESS;
}
static int
brown1_df (CBLAS_TRANSPOSE_t TransJ, const gsl_vector * x,
const gsl_vector * u, void * params, gsl_vector * v,
gsl_matrix * JTJ)
{
gsl_matrix_view J = gsl_matrix_view_array(brown1_J, brown1_N, brown1_P);
double x0 = gsl_vector_get (x, 0);
double x1 = gsl_vector_get (x, 1);
double x2 = gsl_vector_get (x, 2);
double x3 = gsl_vector_get (x, 3);
size_t i;
for (i = 0; i < brown1_N; i++)
{
double ti = 0.2 * (i + 1);
double ui = x0 + x1 * ti - exp (ti);
double vi = x2 + x3 * sin (ti) - cos (ti);
gsl_matrix_set (&J.matrix, i, 0, 2 * ui);
gsl_matrix_set (&J.matrix, i, 1, 2 * ui * ti);
gsl_matrix_set (&J.matrix, i, 2, 2 * vi);
gsl_matrix_set (&J.matrix, i, 3, 2 * vi * sin (ti));
}
if (v)
gsl_blas_dgemv(TransJ, 1.0, &J.matrix, u, 0.0, v);
if (JTJ)
gsl_blas_dsyrk(CblasLower, CblasTrans, 1.0, &J.matrix, 0.0, JTJ);
(void)params; /* avoid unused parameter warning */
return GSL_SUCCESS;
}
static int
brown1_fvv (const gsl_vector * x, const gsl_vector * v,
void *params, gsl_vector * fvv)
{
double v0 = gsl_vector_get (v, 0);
double v1 = gsl_vector_get (v, 1);
double v2 = gsl_vector_get (v, 2);
double v3 = gsl_vector_get (v, 3);
size_t i;
for (i = 0; i < brown1_N; i++)
{
double ti = 0.2 * (i + 1);
double term1 = v0 + ti*v1;
double term2 = v3*sin(ti);
gsl_vector_set (fvv, i, 2.0 * (term1*term1 + v2*v2 +
term2 * (2*v2 + term2)));
}
(void)x; /* avoid unused parameter warning */
(void)params; /* avoid unused parameter warning */
return GSL_SUCCESS;
}
static gsl_multilarge_nlinear_fdf brown1_func =
{
brown1_f,
brown1_df,
brown1_fvv,
brown1_N,
brown1_P,
NULL,
0,
0,
0,
0
};
static test_fdf_problem brown1_problem =
{
"brown_dennis",
brown1_x0,
NULL,
&brown1_epsrel,
&brown1_checksol,
&brown1_func
};
|