void MKL_fgmres(d_matrices& d) {

	int N = 2 * d.size_E_reduced;
	MKL_INT* ipar = MKLINT_Vector(128);
	double* dpar = FLOAT_Vector(128);
	double* tmp = FLOAT_Vector(N * (2 * 1000 + 1) + 1000 * (1000 + 9) / 2 + 1);

	//double* expected_solution = FLOAT_Vector(N);
	//expected_solution[0] = 1.0;
	//expected_solution[1] = -1.0;  // What is this step for?

	//double* rhs;
	//rhs = d.J_double;
	//double rhs[78112] = { 0.0 };

	//double rhs[N] = { 0.0 }, b[N] = { 0.0 };
	//double* b = FLOAT_Vector(N);
	double* computed_solution=FLOAT_Vector(N);
	//double* residual=FLOAT_Vector(N);

	MKL_INT itercount;// expected_itercount = 4;
	MKL_INT RCI_request, i, ivar;
	double dvar;

	// Descriptor of main sparse matrix properties
	struct matrix_descr descrA;


	// Structure with sparse matrix stored in CSR format
	sparse_matrix_t       csrA;
	sparse_operation_t    transA = SPARSE_OPERATION_NON_TRANSPOSE;
	descrA.type = SPARSE_MATRIX_TYPE_GENERAL;
	descrA.mode = SPARSE_FILL_MODE_FULL; /// ?????
	descrA.diag = SPARSE_DIAG_NON_UNIT;
	mkl_sparse_d_create_csr(&csrA, SPARSE_INDEX_BASE_ONE, N, N, d.i_MA_double, d.i_MA_double + 1, d.j_MA_double, d.v_MA_double);

	printf("--------------------------------------------------------\n");
	printf("The FULLY ADVANCED example of usage of RCI FGMRES solver\n");
	printf("   to solve a non-symmetric indefinite non-degenerate\n");
	printf("          algebraic system of linear equations\n");
	printf("--------------------------------------------------------\n\n");
	/*---------------------------------------------------------------------------
	* Initialize variables and the right hand side through matrix-vector product
	*---------------------------------------------------------------------------*/
	ivar = N;
	//mkl_sparse_d_mv(transA, 1.0, csrA, descrA, expected_solution, 0.0, rhs);
	/*---------------------------------------------------------------------------
    * Initialize the initial guess
    *---------------------------------------------------------------------------*/
	for (i = 0; i < N; i++)
	{
		computed_solution[i] = 0.5;
	}
	/*---------------------------------------------------------------------------
	* Initialize the solver
	*---------------------------------------------------------------------------*/
	dfgmres_init(&ivar, computed_solution, d.J_double, &RCI_request, ipar, dpar, tmp);
	if (RCI_request != 0)
		goto FAILED;
	/*---------------------------------------------------------------------------
	* Set the desired parameters:
	* LOGICAL parameters:
	* do residual stopping test
	* do not request for the user defined stopping test
	* do the check of the norm of the next generated vector automatically
	* DOUBLE PRECISION parameters
	* set the relative tolerance to 1.0D-3 instead of default value 1.0D-6
	*---------------------------------------------------------------------------*/
	ipar[4] = 1000; // Maximum iteration number.
	ipar[8] = 1;
	ipar[9] = 0;
	ipar[11] = 1;
	dpar[0] = 1.0E-3;
	/*---------------------------------------------------------------------------
	* Check the correctness and consistency of the newly set parameters
	*---------------------------------------------------------------------------*/
	dfgmres_check(&ivar, computed_solution, d.J_double, &RCI_request, ipar, dpar, tmp);
	if (RCI_request != 0 && RCI_request != -1001)
		goto FAILED;
	/*---------------------------------------------------------------------------
 * Print the info about the RCI FGMRES method
 *---------------------------------------------------------------------------*/
	printf("Some info about the current run of RCI FGMRES method:\n\n");
	if (ipar[7])
	{
		printf("As ipar[7]=" IFORMAT ", the automatic test for the maximal number of ", ipar[7]);
		printf("iterations will be\nperformed\n");
	}
	else
	{
		printf("As ipar[7]=" IFORMAT ", the automatic test for the maximal number of ", ipar[7]);
		printf("iterations will be\nskipped\n");
	}
	printf("+++\n");
	if (ipar[8])
	{
		printf("As ipar[8]=" IFORMAT ", the automatic residual test will be performed\n", ipar[8]);
	}
	else
	{
		printf("As ipar[8]=" IFORMAT ", the automatic residual test will be skipped\n", ipar[8]);
	}
	printf("+++\n");
	if (ipar[9])
	{
		printf("As ipar[9]=" IFORMAT ", the user-defined stopping test will be ", ipar[9]);
		printf("requested via\nRCI_request=2\n");
	}
	else
	{
		printf("As ipar[9]=" IFORMAT ", the user-defined stopping test will not be ", ipar[9]);
		printf("requested, thus,\nRCI_request will not take the value 2\n");
	}
	printf("+++\n");
	if (ipar[10])
	{
		printf("As ipar[10]=" IFORMAT ", the Preconditioned FGMRES iterations will be ", ipar[10]);
		printf("performed, thus,\nthe preconditioner action will be requested via");
		printf("RCI_request=3\n");
	}
	else
	{
		printf("As ipar[10]=" IFORMAT ", the Preconditioned FGMRES iterations will not ", ipar[10]);
		printf("be performed,\nthus, RCI_request will not take the value 3\n");
	}
	printf("+++\n");
	if (ipar[11])
	{
		printf("As ipar[11]=" IFORMAT ", the automatic test for the norm of the next ", ipar[11]);
		printf("generated vector is\nnot equal to zero up to rounding and ");
		printf("computational errors will be performed,\nthus, RCI_request will not take the value 4\n");
	}
	else
	{
		printf("As ipar[11]=" IFORMAT ", the automatic test for the norm of the next ", ipar[11]);
		printf("generated vector is\nnot equal to zero up to rounding and ");
		printf("computational errors will be skipped,\nthus, the user-defined test ");
		printf("will be requested via RCI_request=4\n");
	}
	printf("+++\n\n");
	/*---------------------------------------------------------------------------
	* Compute the solution by RCI (P)FGMRES solver without preconditioning
	* Reverse Communication starts here
	*---------------------------------------------------------------------------*/
ONE:dfgmres(&ivar, computed_solution, d.J_double, &RCI_request, ipar, dpar, tmp);
	/*---------------------------------------------------------------------------
	* If RCI_request=0, then the solution was found with the required precision
	*---------------------------------------------------------------------------*/
	if (RCI_request == 0)
		goto COMPLETE;
	/*---------------------------------------------------------------------------
	* If RCI_request=1, then compute the vector A*tmp[ipar[21]-1]
	* and put the result in vector tmp[ipar[22]-1]
	*---------------------------------------------------------------------------
	* NOTE that ipar[21] and ipar[22] contain FORTRAN style addresses,
	* therefore, in C code it is required to subtract 1 from them to get C style
	* addresses
	*---------------------------------------------------------------------------*/
	if (RCI_request == 1)
	{
		mkl_sparse_d_mv(transA, 1.0, csrA, descrA, &tmp[ipar[21] - 1], 0.0, &tmp[ipar[22] - 1]);
		goto ONE;
	}
	/*---------------------------------------------------------------------------
	* If RCI_request=anything else, then dfgmres subroutine failed
	* to compute the solution vector: computed_solution[N]
	*---------------------------------------------------------------------------*/
	else
	{
		goto FAILED;
	}
	/*---------------------------------------------------------------------------
	* Reverse Communication ends here
	* Get the current iteration number and the FGMRES solution (DO NOT FORGET to
	* call dfgmres_get routine as computed_solution is still containing
	* the initial guess!)
	*---------------------------------------------------------------------------*/
COMPLETE:dfgmres_get(&ivar, computed_solution, d.J_double, &RCI_request, ipar, dpar, tmp, &itercount);
	/*---------------------------------------------------------------------------
	* Print solution vector: computed_solution[N] and the number of iterations: itercount
	*--------------------------------------------------------------------------- */
	printf(" The system has been solved \n");
	printf("\n The following solution has been obtained: \n");

	d.A = MKLCMP_Vector(2 * d.size_E_reduced);

	for (i = 0; i < d.size_E_reduced; i++)
	{
		d.A[i].real = computed_solution[i];
		d.A[i].imag = computed_solution[i + d.size_E_reduced];
		//printf("computed_solution[" IFORMAT "]=", i);
		//printf("%e\n", computed_solution[i]);
	}
	/*printf("\n The expected solution is: \n");
	for (i = 0; i < N; i++)
	{
		printf("expected_solution[" IFORMAT "]=", i);
		printf("%e\n", expected_solution[i]);
		expected_solution[i] -= computed_solution[i];
	}*/
	printf("\n Number of iterations: " IFORMAT "\n", itercount);
	i = 1;
	//dvar = dnrm2(&ivar, expected_solution, &i);

	/*-------------------------------------------------------------------------*/
	/* Release internal  Intel oneMKL memory that might be used for computations         */
	/* NOTE: It is important to call the routine below to avoid memory leaks   */
	/* unless you disable  Intel oneMKL Memory Manager                                   */
	/*-------------------------------------------------------------------------*/
	MKL_Free_Buffers();

	/*if (itercount == expected_itercount && dvar < 1.0e-14)
	{
		printf("\nThis example has successfully PASSED through all steps of computation!\n");
		return 0;
	}
	else
	{
		printf("\nThis example may have FAILED as either the number of iterations differs\n");
		printf("from the expected number of iterations " IFORMAT ", ", expected_itercount);
		printf("or the computed solution\n");
		printf("differs much from the expected solution (Euclidean norm is %e), or both.\n", dvar);
		return 1;
	}*/

	/*-------------------------------------------------------------------------*/
	/* Release internal  Intel oneMKL memory that might be used for computations         */
	/* NOTE: It is important to call the routine below to avoid memory leaks   */
	/* unless you disable  Intel oneMKL Memory Manager                                   */
	/*-------------------------------------------------------------------------*/
FAILED:printf("\nThis example FAILED as the solver has returned the ERROR code " IFORMAT, RCI_request);
	mkl_sparse_destroy(csrA);
	MKL_Free_Buffers();

	//int kkk = 0;
}