Solvers not working

General CULA Dense (LAPACK & BLAS) support and troubleshooting. Use this forum if you are having a general problem or have encountered a bug.

Solvers not working

Postby deveight » Wed Jul 03, 2013 6:23 am


The solver is meant to be part of an iterative algorithm where the (square) A and b in the Ax = b equations are slightly changed (only values, structure is exactly the same) between the iterations. The first iteration is solved without any problems, but during the second the solver will not provide an even remotely good answer (and no errors are reported).

I've reduced it to a minimal example, as the following (I'm also using the CUSP library here to simplify handling of arrays and matrices, but I think you should be able to understand what's going on):
Code: Select all
#include "cula.hpp"
#include <cusp/multiply.h>
#include <cusp/transpose.h>
#include <cusp/io/matrix_market.h>
#include <cusp/print.h>

using namespace std;

typedef double valueType;
typedef cusp::array1d<valueType, cusp::host_memory> Array1H;
typedef cusp::array1d<int, cusp::host_memory> iArray1H;
typedef cusp::array2d<valueType, cusp::host_memory> Array2H;

int main() {

   Array2H A2_h;
   Array1H A_h;
   Array1H originalb_h;

   cusp::io::read_matrix_market_file(A2_h, "AATVA.mtx");
   const int n = A2_h.num_cols;
   A_h = A2_h.values;

   Array1H b_h(n);
   Array1H temp_h(n);
   cusp::io::read_matrix_market_file(originalb_h, "baffTVA.mtx");
   cusp::blas::copy(originalb_h, b_h);

   iArray1H ipiv_h(n);
   culaStatus s = culaGesv(n, 1, &*A_h.begin(), n, &*ipiv_h.begin(), &*b_h.begin(), n); // b = A\b

   if(s != culaNoError) {
      if(s == culaArgumentError)
           printf("Argument %d has an illegal value\n", culaGetErrorInfo());
           printf("%s\n", culaGetStatusString(s));

   cusp::multiply(A2_h, b_h, temp_h); // temp = A*b
   cusp::blas::axpy(originalb_h, temp_h, -1.0f); // temp = temp - originalb
   cusp::print(temp_h); // prints the error, that is b - A*x

   return 0;

OS: latest Redhat
CUDA: 5.0
GPU: Telsa K20c

And here are the two A:s and b:s I'm using:
Iteration 1: (solves fine)

Iteration 2: (does not solve)

Does anyone have any idea what's going on?
Posts: 1
Joined: Fri Jun 28, 2013 9:32 am

Re: Solvers not working

Postby john » Wed Jul 03, 2013 10:02 am

These matrices have definitely been made dense, correct? And the symmetric portions have been populated? GESV assumes dense and fully filled, and the sparse MTX files are neither.

Keep in mind that GESV modifies the matrix and leaves the LU factorized matrix in its place, so if you call it twice with the same matrix it will not solve the second one properly. You need to use GETRS for subsequent solves against the same matrix.
Posts: 587
Joined: Thu Jul 23, 2009 2:31 pm

Return to CULA Dense Support

Who is online

Users browsing this forum: No registered users and 4 guests