Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.
6539 Discussions

Re: tridiagonal solver in oneMKL for dpc++ (linear solver routine getrs())

Oluwatosin
New Contributor I
646 Views

Hello,

I have been trying to use the getrs() routine, but I can't seem to get pass, size = 50 when executing the code, I need to be able to measure for size close to 100,000,000. I have to call first the getrf() routine according to the documentation, and I think this is where I miss it. I would appreciate any suggestions to improve the getrf() section of program; thank you. 

I am also following the matrix_multiplication with gemm code sample to generate my code;

The basic idea of what I am working on is:  
1. Create a matrix, A, of size n * n;
2. Initialize the matrix A, to zero;
3. Populate the sub, super and main diagonals with -1.0, -1,0 and 2.0  or random numbers;
4. Create a vector, for the right hand side and initialize it randomly.
4. Call the getrf() routine to factorize the system (device code);
5. Call the getrs() routine to solve the system (device code);
6. Compare the results of getrs() with results of host tridiagonal solver algorithm.

Here's what I have for the first four steps (getrf() )


int main(int argc, char* argv[])

{

size_t size = atoi(argv[1]);

 

sycl::queue device_queue{sycl::default_selector{}};

sycl::context context  = device_queue.get_context();

sycl::device device = device_queue.get_device();

 

int_64 i , j;

int64_t m = size, n = size, nrhs = size, lda = size;

 

// Allocate memory
double **A = sycl::malloc_shared<double*>(m * n, device_queue);

for(i = 0; i < lda; i++){

A[i] = sycl::malloc_shared<double>(lda, device_queue);}

 

// Initialize array to zero;
for(i = 0; i < m; i++){

for(j = 0; j < n; j++){
A[i][j] = 0.0;}}

// Populate sub,super, main diagonals of A

for(i = 0; i < lda; i++){

A[i][i] = 2.0; // Main diagonal

if(i < lda){ 

A[i][i+1] = -1.0;} // Super diagonal

if(i > 0){

A[i][i-1] = -1.0;} // Sub diagonal

}

// Other arrays

int64_t *ipiv = sycl::malloc_shared<std::int64_t>(ldb, device_queue);

auto b = sycl::malloc_shared<double>(ldb, device_queue);

 

//Initialize result array, b

for(i = 0; i < ldb; i++){

b[i] = 1.0;}

 

int64_t info = 0; 

 

try

{

int64_t scratchpad_size_getrf = oneapi::mkl::lapack::getrf_scratchpad_size<double>(device_queue, m, n, lda);

double* scratchpad_getrf = sycl::malloc_shared<double>(scratchpad_size_getrf, device, context);

auto event1 = oneapi::mkl::lapack::getrf(device_queue, m, n, *A, lda, ipiv, scratchpad_getrf, scratchpad_getrf, scratchpad_size_getrf);

event1.wait_and_throw();

cl::sycl::fre(scartchpad_getrf, context);

} catch(oneapi::mkl::lapack::exception const& e1){

std::cout << ............

if(e1.info() > 0){

info = e1.info();

}

return info;}

}


I would appreciate any suggestions to improve this section of the program.Thank you!

Labels (1)
0 Kudos
1 Solution
MRajesh_intel
Moderator
638 Views

Hi,


>>I would appreciate any suggestions to improve the getrf() section of program

Please refer to the sample(getri.cpp) present in $MKL_ROOT/examples/dpcpp/lapack/source folder.


Regards

Rajesh.


View solution in original post

7 Replies
MRajesh_intel
Moderator
639 Views

Hi,


>>I would appreciate any suggestions to improve the getrf() section of program

Please refer to the sample(getri.cpp) present in $MKL_ROOT/examples/dpcpp/lapack/source folder.


Regards

Rajesh.


Oluwatosin
New Contributor I
623 Views
Hi, thanks. I am working on the DevCloud. Can I get an external link, like a git repository, to download this example? Thank you.
MRajesh_intel
Moderator
612 Views

Hi,

 

You can view the MKL samples in the dir:

/opt/intel/oneapi/mkl/latest/examples

 

You can unzip the examples_dpcpp.tgz and check the getri.cpp example for reference.

 

Regards

Rajesh.

 

 

 

Oluwatosin
New Contributor I
604 Views
MRajesh_intel
Moderator
586 Views

Hi,

 

Can you please let us know if your issue has been resolved?

 

Regards

Rajesh.

 

Oluwatosin
New Contributor I
578 Views
Hello,
Yes, thank you!
MRajesh_intel
Moderator
574 Views

Hi,


Thanks for the confirmation!


As this issue has been resolved, we will no longer respond to this thread. If you require any additional assistance from Intel, please start a new thread. Any further interaction in this thread will be considered community only.


Have a Good day.


Regards

Rajesh



Reply