topic In general, nonlinear in Intel® oneAPI Math Kernel Library

dtrnlspbc_solve yields NaN

Dimitris_Dakopoulos — Thu, 19 Sep 2013 13:31:39 GMT

Dear all,

For my application I use the bounded nonlinear solver almost identical to the examples provided in MKL.

I show sample of my code (unfortunately it's not possible to presend a more detailed compilable version).

So, in the main solver loop, after a couple of iterations I notice that after running dtrnlsbc_solve optimization_vector

is all nans. I am certain that no-one else messes with it in the meantime.

Also, fvec and fjac do not have any nans, infs or very small or large values. Does anybody have any idea why could this

happen? Maybe an improper initialization?

Note also, thatthis problematic case is unique. It runs fine in all other cases.

Thank you in advance

Dimitris

System: Linux CentOS 4, gcc-4.7.2, MKL 11.0

[cpp]

std::vector<double> eps = {1e-5, 1e-5, 1e-5, 1e-5, 1e-5, 1e-5};
std::vector<double> fvec(m, 0.);
std::vector<double> fjac(m * n, 0.);
_TRNSP_HANDLE_t handle;

int iter1 = 1000, iter2 = 100;
double rs = 0.0;

if(dtrnlspbc_init(&handle, &n, &m,
                  optimization_vector, lower_bounds
                  upper_bounds, eps.data(), &iter1, &iter2, &rs)
   != TR_SUCCESS)
{
    exit(1);
}

int info[6];
if(dtrnlspbc_check(&handle, &n, &m, fjac.data(), fvec.data(),
                   lower_bounds, upper_bounds, eps.data(),
                   info) != TR_SUCCESS)
{
    exit(1);
}
else if(info[0] || info[1] || info[2] || info[3])
{
    exit(1);
}

int RCI_Request = 0;

// main solver loop
while(1)
{
    // THE PROBLEM IS HERE
    // I print optimization_vector, fvec, fjac and they look fine but afterwards

    // optimization_vector is made into nans
    if(dtrnlspbc_solve(&handle, fvec.data(), fjac.data(), &RCI_Request)
       != TR_SUCCESS)
    {
        exit(1);
    }

    if(RCI_Request == -1 || RCI_Request == -2 || RCI_Request == -3
       || RCI_Request == -4 || RCI_Request == -5 || RCI_Request == -6)
    {
        break; // success!
    }
    else if(RCI_Request == 1)
    {
        objective(&m, &n, optimization_vector, fvec.data(), this);
    }
    else if(RCI_Request == 2)
    {
        if(djacobix(objective, &n, &m, fjac.data(),
                    optimization_vector, eps.data(), this)
           != TR_SUCCESS)
        {
            exit(1);
        }
    }
}

double r1, r2;
int iter, st_cr;
if(dtrnlspbc_get(&handle, &iter, &st_cr, &r1, &r2) != TR_SUCCESS)
{
    exit(1);
}

if(dtrnlspbc_delete(&handle) != TR_SUCCESS)
{
    exit(1);
}

mkl_free_buffers();

[/cpp]

You have left out important

mecej4 — Fri, 20 Sep 2013 12:41:54 GMT

You have left out important information, such as declarations and #include statements, and it would be useful to see complete source code for the failed run. Without that information, you are essentially asking us to give our opinions on your speculations.

Exactly. I would be glad if I

Dimitris_Dakopoulos — Fri, 20 Sep 2013 12:53:24 GMT

Exactly. I would be glad if I got comments:

* on the structure of the flow I am following

* by users having dealt with something similar

* theoritical hints on whether the solver is possible to fail and in which cases.

Unfortunately I cannot strip down the case from my application and/or provide more "usable" code.

UPDATE: After lots of debugging seems that there is a problem with my objective function

UPDATE: For some reason, the

Dimitris_Dakopoulos — Fri, 20 Sep 2013 13:49:06 GMT

UPDATE:
For some reason, the solver produces a solution vector value outside the boundaries I provided (0, inf), more precicely a negative one. My objective function contained a square root of that so it produced NaN.

Now, I guess, there is another question. Why did this happen?

In general, nonlinear

mecej4 — Fri, 20 Sep 2013 17:42:07 GMT

In general, nonlinear constrained optimization methods do not visit only feasible points on the way to the optimum.

If you displayed the code for your objective function, at least, perhaps some progress could be made. Keep in mind that MKL is not open source, and giving answers to the questions you asked about the methodology may involve divulging proprietary information.

I think that it would be more productive to ensure that you have formulated the problem correctly and are invoking the MKL routines properly.