Developer Reference for Intel® oneAPI Math Kernel Library for C

ID 766684
Date 3/22/2024
Public

A newer version of this document is available. Customers should click here to go to the newest version.

Document Table of Contents

dfgmres

Makes the FGMRES iterations.

Syntax

void dfgmres (const MKL_INT *n, double *x, double *b, MKL_INT *RCI_request, MKL_INT *ipar, double *dpar, double *tmp );

Include Files

  • mkl.h

Description

The routine dfgmres performs the FGMRES iterations [Saad03], using the value that was in the array x before the first call as an initial approximation of the solution vector. To update the current approximation to the solution, the dfgmres_get routine must be called. The RCI FGMRES iterations can be continued after the call to the dfgmres_get routine only if the value of the parameter ipar[12] is not equal to 0 (default value). Note that the updated solution overwrites the right-hand side in the vector b if the parameter ipar[12] is positive, and the restarted version of the FGMRES method can not be run. If you want to keep the right-hand side, you must be save it in a different memory location before the first call to the dfgmres_get routine with a positive ipar[12].

The parameter RCI_request gives information about the task completion and requests results of certain operations that the solver requires.

The lengths of all the vectors must be defined in a previous call to the dfgmres_init routine.

Input Parameters

n

Sets the size of the problem.

x

Array of size n. Contains the initial approximation to the solution vector.

b

Array of size n. Contains the right-hand side vector.

tmp

Array of size [12]. Refer to the FGMRES Common Parameters.

Output Parameters

RCI_request

Informs about result of work of the routine.

ipar

Array of size 128. Refer to the FGMRES Common Parameters.

dpar

Array of size 128. Refer to the FGMRES Common Parameters.

tmp

Array of size ((2*ipar[14]+1)*n+ipar[14]*ipar[14]+9)/2 + 1. Refer to the FGMRES Common Parameters.

Return Values

RCI_request=0

Indicates that the task completed normally and the solution is found and stored in the vector x. This occurs only if the stopping tests are fully automatic. For the user defined stopping tests, see the description of the RCI_request= 2 or 4.

RCI_request=-1

Indicates that the routine was interrupted because the maximum number of iterations was reached, but the relative stopping criterion was not met. This situation occurs only if you request both tests.

RCI_request= -10

Indicates that the routine was interrupted because of an attempt to divide by zero. Usually this happens if the matrix is degenerate or almost degenerate. However, it may happen if the parameter dpar is altered, or if the method is not stopped when the solution is found.

RCI_request= -11

Indicates that the routine was interrupted because it entered an infinite cycle. Usually this happens because the values ipar[7], ipar[8], ipar[9] were altered outside of the routine, or the dfgmres_check routine was not called.

RCI_request= -12

Indicates that the routine was interrupted because errors were found in the method parameters. Usually this happens if the parameters ipar and dpar were altered by mistake outside the routine.

RCI_request= 1

Indicates that you must multiply the matrix by tmp[ipar[21] - 1:ipar[21] + n - 2], put the result in the tmp[ipar[22] - 1:ipar[22] + n - 2], and return control back to the routine dfgmres.

RCI_request= 2

Indicates that you must perform the stopping tests. If they fail, return control to the dfgmres routine. Otherwise, the FGMRES solution is found, and you can run the fgmres_get routine to update the computed solution in the vector x.

RCI_request= 3

Indicates that you must apply the inverse preconditioner to tmp[ipar[21] - 1:ipar[21] + n - 2], put the result in the tmp[ipar[22] - 1:ipar[22] + n - 2], and return control back to the routine dfgmres.

RCI_request= 4

Indicates that you must check the norm of the currently generated vector. If it is not zero within the computational/rounding errors, return control to the dfgmres routine. Otherwise, the FGMRES solution is found, and you can run the dfgmres_get routine to update the computed solution in the vector x.