//      f08be Example Program Text
//      C# version, NAG Copyright 2008
using System;
using NagLibrary;
using System.Globalization;
namespace NagDotNetExamples
{
  public class F08BEE
  {
    static string datafile = "ExampleData/f08bee.d";
    static void Main(String[] args)
    {
      if (args.Length == 1)
      {
        datafile = args[0];
      }
      StartExample();
    }
    public static void StartExample()
    {
      try
      {
        DataReader sr = new DataReader(datafile);
  
        double tol=0.0;
        int i,  info,  j,  k,  m,  n,  nrhs;
        int ifail;
        Console.WriteLine("f08be Example Program Results");
        //      Skip heading in data file
        sr.Reset();
        sr.Reset();
        m = int.Parse(sr.Next());
        n = int.Parse(sr.Next());
        nrhs = int.Parse(sr.Next());
        double[,] a = new double[m, n];
        double[,] b = new double[m, nrhs];
        double[] tau = new double[Math.Min(m,n)];
        double[] work = new double[512];
        double[,] x = new double[n, nrhs];
        int[] jpvt = new int[n];
        if ((m >= 1 && n >= 1 && m >= n) && nrhs >= 1)
        {
          // 
          //         Read A and B from data file
          // 
          sr.Reset();
          for (i = 1 ; i <= m ; i++)
          {
            for (j = 1 ; j <= n ; j++)
            {
              a[i - 1 , j - 1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
            }
          }
          sr.Reset();
          for (i = 1 ; i <= m ; i++)
          {
            for (j = 1 ; j <= nrhs ; j++)
            {
              b[i - 1 , j - 1] = double.Parse(sr.Next(), CultureInfo.InvariantCulture);
            }
          }
          // 
          //         Initialize jpvt to be zero so that all columns are free
          // 
          F06.f06db(n, 0, jpvt, 1, out ifail);
          // 
          //         Compute the QR factorization of a
          // 
          F08.f08be(m, n, a, jpvt, tau, out info);

          if (info < 0)
          {
              return;
          }
          // 
          //         Choose tol to reflect the relative accuracy of the input data
          // 
          tol = 0.010e0;
          // 
          //         Determine which columns of R to use
          // 
          for (k = 1 ; k <= n ; k++)
          {
            if (Math.Abs(a[k - 1 , k - 1]) <= tol * Math.Abs(a[0 , 0]))
            {
              goto L40;
            }
          }
          // 
          //         Compute C = (Q**T)*B, storing the result in B
          // 
          L40: ;
          k = k - 1;
          // 
          F08.f08ag("Left", "Transpose", m, nrhs, n, a, tau, b, out info);
          // 
          //         Compute least-squares solution by backsubstitution in R*B = C
          //
          double [] bcolumn = new double [m];
          for (i = 1 ; i <= nrhs ; i++)
          {
            // 
            for (j = 0; j < k; j++)
            bcolumn[j] = b[j,i-1];
            F06.f06pj("Upper", "No transpose", "Non-Unit", k, a, bcolumn, 1, out ifail);
            // 
            //            Copy bcolumn to b's column including zero elements
            for (j = 0; j < m; j++)
            b[j,i-1] = bcolumn[j];
            // 
          }
          // 
          //         Unscramble the least-squares solution stored in B
          // 
          for (i = 1 ; i <= n ; i++)
          {
            for (j = 1 ; j <= nrhs ; j++)
            {
              x[jpvt[i - 1] - 1 , j - 1] = b[i - 1 , j - 1];
            }
          }
          // 
          //         Print least-squares solution
          // 
          Console.WriteLine("");
          // 
          X04.x04ca("General", "", n, nrhs, x, "Least-squares solution", out ifail);
          // 
        }
      }
      catch (Exception e)
      {
        Console.WriteLine(e.Message);
        Console.WriteLine("Exception Raised");
      }
    }
  }
}