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

Problem using FGMRES

zhang__boyuan
New Contributor I
1,676 Views

Hi,

I'm trying to use the RCI FGMRES function to solver a large sparse matrix linear equation. Essentially I used the example code fgmres_no_precon.c by modifing some of the input data. Specifically, I modified ia, ja, a and the rhs required by the algorithm, as well as the size of the matrix, n.

Here, I have two questions:

1). About the array tmp, I'm not sure what should be the size of array tmp. In the webpage

https://software.intel.com/content/www/us/en/develop/documentation/onemkl-developer-reference-c/top/sparse-solver-routines/iterative-sparse-solvers-based-on-reverse-communication-interface-rci-iss/rci-iss-routines/dfgmres.html

it mentioned that tmp should be an array of size ((2*ipar(15)+1)*n+ipar(15)*ipar(15)+9)/2 + 1).

However, in another webpage

https://software.intel.com/content/www/us/en/develop/documentation/onemkl-developer-reference-c/top/sparse-solver-routines/iterative-sparse-solvers-based-on-reverse-communication-interface-rci-iss/fgmres-interface-description.html

it mentioned here tmp should have size ((2*ipar(14)+1)*n+ipar(14)*ipar(14)+9)/2 + 1).

Which one of them is correct?

2). I cannot get the correct solution with my input data using FGMRES. I used MATLAB GMRES solver to solve the same matrix equation and got the solution I expected, where most elements in the solution vector x are of order 10^(-25), but the FGMRES solution are to the order 10^(-5) or so. 

I'm attaching a, ia, ja and rhs data for reference. The size of my problem is n = 28028. I also attached the result I got from MATLAB in the x.csv file, and the FGMRES part in my code for your convenience. 

Any advice on why the FGMRES solver is not working as expected will be greatly appreciated.

 

Best regards,

Boyuan

0 Kudos
1 Solution
Khang_N_Intel
Employee
1,346 Views
0 Kudos
8 Replies
mecej4
Honored Contributor III
1,636 Views

The expression for the size of tmp given in your first link appears to have been copied from the Fortran documentation for dfgmres; it has unbalanced parentheses. Unfortunately, such errors are not infrequent in the MKL documentation.

The second link that you gave contains the correct C expression that matches the expression in the MKL-Fortran documentation for dfgmres. However, when you copied the expression into your post you probably replaced ipar[14] by ipar(14) and left out a '('. The default lower bound for arrays is 1 in Fortran and 0 in C; subscripts are placed between parentheses in Fortran and between brackets in C.

Not much can be done with the code or the data in your attachments. The code is incomplete and cannot be compiled. The data files cannot be read without your telling us the format, structure and content of the files. As of now, you have a C code that is murky, and Matlab results for an unstated problem. The simplest guess for an explanation is that they solve different problems.

I suggest that you present complete evidence (code, data, instructions to compile and run) and commentary for a smaller value of n, say, 10 or 20.

0 Kudos
zhang__boyuan
New Contributor I
1,632 Views

Hi,

Thanks for the reply. Actually, what confused me about the size of array tmp is that in the first link it is ipar[14] in the expression, whereas in the second link it is ipar[15]. I think in general ipar[14] is not equal to ipar[15].

As for the unexpected result, I'm pretty sure that MATLAB and MKL FGMRES solved the same problem, since I just exported the coefficient matrix and the rhs from C code to MATLAB. I will think about how I can make it easier for others to look into my problem though.

 

Best,

Boyuan

0 Kudos
MRajesh_intel
Moderator
1,599 Views

Hi,

 

Can you please let us know what ipar value you were using in your code and did the results change by replacing iparm parameter from iparm[14] to iparm[15] or vice-versa.

 

Also, can you please share your MKL version, OS version?

 

Thanks

Rajesh.

 

0 Kudos
zhang__boyuan
New Contributor I
1,589 Views

Hi Rajesh,

 

I set ipar[4]=200 and ipar[14] = 150. I didn't assign a value to ipar[15]. When the iteration is done (satisfied the user-defined stopping criteria), ipar[15]=4232229.

 

The MKL I'm using is from the latest version of oneAPI, and my OS is Windows 10.

 

I think I solved the problem by introducing user-defined stopping criteria in the example code fgmres_full_funct.c. Another thing I did is to give the correct order (~10^-24) to the initial guess which helps a lot in achieving convergence.

 

In the other example, fgmres_no_precon.c, the user-defined stopping criteria is not used. At the begining I followed this example where I thought the relative error parameter dpar[0] should be able to have convergence control and give me the expected solution, but it seems not the case. 

 

Thanks again for following up this thread.

Best,

Boyuan

 

 

0 Kudos
Khang_N_Intel
Employee
1,582 Views

Hi Boyuan,


Mecej4 is right! Fortran array starts with index 1 while C array starts with 0.

So, use ipar[14] if you are using C. This is an error in the document. We will fix it.


Best,

Khang


0 Kudos
Khang_N_Intel
Employee
1,576 Views

Hi Boyuan,


The owner of the document is done with fixing the info. It might take a while before you can see the changes onlie.


Best,

Khang


0 Kudos
Khang_N_Intel
Employee
1,482 Views

Hi Boyuan,


The document will be updated in the next release of oneMKL, 2021.4


Best regards,

Khang


0 Kudos
Khang_N_Intel
Employee
1,347 Views

Hi Boyuan,


The document has been updated:

https://software.intel.com/content/www/us/en/develop/documentation/onemkl-developer-reference-c/top/sparse-solver-routines/iterative-sp-solvers-reverse-comm-iface-rci-iss/rci-iss-routines/dfgmres.html


Since the solution has been provided, there will be no more communication on this thread.


Best regards,

Khang



0 Kudos
Reply