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

Problem using FGMRES

zhang__boyuan
新分销商 I
3,338 次查看

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 项奖励
1 解答
Khang_N_Intel
员工
3,008 次查看
0 项奖励
8 回复数
mecej4
名誉分销商 III
3,298 次查看

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 项奖励
zhang__boyuan
新分销商 I
3,294 次查看

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 项奖励
MRajesh_intel
主持人
3,261 次查看

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 项奖励
zhang__boyuan
新分销商 I
3,251 次查看

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 项奖励
Khang_N_Intel
员工
3,244 次查看

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 项奖励
Khang_N_Intel
员工
3,238 次查看

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 项奖励
Khang_N_Intel
员工
3,144 次查看

Hi Boyuan,


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


Best regards,

Khang


0 项奖励
Khang_N_Intel
员工
3,009 次查看

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 项奖励
回复