Community
cancel
Showing results for 
Search instead for 
Did you mean: 
Sh__Rustem
Beginner
80 Views

dtrnlspbc_init fails with TR_INVALID_OPTION

It appears that the trust region optimizer can only handle cases when number of function arguments is equal to number of function values (square Jacobian).

In particular, in the code below (this is a slightly modified fragment of optimization example from the MKL library package), the initialization function  dtrnlspbc_init fails (returns TR_INVALID_OPTION(1502)) whenever m != n. In particular, if n = 3, m = 1 (a scalar function of 3 variables). Note, that this initialization function does not know anything about objective function, it succeeds whenever m = n, like m = 3, n =3 or m = 5, n = 5.  It fails if n != n. 

Is this expected behavior?

If it is not, what am I doing wrong?

Thanks.  

 

int main()
{
    /* user’s objective function */
    extern void extended_powell(MKL_INT *, MKL_INT *, double *, double *, void *);
    /* n - number of function variables
    m - dimension of function value */
    MKL_INT n = 3, m = 1;
    /* precisions for stop-criteria (see manual for more details) */
    double eps[6];
    /* precision of the Jacobian matrix calculation */
    double jac_eps;
    /* solution vector. contains values x for f(x) */
    double *x = NULL;
    /* iter1 - maximum number of iterations
    iter2 - maximum number of iterations of calculation of trial-step */
    MKL_INT iter1 = 1000, iter2 = 100;
    /* initial step bound */
    double rs = 0.0;
    /* reverse communication interface parameter */
    MKL_INT RCI_Request;
    /* controls of rci cycle */
    MKL_INT successful;
    /* function (f(x)) value vector */
    double *fvec = NULL;
    /* jacobi matrix */
    double *fjac = NULL;
    /* lower and upper bounds */
    double *LW = NULL, *UP = NULL;
    /* number of iterations */
    MKL_INT iter;
    /* number of stop-criterion */
    MKL_INT st_cr;
    /* initial and final residuals */
    double r1, r2;
    /* TR solver handle */
    _TRNSPBC_HANDLE_t handle;
    /* cycle’s counter */
    MKL_INT i;
    /* results of input parameter checking */
    MKL_INT info[6];
    /* memory allocation flags */
    MKL_INT mem_error, error;

    /*Additional users data */
    u_data m_data;
    m_data.a = 1;
    m_data.sum = 0;

    error = 0;
    /* memory allocation */
    mem_error = 1;
    x = (double *)mkl_malloc(sizeof(double) * n, 64);
    if (x == NULL) goto end;
    fvec = (double *)mkl_malloc(sizeof(double) * m, 64);
    if (fvec == NULL) goto end;
    fjac = (double *)mkl_malloc(sizeof(double) * m * n, 64);
    if (fjac == NULL) goto end;
    LW = (double *)mkl_malloc(sizeof(double) * n, 64);
    if (LW == NULL) goto end;
    UP = (double *)mkl_malloc(sizeof(double) * n, 64);
    if (UP == NULL) goto end;
    /* memory allocated correctly */
    mem_error = 0;
    /* set precisions for stop-criteria */
    for (i = 0; i < 6; i++)
    {
        eps = 0.00001;
    }
    /* set precision of the Jacobian matrix calculation */
    jac_eps = 0.00000001;
    /* set the initial guess */
    for (i = 0; i < n; i++)
        x = 0;
    /* set initial values */
    for (i = 0; i < m; i++)
        fvec = 0.0;
    for (i = 0; i < m * n; i++)
        fjac = 0.0;
    /* set bounds */
    for (i = 0; i < n; i++)
    {
        LW = -1.;
        UP = 1.0;
    }

    /* initialize solver (allocate memory, set initial values)
    handle       in/out: TR solver handle
    n       in:     number of function variables
    m       in:     dimension of function value
    x       in:     solution vector. contains values x for f(x)
    LW           in:             lower bound
    UP           in:             upper bound
    eps     in:     precisions for stop-criteria
    iter1   in:     maximum number of iterations
    iter2   in:     maximum number of iterations of calculation of trial-step
    rs      in:     initial step bound */
    MKL_INT st = dtrnlspbc_init(&handle, &n, &m, x, LW, UP, eps, &iter1, &iter2, &rs);

0 Kudos
1 Reply
mecej4
Black Belt
80 Views

I have used DTRNLSP for dozens of NLS (nonlinear least squares) problems where the number of observations exceeds the number of model coefficients. For example, the NIST Hahn1 test problem, with m = 236 observations and n = 7 coefficients.

Note, however, that the trust region method is designed for minimizing F(x), where the dimension m of F is ≥ the dimension n of x. See, for example, https://software.intel.com/en-us/mkl-developer-reference-c-nonlinear-least-squares-problem-without-constraints#4EA75446-AA99-423C-B497-00D86E1CDEFE . Thus, if you set m = 1 and n = 3, this condition is violated.

If you want help finding more about this or related bugs, please provide the complete source code as an attachment (include data files, if any, and instructions to build, and state the compiler version and options used to build).

Reply