64-bit interface PARDISO and PARDISO_64

[tlcfem]

The signature of `pardiso` in `mkl_pardiso.h` reads

```c

void pardiso( _MKL_DSS_HANDLE_t pt, const MKL_INT *maxfct, const MKL_INT *mnum,
const MKL_INT *mtype, const MKL_INT *phase, const MKL_INT *n,
const void *a, const MKL_INT *ia, const MKL_INT *ja,
MKL_INT *perm, const MKL_INT *nrhs, MKL_INT *iparm,
const MKL_INT *msglvl, void *b, void *x, MKL_INT *error );

void pardiso_64( _MKL_DSS_HANDLE_t pt, const long long int *maxfct,
const long long int *mnum, const long long int *mtype,
const long long int *phase, const long long int *n,
const void *a, const long long int *ia,
const long long int *ja, long long int *perm,
const long long int *nrhs, long long int *iparm,
const long long int *msglvl, void *b, void *x,
long long int *error );

```

When using ilp64 interface, `MKL_INT` is `long long`, making the signatures indifferent.

Is this intended?

It makes more sense that `pardiso` is using 32-bit int, while `pardiso_64` uses 64-bit int, otherwise why two functions?

If it is intended and correct, does it mean `pardiso` can be called when ilp64 is used?

If it is not intended, shall `MKL_INT` be replaced by `int`?

[Spencer_P_Intel]

Howdy,

In the last few years, we have been introducing the xyz_64() interfaces, to be compared with the xyz() interface (xyz = pardiso in your case), with the intent that MKL_INT is used for integer arguments in xyz() APIs and MKL_INT64 is used for xyz_64().  Practically this means that for ilp64 linking, the arguments are exactly the same, as MKL_INT is a signed 64 bit integer with ilp64 and a signed 32 bit integer with lp64 and MKL_ILP64 is always a signed 64 bit integer. 

If you happen to notices that the arguments are not the same for ilp64, then please let us know and it may be a bug on our part....  

The real useful intent was to allow users linking with oneMKL using lp64 to use 32 bit integers for calls to xyz and 64 bit integers for calls to xyz_64 in the same application.  This is especially important in sparse linear algebra (like pardiso) where the dimension and number of elements stored in the sparse matrix determines which integer type makes the most sense to use to represent it, with the recommendation to pick the smallest available integer type that can represent your matrix workload for get better performance.

So yes, in ilp64, they are essentially the same, but in lp64, it opens up a wider range of solutions.

Best Regards,

Spencer

[tlcfem]

Thank you for the explanation.

Two follow-up questions.

1. I presume it is then logical to use one single function in CPP templates, for example,

```cpp

template<int_t T> void solver_wrapper(...) {

T mtype;

T msglvl;

...

pardiso(...);

}

```

For T being either `int` or `long long`, assuming linking proper interface.

2. The documentation states: `The pardisoinit routine cannot be used together with the pardiso_64 routine.`, See: https://www.intel.com/content/www/us/en/docs/onemkl/developer-reference-c/2025-1/pardisoinit.html#GUID-C16BFBCA-EF1C-4CDB-BC73-41655B6DD8F5

If pardiso and pardiso_64 are identical with ilp64, what makes `pardisoinit` not accessible?

The signature is

```c

void pardisoinit( _MKL_DSS_HANDLE_t pt, const MKL_INT *mtype, MKL_INT *iparm );

```

[tlcfem]

Or maybe it means `pardisoinit` cannot be used to initialize buffers for `pardiso_64` with lp64 interface? (Guess it will only allocate initialize half the size required?)

For ilp64 it is just fine?

[c_sim · Accepted solution]

Hi,

yes, you are right. You can use pardisoinit with pardiso_64, if using ILP64 linking. However, with LP64 linking, using pardisoinit and pardiso_64 together is incorrect and might even result in a segfault. In the future we might add pardisoinit_64, but currently we don't have it. As a workaround you could initialize the handles and iparm manually in this case, i.e., instead of pardisoinit you can do something like:

for ( i = 0; i < 64; i++ )
{
     pt[i] = 0;
}
for ( i = 0; i < 64; i++ )
{
     iparm[i] = 0;
}
iparam[0] = 1;
//Fill iparm values
 iparam[1] = 2;
 ...

For default values of iparm, please refer to pardiso iparm Parameter. For examples on this kind of initialization please refer to oneMKL examples (in share/doc/mkl/examples/examples_core_c.tgz)

Regarding question 1 on templating and pardiso() function: Yes, you can use templating strategy for integer type and be careful with linking. However, since pardiso() function just uses one integer type depending on linking, it might not be that useful. It might be better to use MKL_INT instead.

The other possibility is to use pardiso() and pardiso_64() together with lp64 interface. With this approach, you could template for integer type and select between pardiso() and pardiso_64() variant depending on the integer type.

Hope this answers your questions.

Kind Regards,

Chris

[tlcfem]

Thank you both Spencer and Chris!