topic Interesting thank you - I in Intel® oneAPI Math Kernel Library

PARDISO b array

JohnNichols — Fri, 12 Sep 2014 15:33:13 GMT

In most problems involving a matrix inversion, one only looks at a single b array in the Ax=b problem.

PARDISO has however an arrangement, where b is actually b[*,nhrs] in both the Geneva and the MKL versions. I checked the manuals.

But in the sample, just supplied by Intel, the b matrix in the C# example is called as a vector b

and passed to mkl as such. Clearly it runs and gives the correct answer, but why does the compiler not flag an error? This occurs in the other samples from other places, including the latest GENEVA sample

int mtype = 11; /* Real unsymmetric matrix */
        /* RHS and solution vectors. */
        double[] bData = new double;
        bData[101] = 450;
        AlignedDoubleArrayPointer b = new AlignedDoubleArrayPointer(bData);

Clearly not an ideal solution that one misses the second index on the array?

JMN

What is passed from Fortran

mecej4 — Fri, 12 Sep 2014 16:30:00 GMT

What is passed from Fortran is the base address of the array, and any of the dimensions of the array that are needed by the called subroutine to compute addresses from a subscripted array reference (1-D, 2-D, etc.). Thus, 'x' as a Fortran call argument is the same as another call with x(1) or x(1,1), etc. In the specific case of solving with multiple right hand sides, following a successful factorization, each r.h.s. is a column in the matrix, and that should throw some light on the reason for using column storage mode for matrices in Fortran.

Unless you force strict interface checking (by providing an interface or asking for interface checking), most library subroutines designed to be called from Fortran work correctly when called using the calling sequence in the documentation.

Some of the things I wrote above are no longer valid if you use assumed shape arrays or have VALUE arguments -- features that did not exist in Fortran 77 and earlier.

Interesting thank you - I

JohnNichols — Fri, 12 Sep 2014 19:26:16 GMT

Interesting thank you - I tried adding a second column and the program would not compile, I tried adding extra elements to the single dimension vector b and telling the program I had 2 b sets, and got an exception on the second one.

Thankfully I do not need it - I will leave well enough alone.

JMN

Quote:John Nichols wrote: I

mecej4 — Sun, 14 Sep 2014 13:30:01 GMT

John Nichols wrote:
I tried adding a second column and the program would not compile

If you did that in C or C#, what I said about Fortran conventions will not apply without suitable modifications.

Hi John,

Ying_H_Intel — Mon, 15 Sep 2014 08:45:55 GMT

Hi John,

Not sure how you add a second column and the program would not compile, I tried adding extra elements to the single dimension vector b and telling the program I had 2 b sets.

The function pardiso support both Fortran and C interface, so it can accept that 2-D array as 1-D memory-continuous array and nrhs is used for second index. for example you can handle the case 2b sets as below

nrhs= 2

	
double[] bData = new double[n *nrhs];

	
				4
				        bData[101] = 450;

		
				5
				        AlignedDoubleArrayPointer b = new AlignedDoubleArrayPointer(bData);
			
Best Regards,
Ying