### Solvers not working

Posted:

**Wed Jul 03, 2013 6:23 am**Hi,

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 https://code.google.com/p/cusp-library/ library here to simplify handling of arrays and matrices, but I think you should be able to understand what's going on):

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)

https://www.dropbox.com/s/84qf138lq8p240w/AAETT.mtx

https://www.dropbox.com/s/gbjvutqx21q3gw2/baffETT.mtx

Iteration 2: (does not solve)

https://www.dropbox.com/s/53638a6wtacl4h5/AATVA.mtx

https://www.dropbox.com/s/p1ixl70eo0m5xat/baffTVA.mtx

Does anyone have any idea what's going on?

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 https://code.google.com/p/cusp-library/ 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() {

culaInitialize();

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());

else

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

culaShutdown();

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)

https://www.dropbox.com/s/84qf138lq8p240w/AAETT.mtx

https://www.dropbox.com/s/gbjvutqx21q3gw2/baffETT.mtx

Iteration 2: (does not solve)

https://www.dropbox.com/s/53638a6wtacl4h5/AATVA.mtx

https://www.dropbox.com/s/p1ixl70eo0m5xat/baffTVA.mtx

Does anyone have any idea what's going on?