/* Example Program, with memory declared, based on:
 *
 * nag_dorgqr (f08afc) Example Program.
 *
 * Copyright 2001 Numerical Algorithms Group.
 *
 * Mark 7, 2001.
 */

#include <stdio.h>
#include <string.h>
#include <nag.h>
#include <nag_stdlib.h>
#include <nagf08.h>
#include <nagx04.h>

#define MMAX 10
#define NMAX 8

int main(void)
{
  /* Scalars */
  Integer  i, j, m, n, pda, tau_len;
  Integer  exit_status=0;
  NagError fail;

  /* Arrays */
  char   title[30];
  double a_row[MMAX][NMAX], a_column[MMAX][NMAX], tau[NMAX];
  char matrx_data [] = {
    " -0.57  -1.28  -0.39   0.25 "
    " -1.93   1.08  -0.31  -2.14 "
    "  2.30   0.24   0.40  -0.35 "
    " -1.93   0.64  -0.66   0.08 "
    "  0.15   0.30   0.15  -2.13 "
    " -0.02   1.03  -1.43   0.50 "
  }, *matrix_data_ptr = matrx_data;

  /* Initialise strtok */
  matrix_data_ptr = strtok(matrix_data_ptr, " \t\n");


  INIT_FAIL(fail);
  
  m = 6;
  n = 4;;
  pda = NMAX;
  tau_len = MIN(m, n);

  /* Read A from data above */
#ifdef NAG_ROW_MAJOR
  for (i = 0; i < m; ++i)
    {
      for (j = 0; j< n; j++)
        {
          sscanf(matrix_data_ptr, "%lf", &a_row[i][j]);
          matrix_data_ptr = strtok(0, " \t\n");
        }
    }

  /* Compute the QR factorization of A */
  Vprintf("Using row major storage, declared memory\n");
  f08aec(Nag_RowMajor, m, n, &a_row[0][0], pda, tau, &fail);
  if (fail.code != NE_NOERROR)
    {
      Vprintf("Error from f08aec.\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
  /* Form the leading N columns of Q explicitly */
  f08afc(Nag_RowMajor, m, n, n, &a_row[0][0], pda, tau, &fail);
  if (fail.code != NE_NOERROR)
    {
      Vprintf("Error from f08afc.\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
  /* Print the leading N columns of Q only */
  Vsprintf(title, "The leading %2ld columns of Q\n", n);
  x04cac(Nag_RowMajor, Nag_GeneralMatrix, Nag_NonUnitDiag, m, n, 
	 &a_row[0][0], pda, title, 0, &fail);
  if (fail.code != NE_NOERROR)
    {
      Vprintf("Error from x04cac.\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
#else
  Vprintf("Using column major storage, declared memory\n");
  for (i = 0; i < m; ++i)
    for (j = 0; j< n; j++)
      {
	/* Note column data is transposed */
	sscanf(matrix_data_ptr, "%lf", &a_column[j][i]);
	matrix_data_ptr = strtok(0, " \t\n");
      }
  
  f08aec(Nag_ColMajor, m, n, &a_column[0][0], pda, tau, &fail);
  if (fail.code != NE_NOERROR)
    {
      Vprintf("Error from f08aec.\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
  /* Form the leading N columns of Q explicitly */
  f08afc(Nag_ColMajor, m, n, n, &a_column[0][0], pda, tau, &fail);
  if (fail.code != NE_NOERROR)
    {
      Vprintf("Error from f08afc.\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
  /* Print the leading N columns of Q only */
  Vsprintf(title, "The leading %2ld columns of Q\n", n);
  x04cac(Nag_ColMajor, Nag_GeneralMatrix, Nag_NonUnitDiag,  m, n,
	 &a_column[0][0], pda, title, 0, &fail);
  if (fail.code != NE_NOERROR)
    {
      Vprintf("Error from x04cac.\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
#endif
 END:
  return exit_status;
}