Solving with multiple RHS

[badlishah]

Greetings,

I am trying to get X for the AX=B type of problem. I have a large banded sparse matrix A (200k x 200k) of single precision complex numbers and multiple RHS as the sparse matrix B. At present I am using the Pardiso direct solver to do this. What would be the best approach to get the solution X if I were to use Magma's sparse iterative solver instead?

AFAIK sparse RHS is not yet supported, would there be a workaround for this problem? I have read about the use of block-Krylov methods but I am not familiar with them.

Any advice is much appreciated. Thanks.

cheers,
Badlishah

[hartwig anzt]

Badlishah,

at the current state, I would suggest to transform the sparse RHS-matrix into a dense matrix, such that you have a block of vectors (dense matrix) RHS.
magma_zmconvert( B, &B2, Magma_CSR, Magma_DENSE, queue );

Then, you use an initial guess for the block solution (X) also as dense matrix.
magma_zvinit(&X, Magma_CPU, B.num-rows, B.num_cols, MAGMA_Z_ZERO, queue );

With this setup you can use any of the Krylov solvers: CG, BiCGSTAB, GMRES...

Let me know whether this helps!

Hartwig

[badlishah]

Hartwig,

Thanks for the tip. I will try this out and post the observations here.

cheers,
Badlishah

[badlishah]

Hartwig,

I have a little problem here. I may have missed something. I keep getting this error:
error: only 1 RHS supported for this solver class.

I am using GEMRES solver.

my code snippet:

Code: Select all

	loadAMatrix(&A_h, queue, argv[1]);
	loadData(&B_h, argv[2], argv[3], argv[4], A_h.num_rows, queue);

	magma_cmtransfer( A_h, &A_d, Magma_CPU, Magma_DEV, queue );
	magma_cmconvert( B_h, &B_D, Magma_CSR, Magma_DENSE, queue );
	magma_cmtransfer( B_D, &B_d, Magma_CPU, Magma_DEV, queue );
	magma_c_vinit( &X_d, Magma_DEV, B_D.num_rows, B_D.num_cols, zero, queue );
	info = magma_c_solver( A_d, B_d, &X_d, &zopts, queue );

Did I missed anything ?

cheers,
Badlishah

[hartwig anzt]

Badlishah,

indeed, not all solvers support block-RHS at this point. Please try with preconditioned Conjugate Gradient (called zbpcg).

Thanks, Hartwig

[badlishah]

Hartwig,

Thanks a lot. I don't get that error anymore. I am using magma_cbpcg and preconditioner set to none. However I am getting solver returned: not positive definite (SPD/HPD) (-203) error instead. I am checking the program whether data is read correctly or not.What would cause such error?


cheers,
Badlishah

[hartwig anzt]

Badlishah,

this error typically occurs if the linear system is not positive definite, which is needed for using CG. In particular, the system needs to be positive ( A=A' ) and all eigenvalues of the system need to be positive. Is this true for your system? What is the origin?

Thanks, Hartwig

[badlishah]

Hi Hartwig,

I am trying to solve a partial differential equation, to be exact the inverse of the wave equation of some seismic data. All this while I am using a direct solver, however direct solvers take a lot of memory and do not capitalize on the GPU capabilities. That is my main motivation to explore the use of iterative solvers instead .. if only I can get it right ;)

[hartwig anzt]

Maybe you can write out one sample matrix, such that I can try to figure out what is happening? If yes, just pass me the URL.

Thanks, Hartwig