Segmentation fault (core dumped) in PARDISO()?

[SAL_AHM]

I am using the pardiso mkl to solve a system of equations A*X=B. My matrix is in CRS and save in vectors. So I changed to array as:

MKL_INT n = f.size();//or   _id.size()-1;
MKL_INT nnz = _sk.size();
    MKL_INT ia[n+1];
    cout<<"Salman Ahmad: 0"<<endl;
    transform(_id.begin(), _id.end(), ia,[](const int & elem){return elem;});
    MKL_INT ja[nnz];
    cout<<"Salman Ahmad: 1"<<endl;
    transform(_ik.begin(), _ik.end(), ja,[](const int & elem){return elem;});
    double a[nnz];
    cout<<"Salman Ahmad: 2"<<endl;
    transform(_sk.begin(), _sk.end(), a,[](const double & elem){return elem;});

But now they give me the segmentation fault on the function:

PARDISO (pt, &maxfct, &mnum, &mtype, &phase,
		&n, a, ia, ja, &idum, &nrhs,
		iparm, &msglvl, &ddum, &ddum, &error);

Please help me, how to handle this problem?

[VarshaS_Intel]

Hi,

Thanks for posting in Intel Communities.

Could you please cross-verify if the linking of your application is done as expected according to the format of integers that were been used by you?

Could you please refer to the below link for more details regarding the usage of link line advisor:

https://www.intel.com/content/www/us/en/developer/tools/oneapi/onemkl-link-line-advisor.html

In addition, we would like to request you to kindly provide us with the environment details, sample reproducer, and complete error you are facing as it helps us reproduce the issue in our environment and assist you accordingly.

It is needed to choose the LP64 interface or ILP64 interface as the difference between them is integer type length i.e. ILP64 interface uses the 64-bit integer type, while LP64 uses the 32-bit integer type.

Thanks & Regards,

Varsha

[SAL_AHM]

ff

[SAL_AHM]

Dear Varsha,

  Thank you so much for your reply. I am the new user so sorry for mistakes in advance. I am trying to explain you my problem as: I am coding in C++, before I used the GCC compiler but it is not working with intel, give me the error of undefined reference. So now I use the ICC compiler. I defined the pardiso function below in which the _id vector is extents of rows, _ik column index vector and _sk is nonzero vector.   #include <stdio.h>
#include <stdlib.h>
#include <math.h>
// PARDISO prototype.
#if defined(_WIN32) || defined(_WIN64)
#define pardiso_ PARDISO
#else
#define PARDISO pardiso_
#endif
#if defined(MKL_ILP64)
#define MKL_INT long long
#else
#define MKL_INT int
#endif
extern "C" MKL_INT PARDISO
(void *, MKL_INT *, MKL_INT *, MKL_INT *, MKL_INT *, MKL_INT *,
double *, MKL_INT *, MKL_INT *, MKL_INT *, MKL_INT *, MKL_INT *,
MKL_INT *, double *, double *, MKL_INT *);   template <class ELEM>
std::vector<double> FEM_Matrix_2<ELEM>::Solve_ pardiso_mkl(std::vector< double> const &f) const
{
MKL_INT n = f.size();//or   _id.size()-1; // f is right hand side data vector
MKL_INT nnz = _sk.size();
    MKL_INT ia[n+1];
    std::copy(_id.begin(), _id.end(), ia);
    //for(int i=0; i<n+1;++i){cout<<"ia: "<<ia[i]<<endl;}
   
    MKL_INT ja[nnz];
    std::copy(_ik.begin(), _ik.end(), ja);
    double a[nnz];
    cout<<"Salman Ahmad: 2"<<endl;
    std::copy(_sk.begin(), _sk.end(), a);
    for(int i=0; i<nnz;++i){cout<<"ik: "<<ja[i]<<", sk: "<<a[i]<<endl;}
   
MKL_INT mtype = 11; /* Real unsymmetric matrix */
/* RHS and solution vectors. */
double b[n], x[n];
MKL_INT nrhs = 1; /* Number of right hand sides. */
/* Internal solver memory pointer pt, */
/* 32-bit: int pt[64]; 64-bit: long int pt[64] */
/* or void *pt[64] should be OK on both architectures */
void *pt[64];
/* Pardiso control parameters. */
MKL_INT iparm[64];
MKL_INT maxfct, mnum, phase, error, msglvl;
/* Auxiliary variables. */
MKL_INT i;
double ddum; /* Double dummy */
MKL_INT idum; /* Integer dummy. */
/* ------------------------------ ------------------------------ -------- */
/* .. Setup Pardiso control parameters. */
/* ------------------------------ ------------------------------ -------- */
for (i = 0; i < 64; i++) {
iparm[i] = 0;
}
iparm[0] = 1; /* No solver default */
iparm[1] = 2; /* Fill-in reordering from METIS */
/* Numbers of processors, value of OMP_NUM_THREADS */
iparm[2] = 1;
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; /* Not in use */
iparm[7] = 2; /* Max numbers of iterative refinement steps */
iparm[8] = 0; /* Not in use */
iparm[9] = 13; /* Perturb the pivot elements with 1E-13 */
iparm[10] = 1; /* Use nonsymmetric permutation and scaling MPS */
iparm[11] = 0; /* Not in use */
iparm[12] = 0; /* Not in use */
iparm[13] = 0; /* Output: Number of perturbed pivots */
iparm[14] = 0; /* Not in use */
iparm[15] = 0; /* Not in use */
iparm[16] = 0; /* Not in use */
iparm[17] = -1; /* Output: Number of nonzeros in the factor LU */
iparm[18] = -1; /* Output: Mflops for LU factorization */
iparm[19] = 0; /* Output: Numbers of CG Iterations */
maxfct = 1; /* Maximum number of numerical factorizations. */
mnum = 1; /* Which factorization to use. */
msglvl = 1; /* Print statistical information in file */
error = 0; /* Initialize error flag */
/* ------------------------------ ------------------------------ -------- */
/* .. Initialize the internal solver memory pointer. This is only */
/* necessary for the FIRST call of the PARDISO solver. */
/* ------------------------------ ------------------------------ -------- */
for (i = 0; i < 64; i++) {
pt[i] = 0;
}
/* ------------------------------ ------------------------------ -------- */
/* .. Reordering and Symbolic Factorization. This step also allocates */
/* all memory that is necessary for the factorization. */
/* ------------------------------ ------------------------------ -------- */

phase = 11;
PARDISO(pt, &maxfct, &mnum, &mtype, &phase,
&n, a, ia, ja, &idum, &nrhs,
iparm, &msglvl, &ddum, &ddum, &error);
if (error != 0) {
printf("\nERROR during symbolic factorization: %d", error);
exit(1);
}
printf("\nReordering completed ... ");
printf("\nNumber of nonzeros in factors = %d", iparm[17]);
printf("\nNumber of factorization MFLOPS = %d", iparm[18]);
/* ------------------------------ ------------------------------ -------- */
/* .. Numerical factorization. */
/* ------------------------------ ------------------------------ -------- */
phase = 22;
PARDISO(pt, &maxfct, &mnum, &mtype, &phase,
&n, a, ia, ja, &idum, &nrhs,
iparm, &msglvl, &ddum, &ddum, &error);
if (error != 0) {
printf("\nERROR during numerical factorization: %d", error);
exit(2);
}
printf("\nFactorization completed ... ");
/* ------------------------------ ------------------------------ -------- */
/* .. Back substitution and iterative refinement. */
/* ------------------------------ ------------------------------ -------- */
phase = 33;
iparm[7] = 2; /* Max numbers of iterative refinement steps. */
/* Set right hand. */
    std::copy(f.begin(), f.end(), b);
//for (i = 0; i < n; i++) {
// b[i] = 1;
//}
PARDISO(pt, &maxfct, &mnum, &mtype, &phase,
&n, a, ia, ja, &idum, &nrhs,
iparm, &msglvl, b, x, &error);
if (error != 0) {
printf("\nERROR during solution: %d", error);
exit(3);
}
printf("\nSolve completed ... ");
printf("\nThe solution of the system is: ");
for (i = 0; i < n; i++) {
printf("\n x [%d] = % f", i, x[i] );
}
printf ("\n");
/* ------------------------------ ------------------------------ -------- */
/* .. Termination and release of memory. */
/* ------------------------------ ------------------------------ -------- */
phase = -1; /* Release internal memory. */
PARDISO(pt, &maxfct, &mnum, &mtype, &phase,
&n, &ddum, ia, ja, &idum, &nrhs,
iparm, &msglvl, &ddum, &ddum, &error);



std::vector<double> sol(f.size());// I will save the solution in this

    return sol;
}   This is working for matrix data in ( http://sep.stanford.edu/sep/claudio/Research/Prst_ExpRefl/ShtPSPI/intel/mkl/10.0.3.020/examples/solver/source/pardiso_unsym_c.c)   

But if I change to my data like above it gives error:  

Segmentation fault (core dumped).    

Best Regards,

Sal_Ahm

[VarshaS_Intel]

Hi,

Thanks for the information.

When we tried to compile your code at our end, we are getting the error below:

sample.cpp:21:43: error: no template named 'FEM_Matrix_2'

template <class ELEM> std::vector<double> FEM_Matrix_2<ELEM>::Solve_ pardiso_mkl(std::vector< double> const &f) const

                     ^

1 error generated.

which means there will be other code where you have defines the "FEM_Matrix_2". Could you please provide us with the complete sample reproducer code?

And also, could you please let us know the OS details and Intel oneAPI version you are using? 

You can also find the example code for Intel MKL(Linux- /opt/intel/oneapi/mkl/2023.0.0/examples/examples_core_c.zip , Windows- C:\Program Files (x86)\Intel\oneAPI\mkl\2023.0.0\examples\examples_core_c.zip\c)

Thanks & Regards,

Varsha

[SAL_AHM]

Dear Vaesha,

Thank you so much for your support.  I am using 2023.0.0 oneAPI and my OS is linux with specification is 768 GB ram & dual Xeon 24-core CPUs.

The code is in the attachment. The Solve_parsido_mkl(...) is in element.h (line: 2199-2353) and call in main.cpp (line: 182). 

I also did a separate post for reordering problem but I this this is the same issue. 

Best Regards,

SAL_AHM

[SAL_AHM]

Dear Vaesha,

I changed the compiler to dpcpp the detail given in ONEAPI_default.mk file. But still they give me the error:

ERROR during symbolic factorization: -1

The source code is in the attachment.   Please you can check this where is the problem?

Best Regards,

Salman Ahmad

[SAL_AHM]

Hi Varsha,

I fixed the problems now but still have one problem. If I increase the number of equations such that matrix size then it gives me the Segmentation fault (core dumped).

Please can you check the attached code in which the number of equations is 17011?

Best Regards,

Sal_Ahm

[VarshaS_Intel]

Hi,

Thanks for the details.

Could you please let us know the other libraries that are included

SUPERLU_INC= -I/usr/include/superlu -I/usr/include/superlu-dist -I/home/asalman/cuhk_ahmad_2D/navier_stokes_2D/SRC

As we are facing an error here searched for the Navier stokes 2D but we could not find the exact one.

Could you please help us with how we can get the libraries? Also, could you please share with us the complete steps you have followed in order to investigate more on our end?

Thanks & Regards,

Varsha

[SAL_AHM]

Dear Versha,

Thank you so much for your support. 

SUPERLU_INC= -I/usr/include/superlu -I/usr/include/superlu-dist -I/home/asalman/cuhk_ahmad_2D/navier_stokes_2D/SRC

is a SuperLU library used for solutions but now it is closed, no need anymore. Find the attached code without this library. Maybe the problem is the size of the array. Also I tried to save the matrix data in a .CSR file (in elements.h line: 2207-2237) and call this function in main.cpp line 173. But unfortunately I don't understand how to call ia[], ja[] and a[] from the .CRS file. 

Best Regards,

SAL_MAN

[VarshaS_Intel]

Hi,

We could see similar behavior at our end, if possible could you please isolate your issue with minimal sample reproducer code so that we can investigate more on the root cause?

Thanks & Regards,

Varsha

[SAL_AHM]

Dear Versha,

Thank you so much for your support.

I saved my data in .txt files and call in .cpp file for ia[], ja[], a[] and b[]. So please find the attachment.

Best Regards,

Sal_Ahm

[VarshaS_Intel]

Hi,

Thank you for providing the reproducer.

Based on the error you encountered, it seems like there is a mistake in line 135 of the sample code. To fix this error, please replace line 135 with the sample code:

"b = (double*)MKL_malloc(n*sizeof(double),128);"

We did not get a "Segmentation Fault" error. Could you please try it and let us know if you face any further issues?

Thanks & Regards,

Varsha

[SAL_AHM]

Dear Versha,

Thank you so much but there is still a problem. It runs without error but does not give the solution. I noted as in:

Statistics:
===========
Parallel Direct Factorization is running on 48 OpenMP

< Linear system Ax = b >
number of equations: 155443
number of non-zeros in A: 0
number of non-zeros in A (%): 0.000000

number of right-hand sides: 1

< Factors L and U >
number of columns for each panel: 72
number of independent subgraphs: 0
< Preprocessing with state of the art partitioning metis>
number of supernodes: 155443
size of largest supernode: 1
number of non-zeros in L: 155443
number of non-zeros in U: 1
number of non-zeros in L+U: 155444
gflop for the numerical factorization: 0.000000

gflop/s for the numerical factorization: 0.000000

Salman Ahmad: 5

Factorization completed ...
=== PARDISO: solving a real nonsymmetric system ===

Summary: ( solution phase )
================

Times:
======
Time spent in direct solver at solve step (solve) : 0.120024 s
Time spent in additional calculations : 0.063136 s
Total time spent : 0.183160 s

Statistics:
===========
Parallel Direct Factorization is running on 48 OpenMP

< Linear system Ax = b >
number of equations: 155443
number of non-zeros in A: 0
number of non-zeros in A (%): 0.000000

number of right-hand sides: 1

< Factors L and U >
number of columns for each panel: 72
number of independent subgraphs: 0
< Preprocessing with state of the art partitioning metis>
number of supernodes: 155443
size of largest supernode: 1
number of non-zeros in L: 155443
number of non-zeros in U: 1
number of non-zeros in L+U: 155444
gflop for the numerical factorization: 0.000000

gflop/s for the numerical factorization: 0.000000

Solve completed ...

It gives me the number of non-zeros in A: 0. Why this is not taking the input data from .txt files?

Best Regards,

Salm_Ahm

[VarshaS_Intel]

Hi,

Could you please let us know the significance of using "aman" in the code(fopen("file_ia.txt","aman"))? As it helps us in analyzing the issue better.

Thanks & Regards,

Varsha

[SAL_AHM]

Hello,

 there is no specific significance for "aman" , just only want to read the .txt file. Instead of it only can use "r".  

Best Regards,

Sal_Ahm

[VarshaS_Intel]

Hi,

>>It gives me the number of non-zeros in A: 0.

Could you please confirm whether you are observing the same behavior with any other datasets?

Thanks & Regards,

Varsha

[VarshaS_Intel · Accepted solution]

Hi,

Please find the attached code where we are able to get the correct results without zeros. Could you please try and let us know if you are still facing any issues?

Thanks & Regards,

Varsha

[SAL_AHM]

Hello,

Thank you so much for your support. Now it is run correctly but there is still one problem with the residual error. It gives me very high residual error given below:

Input and solution norms:
||A|| = 3.16228e+07
||b|| = 6.65132e-05
||x|| = 0.397544
||Y|| = 171350
||Ax-b|| = 171350

Can you check it and suggest some solution? 

Best Regards,

Sal_Ahm

[VarshaS_Intel]

Hi,

Since your initial query was resolved and we could see that you have raised a new thread https://community.intel.com/t5/Intel-oneAPI-Math-Kernel-Library/Why-residual-error-is-so-high-in-PARDISO/m-p/1479516. We will address your queries in that thread. So, we are going ahead and closing this thread.

Thanks & Regards,

Varsha