Sub-Array Manipulations in Fortran

• Array Sub-objects (Array sections)

• Fortran has built-in support for array manipulation in support of Scientific Computing needs

• Some common operations in Numerical Linear Algebra are:

 matrix * vector vector * matrix

• To compute matrix * vector or vector * matrix , we must multiple:

 A slice of a matrix with a vector

Example:

• To implement these operations, we must be able to select slices in a matrix

• F90 has a number of specialized array selection methods (array "sub-objects") that will be discussed next....

• Selecting elements in one-dimensional arrays

• A subset of the elements in a one-dimensional array can be selected using one of these techniques:

• Explicit indices

Example:

 ``` REAL, DIMENSION(5) :: A A((/1, 3, 5/)) // Elements selected: A(1), A(3), A(5) ```

• Implicit DO loop

Example:

 ``` REAL, DIMENSION(5) :: A A(1:3) // Elements selected: A(1), A(2), A(3) A(1:5:2) // Elements selected: A(1), A(3), A(5) A(:) // Elements selected: A(1), A(2), A(3), A(4), A(5) ```

• Array as index

Example:

 ``` REAL, DIMENSION(5) :: A INTEGER, DIMENSION(3) :: K = (/ 1, 3, 5 /) A(K) // Elements selected: A(1), A(3), A(5) ```

• Example Program: (Demo above code)

• Selecting elements in multi-dimensional arrays

• You can use any of the above methods in each dimension in a multi-dimentional array

• Examples: selecting a slice in a matrix

 ``` REAL, DIMENSION(5,5) :: A INTEGER, DIMENSION(3) :: K = (/ 1, 3 /) A(2, (/1, 3, 5/)) --> A(2,1) A(2,3) A(2,5) A(2, 1:3) --> A(2,1) A(2,2) A(2,3) A(2, 1:3:2) --> A(2,1) A(2,3) A(2,5) A(2, :) --> A(2,1) A(2,2) A(2,3) A(2,4) A(2,5) A(2, K) --> A(2,1) A(2,3) ```

• Examples: selecting a blocks in a matrix

 ``` REAL, DIMENSION(5,5) :: A INTEGER, DIMENSION(3) :: K = (/ 1, 3 /) A((/1,3/), (/1,3/)) --> A(1,1) A(1,3) (2-dim. array !!!) A(3,1) A(3,3) A((/1,3/), 1:3) --> A(1,1) A(1,2) A(1,3) A(3,1) A(3,2) A(3,3) ```

• Implementing the Matrix*Vector operation with array operations

• The matrix-vector multiplication consists of a number of

 array multiplications , and summarion of array elements

• Consider the Matrix-vector multiplication:

The first element in the result vector (w(1)) is obtained by:

 ``` 1. multiply first row of matrix and the input vector: A(1, :) * v 2. sum all elements of the result of the first step: SUM ( A(1, :) * v ) ```

• Using F90 slice selection techniques , we can write the Matrix-vector multiplication as:

 ``` INTEGER N = ... REAL, DIMENSION(N,N) :: A INTEGER, DIMENSION(N) :: v, w w(1) = SUM ( A(1, : ) * v ) w(2) = SUM ( A(2, : ) * v ) .. w(N) = SUM ( A(N, : ) * v ) or: DO i = 1, N w(i) = SUM ( A(i, : ) * v ) END DO ```

• Example Program: (Demo above code)

• Alternative implementation of the Matrix*vector operation with array operations

• There is a better way to implement the Matrix-vector multiplication

Consider the following example:

We can see that the Matrix-vector multiplication can be computed as:

 ``` w = A(:,1)*v(1) + A(:,2)*v(2) + A(:,3)*v(3) ```

• The following DO-loop will perform a Matrix-vector multiplication:

 ``` real, dimension(3,3) :: A real, dimension(3) :: v real, dimension(3) :: w integer :: i integer :: N N = size(v) w = 0.0 !! clear whole vector DO i = 1, N w = w + A( :, i ) * v(i) END DO ```

• Hint: use array operation whenever possible

• When running on multi-processor computers (computers with multiple CPUs), the Fortran compiler can usually parallel array operations

• The work to perform an array operations is distributed among the processors

• You will achieve a significant speed increase using array operations whenever possible (if you are using a multip-processor system and a parallizing compiler).

• Reshaping arrays...

• A number of F90 intrinsic functions allow the programmer to reshape arrays (matrices) :

 Transpose Pack Reshape

• Transposing arrays (matrices)...

• Function call to transpose a matrix:

 ``` REAL, DIMENSION(3, 3) :: A, B B = TRANSPOSE(A) ```

• Example Program: (Demo above code)

• Packing arrays (matrices)...

• Packing an array is to collapse a multi-dimensional array into a one-dimensional array

• The PACK intrinsic function uses a logical array to select elements for "packing"

Example:

 ``` REAL, DIMENSION(3, 3) :: A LOGICAL, DIMENSION(3, 3) :: MASK REAL, DIMENSION(4) :: B DATA MASK / .true., .true., .false., & .false., .false., .false., & .false., .true., .true. / B = PACK(A, MASK) Result: B(1) = A(1,1) ! MASK(1) == .true. B(2) = A(1,2) ! MASK(2) == .true. B(3) = A(3,2) ! MASK(8) == .true. B(4) = A(3,3) ! MASK(9) == .true. ```

• Example Program: (Demo above code)

• NOTE:

To select every element for packing, use:

 ``` B = PACK(A, .TRUE.) ```

F90 will automatically convert the logical value .true. into a conformable logical array

• Reshaping one-dimensional arrays into multi-dimensional arrays

• Reshape: re-distribute the elements in a one-dimensional array into a multi-dimensional array

• The RESHAPE function uses a shape array to construct the multi-dimensional array

• Example:

 ``` REAL, DIMENSION(10) :: A ! Source array REAL, DIMENSION(2, 5) :: B ! Destination INTEGER, DIMENSION(2) :: X = (/ 2, 5 /) ! Shape: 2 rows, 5 columns B = RESHAPE(A, X) Result: B(1,1) = A(1) B(1,2) = A(3) ... B(2,1) = A(2) B(2,2) = A(4) ... ```

• Example Program: (Demo above code)