//set the integer pointer to zero for (i = 0; i < 64; i++) { pt[i] = 0; } mtype = -2; //number of right hand side vector nrhs = 1; //Print statistical information in file. 0 no, 1 yes msglvl = 0; //Which factorization to use mnum = 1; //Maximum number of numerical factorization maxfct = 1; for (i = 0; i < 64; i++) { iparm[i] = 0; } pardisoinit (pt, &mtype, iparm); iparm[7] = 15; // Max numbers of iterative refinement steps // iparm[10] = 1; // Use nonsymmetric permutation and scaling MPS // used to gether with iparm[12] iparm[12] = 1; // Maximum weighted matching algorithm is switched-off (default for symmetric). Try iparm[12] = 1 in case of inappropriate accuracy // iparm[0] = 1; // No solver default // iparm[1] = 2; //0 The minimum degree algorithm [Li99]. //2* The nested dissection algorithm from the METIS package [Karypis98]. //3 The parallel (OpenMP) version of the nested dissection algorithm. It can decrease the time of computations on multi-core computers, especially when Intel MKL PARDISO Phase 1 takes significant time. // Numbers of processors, value of OMP_NUM_THREADS // //iparm[2] = 8; iparm[3] = 0; // No iterative-direct algorithm // iparm[4] = 0; // No user fill-in reducing permutation // iparm[5] = 0; // Write solution into x // iparm[6] = 0; // output: number of perfomred iterative refinement factorization iparm[9] = 13; // Perturb the pivot elements with 1E-13 // //iparm[10] = 1; // scaling amtrix //iparm[12] = 1; // scaling amtrix iparm[11] = 0; // solving with transposed or conjugate transposed matrix. //If = 0, PARDISO solves a linear system Ax = b (default value). iparm[13] = 0; // Output: Number of perturbed pivots // //set output <0 to get the output iparm[17] = -1; // Output: Number of nonzeros in the factor LU // iparm[18] = 0; // Output: Mflops for LU factorization // iparm[19] = 0; // Output: Numbers of CG Iterations // iparm[20] = 0; // control the pivoting method for sparse symmetric indefinite // 0: 1x1 diagonal pivoting //1: 1x1 or 2x2 Bunck and Kaufman pivoting //iparm[23] = 1; //parallel factorization control. //0 (default value), then PARDISO uses the previous parallel factorization. //1, then PARDISO uses new two-level scheduling algorithm. This algorithm generally improves scalability in case of parallel factorization on many threads (more than eight). //iparm[24] = 1;parallel forward/backward solve control. // 0 (default value), then PARDISO uses a parallel algorithm for the solve step. // 1, then PARDISO uses sequential forward and backward solve. // This feature is available only for in-core version. iparm[26] = 1; //1. check matrix. 0. no check// iparm[27] = 0; //0. double 1. single precision iparm[34] = 1; //0. fortran array 1. C array iparm[59] = 0; // Depend on the environement variables . FIXME Need to get that variables.// //0: incore; 1: depend on the evironment variable, 2: OOC solution numRows = A.GetSizeM(); assert(numRows == (_INTEGER_t) (b.size())); _INTEGER_t phase,i; // Initialize error flag CString str1; if (x.size() != numRows){ x.resize(numRows); } // .. Reordering and Symbolic Factorization. This step also allocates all memory that is necessary for the factorization. phase = 11; status(INFO_INFO, "Reordering and Symbolic Factorization..."); PARDISO_64 (pt, &maxfct, &mnum, &mtype, &phase,&numRows, A.GetValues(), A.GetRowIndex(), A.GetColumns(), &idum, &nrhs,iparm, &msglvl, &ddum, &ddum, &error); if (error!=0) { switch(error) { case -1: status(INFO_INFO, "Input inconsistent"); //std::cerr << "Input inconsistent" << std::endl; break; case -2: status(INFO_INFO, "Not enough memory"); //std::cerr << "Not enough memory" << std::endl; break; case -3: status(INFO_INFO, "Reordering problem"); break; case -5: status(INFO_INFO,"unclassified error"); break; case -6: status(INFO_INFO, "preordering failed (matrix types 11, 13 only)"); break; case -7: status(INFO_INFO,"diagonal matrix problem"); break; case -8: status(INFO_INFO,"32 bit integer overflow problem"); break; case -9: status(INFO_INFO,"Not enough memory for OOC"); break; case -10: status(INFO_INFO,"Problems with the OOC data file"); break; case -11: status(INFO_INFO,"Read and write problems with OCC data"); break; } //cerr<<" Matrix is not positive definite. Retrying with non-positive definite solver in Reordering step"<