/* Example Program, with memory allocated, 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>

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

  /* Arrays */
  char   *title=0;
  double *a_row=0, *a_column=0, *tau=0;
  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");


#define A_COLUMN(I,J) a_column[(J-1)*pda_column + I - 1]
#define A_ROW(I,J) a_row[(I-1)*pda_row + J - 1]

  INIT_FAIL(fail);
  
  m = 6;
  n = 4;;
  pda_column = m;
  pda_row = n;
  tau_len = MIN(m, n);

  /* Allocate memory */
  if ( !(title = NAG_ALLOC(31, char)) ||
       !(a_row = NAG_ALLOC(m * n, double)) ||
       !(a_column = NAG_ALLOC(m * n, double)) ||
       !(tau = NAG_ALLOC(tau_len, double)) )
    {
      Vprintf("Allocation failure\n");
      exit_status = -1;
      goto END;
    }

#ifdef NAG_ROW_MAJOR
  Vprintf("Using row major storage, allocated memory\n");
  /* Read A from data above */
  for (i = 1; i <= m; ++i)
    {
      for (j = 1; 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 */
  f08aec(Nag_RowMajor, m, n, a_row, pda_row, 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, pda_row, 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, pda_row, 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, allocated memory\n");
  /* Read A from data above */
  for (i = 1; i <= m; ++i)
    {
      for (j = 1; j<= n; j++)
        {
          sscanf(matrix_data_ptr, "%lf", &A_COLUMN(i,j));
          matrix_data_ptr = strtok(0, " \t\n");
        }
    }

  f08aec(Nag_ColMajor, m, n, a_column, pda_column, 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, pda_column, 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, pda_column, title, 0, &fail);
  if (fail.code != NE_NOERROR)
    {
      Vprintf("Error from x04cac.\n%s\n", fail.message);
      exit_status = 1;
      goto END;
    }
#endif
 END:
  if (title) NAG_FREE(title);
  if (a_row) NAG_FREE(a_row);
  if (a_column) NAG_FREE(a_column);
  if (tau) NAG_FREE(tau);

  return exit_status;
}