Intel® Fortran Compiler
Build applications that can scale for the future with optimized code designed for Intel® Xeon® and compatible processors.
29274 Discussions

Array Shapes and External Libraries

Dishaw__Jim
Beginner
772 Views
Suppose you have an array with the following shape

REAL :: A(2,512)

I'm calling the dss_solve_real function in the MKL, which has the following interface (defined in mkl_dss.f90)

INTERFACE
FUNCTION DSS_SOLVE_REAL( HANDLE, OPT, RRHSVALUES, NRHS, RSOLVALUES )
USE MKL_DSS_PRIVATE
TYPE(MKL_DSS_HANDLE), INTENT(INOUT) :: HANDLE
INTEGER(KIND=4), INTENT(IN) :: OPT
INTEGER(KIND=4), INTENT(IN) :: NRHS
REAL(KIND=8), INTENT(IN) :: RRHSVALUES( * )
REAL(KIND=8), INTENT(OUT) :: RSOLVALUES( * )
INTEGER(KIND=4) :: DSS_SOLVE_REAL
END FUNCTION DSS_SOLVE_REAL
END INTERFACE

At first I did the proper thing and declared a second array B(1024) and used the RESHAPE to get A into B and then would pass B as RSOLVALUES. Unfortunately the memory and copy penalty I was paying was getting expensive when the array sizes got larger. I decided to take advantage of how A was actually arranged in memory and decided to try to pass A as RSOLVALUES, which suprisingly enough works without even a compiler warning.

While the above method solves my performance and memory penalty problem, it strikes me as a kludge. The only solution that I know of is to define A as a 1-D array, though I loath to do that because the extra dimension simplifies some aspects of my code, specifically the first dimension specifies the direction and the second specifies the cell. I thought about creating a mapping function that would compute the array position based on a direction and cell number.

Is there a way to "reshape" an array without actually copying it or way to tell the compiler "yes I know I am trying to stick a 2-D array into a 1-D array"? Is there a downside to leaving my code as-is (other than it may fail to compile in the future)?
0 Kudos
6 Replies
Steven_L_Intel1
Employee
772 Views
Just pass A(1,1). It will do the right thing due to Fortran's sequence association rules.
0 Kudos
Dishaw__Jim
Beginner
772 Views
I tried your suggestion about passing A(1,1) and get an error. When I run in the IDE or from the command line, I get the following dialog box

Debug Error!

Program:...isual Studio est_transportdebug est_transport.exe

DAMAGE: after Normal block (#179) at 0x00369870

(Press Retry to debug the application)


The code seems to run okay--the error occurs at program termination.
0 Kudos
Steven_L_Intel1
Employee
772 Views
Interesting - sounds as if there's memory corruption. Are you sure you passed the correct size values to the external routine?
0 Kudos
Dishaw__Jim
Beginner
772 Views
I pass the numbers of rows and columns via the dss_define_structure call. Since I have square matrix, the call looks like

error = dss_define_structure(handle, &
MKL_DSS_NON_SYMMETRIC, &
rowIndex, nRows, nRows, columns, &
3 * nRows) ! Number of non-zero elements

Since I have only one right hand side, my call to dss_solve_real is

error = dss_solve_real(handle, MKL_DSS_DEFAULTS, Jextface, 1, Jface)
IF (error .NE. MKL_DSS_SUCCESS) THEN
WRITE(*,*) 'Unable to factor matrix, error = ', error
STOP
END IF

The error is generated at the end of the program when I deallocate Jface.
0 Kudos
Steven_L_Intel1
Employee
772 Views
Something is overwriting memory that doesn't belong to it. You'll have to find out what and where.
0 Kudos
Dishaw__Jim
Beginner
772 Views
I found it--thanks for the pointer.
0 Kudos
Reply