Use RuntimeSettings to enable MAGMA iterative refinement in CLUSTER

  string message = "INFO: running scale iteration " + to_string(jxi488) + " of " + to_string(nxi) + ".\n";

  Logger *logger = new Logger(LOG_DEBG);

  logger->log(message);

  RuntimeSettings rs(gconf, logger);

  chrono::duration<double> elapsed;

  chrono::time_point<chrono::high_resolution_clock> interval_start, interval_end;

  int jer = 0;

#endif

  // we put the accuracygoal in, get the actual accuracy back out

  double actualaccuracy = cid->accuracygoal;

  invert_matrix(cid->am, ndit, jer, cid->maxrefiters, actualaccuracy, cid->refinemode, output_path, jxi488, mxndm, cid->proc_device);

  invert_matrix(cid->am, ndit, jer, cid->maxrefiters, actualaccuracy, cid->refinemode, output_path, jxi488, mxndm, cid->proc_device, rs);

  // in principle, we should check whether the returned actualaccuracy is indeed lower than the accuracygoal, and do something about it if not

  if (gconf->refine_flag) {

    if (cid->refinemode==2) {

      message = "DEBUG: iterative refinement enabled at run-time.\n";

      logger->log(message, LOG_DEBG);

      message = "INFO: calibration obtained accuracy " + to_string(actualaccuracy) + " (" + to_string(cid->accuracygoal) + " requested) in " + to_string(cid->maxrefiters) + " refinement iterations\n";

      logger->log(message);

      if (actualaccuracy > 1e-1) {

	printf("Accuracy worse than 0.1, stopping");

	exit(1);

      }

    }

  }

  cid->refinemode = 0;

#ifdef DEBUG_AM

  VirtualAsciiFile *outam2 = new VirtualAsciiFile();

using namespace std;

/*! \brief Perform in-place matrix inversion.

/**

 * \brief Perform in-place matrix inversion.

 *

 * \param mat: `complex double **` The matrix to be inverted (must be a square matrix).

 * \param size: `np_int` The size of the matrix (i.e. the number of its rows or columns).

 * \param jxi488: `int` Index of the current wavelength calculation.

 * \param max_size: `np_int` The maximum expected size (required by some call-backs, optional, defaults to 0).

 * \param target_device: `int` ID of target GPU, if available (defaults to 0).

 * \param rs: `const RuntimeSettings &` Runtime options instance (optional).

 */

void invert_matrix(dcomplex **mat, np_int size, int &ier, int &maxrefiters, double &accuracygoal, int refinemode, const string& output_path, int jxi488, np_int max_size=0, int target_device=0);

void invert_matrix(dcomplex **mat, np_int size, int &ier, int &maxrefiters, double &accuracygoal, int refinemode, const string& output_path, int jxi488, np_int max_size=0, int target_device=0, const RuntimeSettings& rs=RuntimeSettings());

#endif

   this program in the COPYING file. If not, see: <https://www.gnu.org/licenses/>.

 */

/*! \file magma_calls.h

/**

 * \file magma_calls.h

 *

 * \brief C++ interface to MAGMA calls.

 *

#ifndef INCLUDE_MAGMA_CALLS_H_

#define INCLUDE_MAGMA_CALLS_H_

/*! \brief Invert a complex matrix with double precision elements.

 *

 * call magma_zinvert1() to do the actual inversion

/**

 * \brief Invert a complex matrix with double precision elements.

 *

 * \param mat: Matrix of complex. The matrix to be inverted.

 * \param n: `np_int` The number of rows and columns of the [n x n] matrix.

 * \param jer: `int &` Reference to an integer return flag.

 * \param device_id: `int` ID of the device for matrix inversion offloading.

 * \param rs: `const RuntimeSettings &` Runtime settings instance.

 */

void magma_zinvert(dcomplex **mat, np_int n, int &jer, int device_id=0);

void magma_zinvert(

  dcomplex **mat, np_int n, int &jer, int device_id=0,

  const RuntimeSettings& rs=RuntimeSettings()

);

/*! \brief Auxiliary function to Invert a complex matrix with double precision elements.

 *

 * Use MAGMA to perform an in-place matrix inversion for a complex

 * matrix with double precision elements.

/**

 * \brief Perform Newton-Schulz iterative refinement of matrix inversion.

 *

 * \param inva: reference to the pointer to the first element of the matrix to be inverted on entry, inverted matrix on exit.

 * \param n: `np_int` The number of rows and columns of the [n x n] matrix.

 * \param jer: `int &` Reference to an integer return flag.

 * \param device_id: `int` ID of the device for matrix inversion offloading.

 * \param rs: `const RuntimeSettings &` Runtime settings instance. [IN]

 * \param a: `magmaDoubleComplex *` Pointer to the first element of the non-inverted matrix on host. [IN]

 * \param m: `const magma_int_t` Number of rows / columns in a. [IN]

 * \param d_a: `magmaDoubleComplex *` Pointer to the matrix on the GPU. [IN/OUT]

 * \param queue: `magma_queue_t` GPU communication queue. [IN]

 */

void magma_zinvert1(dcomplex * &inva, np_int n, int &jer, int device_id);

magma_int_t magma_newton(

  const RuntimeSettings& rs, magmaDoubleComplex* a, const magma_int_t m,

  magmaDoubleComplex* d_a, magma_queue_t queue

);

/*! \brief Invert a complex matrix with double precision elements, applying iterative refinement of the solution

/* \brief Invert a complex matrix with double precision elements, applying iterative refinement of the solution

 *

 * call magma_zinvert1() to perform the first matrix inversion, then magma_refine() to do the refinement (only if maxrefiters is >0)

 *

 * \param output_path: `const string &` Path where the output needs to be placed.

 * \param jxi488: `int` Index of the current wavelength calculation.

 */

void magma_zinvert_and_refine(dcomplex **mat, np_int n, int &jer, int &maxrefiters, double &accuracygoal, int refinemode, int device_id, const string& output_path, int jxi488);

//void magma_zinvert_and_refine(dcomplex **mat, np_int n, int &jer, int &maxrefiters, double &accuracygoal, int refinemode, int device_id, const string& output_path, int jxi488);

/*! \brief Apply iterative refinement of the solution of a matrix inversion.

/* \brief Apply iterative refinement of the solution of a matrix inversion.

 *

 * Iteratively compute and apply a correction to the inverse `inva` of the complex

 * matrix `aorig`, for a maximum number of `maxiters` times, or until achieving a

 * \param output_path: `const string &` Path where the output needs to be placed.

 * \param jxi488: `int` Index of the current wavelength calculation.

 */

void magma_refine(dcomplex *aorig, dcomplex *inva, np_int n, int &jer, int &maxrefiters, double &accuracygoal, int refinemode, int device_id, const string& output_path, int jxi488);

// void magma_refine(dcomplex *aorig, dcomplex *inva, np_int n, int &jer, int &maxrefiters, double &accuracygoal, int refinemode, int device_id, const std::string& output_path, int jxi488);

#endif