Solved: dtrnlsp example won't run

pivello · ‎02-10-2011

Hi all. I'm trying to use dtrnlsp but it fails even for the example provided by the documentation, as dtrnlsp_init returns TR_OUT_OF_MEMORY error.

Linking was done according to intel mkl link line advisor:

MKL_FLAGS = -L$(MKL_LIB_DIR) $(MKL_LIB_DIR)/libmkl_solver_ilp64_sequential.a -Wl,--start-group -lmkl_intel_ilp64 -lmkl_sequential -lmkl_core -Wl,--end-group -lpthread

I'm using MKL 10.3 and the compilers provided by intel parallel studio.

Any hints?

Thanks in advance

Mrcio

PS: The code is inserted bellow, just in case.

[cpp]/*
********************************************************************************
*                              INTEL CONFIDENTIAL
*   Copyright(C) 2004-2010 Intel Corporation. All Rights Reserved.
*   The source code contained  or  described herein and all documents related to
*   the source code ("Material") are owned by Intel Corporation or its suppliers
*   or licensors.  Title to the  Material remains with  Intel Corporation or its
*   suppliers and licensors. The Material contains trade secrets and proprietary
*   and  confidential  information of  Intel or its suppliers and licensors. The
*   Material  is  protected  by  worldwide  copyright  and trade secret laws and
*   treaty  provisions. No part of the Material may be used, copied, reproduced,
*   modified, published, uploaded, posted, transmitted, distributed or disclosed
*   in any way without Intel's prior express written permission.
*   No license  under any  patent, copyright, trade secret or other intellectual
*   property right is granted to or conferred upon you by disclosure or delivery
*   of the Materials,  either expressly, by implication, inducement, estoppel or
*   otherwise.  Any  license  under  such  intellectual property  rights must be
*   express and approved by Intel in writing.
*
********************************************************************************
*   Content : TR Solver C example
*
********************************************************************************
*/

#include 
#include 
#include 

#include 
#include 
#include 

/* nonlinear least square problem without boundary constraints */
int main ()
{
	/* users objective function */
	extern void extendet_powell (MKL_INT *, MKL_INT *, double*, double*);
	/* n - number of function variables
	   m - dimension of function value */
	MKL_INT			n = 4, m = 4;
	/* precisions for stop-criteria (see manual for more detailes) */
	double	eps[6];
	/* solution vector. contains values x for f(x) */
	double	*x;
	/* 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; // reverse communication interface variable
	/* controls of rci cycle */
	MKL_INT			successful;
	/* function (f(x)) value vector */
	double *fvec;
	/* jacobi matrix */
	double *fjac;
	/* number of iterations */
	MKL_INT			iter;
	/* number of stop-criterion */
	MKL_INT			st_cr;
	/* initial and final residauls */
	double r1, r2;
	/* TR solver handle */
	_TRNSP_HANDLE_t handle; // TR solver handle
	/* cycles counter */
	MKL_INT i, ierror;

	/* memory allocation */
	x = (double*) malloc (sizeof (double)*n);
	fvec = (double*) malloc (sizeof (double)*m);
	fjac = (double*) malloc (sizeof (double)*m*n);
	/* set precisions for stop-criteria */
	for (i = 0; i < 6; i++)
	{
		eps  = 0.00001;
	}
	/* set the initial guess */
	for (i = 0; i < n/4; i++)
	{
        x [4*i]		=  3.0;
        x [4*i + 1] = -1.0;
        x [4*i + 2] =  0.0;
        x [4*i + 3] =  1.0;
	}
	/* set initial values */
	for (i = 0; i < m; i++)
		fvec  = 0.0;
	for (i = 0; i < m*n; i++)
		fjac  = 0.0;
	/* initialize solver (allocate mamory, 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)
		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 */
	if ((ierror = dtrnlsp_init (&handle, &n, &m, x, eps, &iter1, &iter2, &rs)) != TR_SUCCESS)
	{
		/* if function does not complete successful then print error message */
		printf ("| error in dtrnlsp_init\n");
      if (ierror == TR_INVALID_OPTION) printf("invalid option on input.\n");
      if (ierror == TR_OUT_OF_MEMORY) printf("memory error.\n");
	  /* Release internal MKL memory that might be used for computations         */
	  /* NOTE: It is important to call the routine below to avoid memory leaks   */
	  /* unless you disable MKL Memory Manager                                   */
	  MKL_FreeBuffers();
		/* and exit */
		return 1;
	}
	/* set initial rci cycle variables */
	RCI_Request = 0;
	successful = 0;
	/* rci cycle */
    while (successful == 0)
	{
		/* call tr solver
			handle		in/out:	tr solver handle
			fvec		in:     vector
			fjac		in:     jacobi matrix
			RCI_request in/out:	return number which denote next step for performing */
		if (dtrnlsp_solve (&handle, fvec, fjac, &RCI_Request) != TR_SUCCESS)
		{
			/* if function does not complete successful then print error message */
			printf ("| error in dtrnlsp_solve\n");
			/* Release internal MKL memory that might be used for computations         */
			/* NOTE: It is important to call the routine below to avoid memory leaks   */
			/* unless you disable MKL Memory Manager                                   */
			MKL_FreeBuffers();
			/* and exit */
			return 1;
		}
		/* according with rci_request value we do next step */
		if (RCI_Request == -1 ||
			RCI_Request == -2 ||
			RCI_Request == -3 ||
			RCI_Request == -4 ||
			RCI_Request == -5 ||
			RCI_Request == -6)
			/* exit rci cycle */
			successful = 1;
	    if (RCI_Request == 1)
		{
			/* recalculate function value
				m		in:     dimension of function value
				n		in:     number of function variables
				x		in:     solution vector
				fvec    out:    function value f(x) */
            extendet_powell (&m, &n, x, fvec);
		}
		if (RCI_Request == 2)
		{
			/* compute jacobi matrix
				extendet_powell	in:     external objective function
				n               in:     number of function variables
				m               in:     dimension of function value
				fjac            out:    jacobi matrix
				x               in:     solution vector
				jac_eps         in:     jacobi calculation precision */
			if (djacobi (extendet_powell, &n, &m, fjac, x, eps) != TR_SUCCESS)
			{
				/* if function does not complete successful then print error message */
				printf ("| error in djacobi\n");
				/* Release internal MKL memory that might be used for computations         */
				/* NOTE: It is important to call the routine below to avoid memory leaks   */
				/* unless you disable MKL Memory Manager                                   */
				MKL_FreeBuffers();
				/* and exit */
				return 1;
			}
		}
	}
	/* get solution statuses
		handle            in:	TR solver handle
		iter              out:	number of iterations
		st_cr             out:	number of stop criterion
		r1                out:	initial residuals
		r2                out:	final residuals */
	if (dtrnlsp_get (&handle, &iter, &st_cr, &r1, &r2) != TR_SUCCESS)
	{
		/* if function does not complete successful then print error message */
		printf ("| error in dtrnlsp_get\n");
	  /* Release internal MKL memory that might be used for computations         */
	  /* NOTE: It is important to call the routine below to avoid memory leaks   */
	  /* unless you disable MKL Memory Manager                                   */
	  MKL_FreeBuffers();
		/* and exit */
		return 1;
	}
	/* free handle memory */
	if (dtrnlsp_delete (&handle) != TR_SUCCESS)
	{
		/* if function does not complete successful then print error message */
		printf ("| error in dtrnlsp_delete\n");
	  /* Release internal MKL memory that might be used for computations         */
	  /* NOTE: It is important to call the routine below to avoid memory leaks   */
	  /* unless you disable MKL Memory Manager                                   */
	  MKL_FreeBuffers();
		/* and exit */
		return 1;
	}
	/* free allocated memory */
	free (x);
	free (fvec);
	free (fjac);
	/* Release internal MKL memory that might be used for computations         */
	/* NOTE: It is important to call the routine below to avoid memory leaks   */
	/* unless you disable MKL Memory Manager                                   */
	MKL_FreeBuffers();
	/* if final residual less then required precision then print pass */
  if (r2 < 0.00001)
  {
    printf ("|         dtrnlsp powell............PASS\n");
    return 0;
  }
  /* else print failed */
  else
  {
    printf ("|         dtrnlsp powell............FAILED\n");
    return 1;
  }
}

/* nonlinear system equations without constraints */
/* routine for extendet powell function calculation
   m     in:     dimension of function value
   n     in:     number of function variables
   x     in:     vector for function calculating
   f     out:    function value f(x) */
void extendet_powell (MKL_INT *m, MKL_INT *n, double *x, double *f)
{
	MKL_INT i;

	for (i = 0; i < (*n)/4; i++)
	{
		f [4*i] = x [4*i] + 10.0*x [4*i + 1];
		f [4*i + 1] = 2.2360679774998*(x [4*i + 2] - x [4*i + 3]);
		f [4*i + 2] = (x [4*i + 1] - 2.0*x [4*i + 2])*(x [4*i + 1] - 2.0*x [4*i + 2]);
		f [4*i + 3] = 3.1622776601684*(x [4*i] - x [4*i + 3])*(x [4*i] - x [4*i + 3]);
	}
	return;
}
[/cpp]

Gennady_F_Intel · ‎02-10-2011

Mrcio,

1) because of you are linking with iLP64 libraries, then please compile the code with -DMKL_ILP64options. Refer to the more details in the MKL Users's Guide.

or

2) link the example with lp64 lib's:

MKL_FLAGS = -L$(MKL_LIB_DIR) $(MKL_LIB_DIR)/libmkl_solver_lp64_sequential.a -Wl,--start-group -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -Wl,--end-group -lpthread

--Gennady

View solution in original post

Gennady_F_Intel · ‎02-10-2011