solve() [3/4]

bool Lapack::solve ( const DMatRR * A, const DMatRR * b, DMatRR * x )

static

Definition at line 240 of file lapack.cpp.

{
  int size = static_cast<int>(A->numRows());
  int bsize = static_cast<int>(b->numColumns());
  int info;
 
  /* make sure matrix is square */
  if (A->numRows() != A->numColumns())
    {
      ERROR("expected a square matrix");
      return false;
    }
 
  /* make sure dimensions of b make sense for Ax=b */
  if (b->numRows() != size)
    {
      ERROR("expected matrices to have same number of rows");
      return false;
      ;
    }
 
  if (size == 0)
    {
      x->resize(size, bsize);
      return true;
    }
 
  int* perm = new int[size];
  std::vector<double> copyA = make_lapack_array(*A);
  std::vector<double> copyb = make_lapack_array(*b);
  
  dgesv_(&size,
         &bsize,
         copyA.data(),
         &size,
         perm,
         copyb.data(),  // also the result
         &size,
         &info);
 
  delete [] perm; // Do we need this for anything??
 
  // Now set x
  x->resize(size, bsize);
  fill_from_lapack_array(copyb, *x);
 
  if (info > 0)
    {
      ERROR("according to dgesv, matrix is singular");
      return false;
    }
  else if (info < 0)
    {
      ERROR("argument passed to dgesv had an illegal value");
      return false;
    }
 
  return true;
}

References dgesv_(), ERROR, fill_from_lapack_array(), make_lapack_array(), DMat< ACoeffRing >::numColumns(), DMat< ACoeffRing >::numRows(), and x.