Matrix Operations with DO Loops in Fortran

In Fortran, matrices are represented as two-dimensional arrays, which are widely used in scientific computing, engineering simulations, numerical analysis, and data processing. Performing operations on matrices often requires nested loops to traverse rows and columns.

This post explores how to access and manipulate 2D arrays (matrices) using DO loops, perform arithmetic operations, implement common matrix algorithms, and follow best practices for readability and efficiency.

1. Introduction to 2D Arrays

A 2D array (matrix) is declared with two dimensions, typically representing rows and columns:

real :: matrix(3,3)
  • The first index refers to the row
  • The second index refers to the column
  • Fortran stores arrays in column-major order, meaning elements of a column are stored contiguously in memory.

Example: Initializing a Matrix

program initialize_matrix
real :: matrix(3,3)
matrix = reshape((/1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0/), (/3,3/))
print *, "Matrix:"
print *, matrix
end program initialize_matrix

Explanation:

  • reshape converts a 1D list into a 2D matrix
  • The matrix now has 3 rows and 3 columns

2. Accessing Matrix Elements

Matrix elements are accessed using row and column indices:

matrix(i,j)
  • i ranges from 1 to the number of rows
  • j ranges from 1 to the number of columns

Example: Print Individual Elements

program print_matrix_elements
integer :: i, j
real :: matrix(3,3)
matrix = reshape((/1,2,3,4,5,6,7,8,9/), (/3,3/))
do i = 1, 3
    do j = 1, 3
        print *, "Element (", i, ",", j, ") =", matrix(i,j)
    end do
end do
end program print_matrix_elements

Explanation:

  • Nested loops traverse all rows and columns
  • Each element is accessed individually using matrix(i,j)

3. Basic Arithmetic Operations on Matrices

DO loops can be used to perform arithmetic operations element-wise.

Example 1: Doubling Each Element

program double_matrix
integer :: i, j
real :: matrix(3,3)
matrix = reshape((/1,2,3,4,5,6,7,8,9/), (/3,3/))
do i = 1, 3
    do j = 1, 3
        matrix(i,j) = matrix(i,j) * 2
    end do
end do
print *, "Doubled matrix:"
do i = 1, 3
    print *, matrix(i,1:3)
end do
end program double_matrix

Output:

Doubled matrix:
2 4 6
8 10 12
14 16 18

Explanation:

  • Outer loop iterates over rows, inner loop iterates over columns
  • Each element is multiplied by 2
  • Result is printed row by row

4. Matrix Addition

You can perform matrix addition by looping over corresponding elements of two matrices.

Example 2: Adding Two Matrices

program matrix_addition
integer :: i, j
real :: A(2,2), B(2,2), C(2,2)
A = reshape((/1,2,3,4/), (/2,2/))
B = reshape((/5,6,7,8/), (/2,2/))
do i = 1, 2
    do j = 1, 2
        C(i,j) = A(i,j) + B(i,j)
    end do
end do
print *, "Matrix A + B:"
do i = 1, 2
    print *, C(i,1:2)
end do
end program matrix_addition

Explanation:

  • Each element of C is the sum of corresponding elements in A and B
  • Nested loops allow row-wise and column-wise access

5. Matrix Subtraction

Matrix subtraction is similar to addition, using loops to subtract element-wise.

Example 3: Subtracting Two Matrices

program matrix_subtraction
integer :: i, j
real :: A(2,2), B(2,2), C(2,2)
A = reshape((/5,6,7,8/), (/2,2/))
B = reshape((/1,2,3,4/), (/2,2/))
do i = 1, 2
    do j = 1, 2
        C(i,j) = A(i,j) - B(i,j)
    end do
end do
print *, "Matrix A - B:"
do i = 1, 2
    print *, C(i,1:2)
end do
end program matrix_subtraction

Explanation:

  • Each element of C contains the difference of corresponding elements
  • Nested loops provide structured element access

6. Matrix Multiplication

Matrix multiplication requires three nested loops:

  1. Outer loop for rows of the first matrix
  2. Middle loop for columns of the second matrix
  3. Inner loop for summing the products

Example 4: Matrix Multiplication

program matrix_multiplication
integer :: i, j, k
real :: A(2,3), B(3,2), C(2,2)
A = reshape((/1,2,3,4,5,6/), (/2,3/))
B = reshape((/7,8,9,10,11,12/), (/3,2/))
C = 0.0
do i = 1, 2
    do j = 1, 2
        do k = 1, 3
            C(i,j) = C(i,j) + A(i,k) * B(k,j)
        end do
    end do
end do
print *, "Matrix A * B:"
do i = 1, 2
    print *, C(i,1:2)
end do
end program matrix_multiplication

Explanation:

  • Uses three nested loops to calculate dot product for each element of C
  • Outer and middle loops select the position in the result matrix
  • Inner loop performs sum of products

7. Transpose of a Matrix

Transposing a matrix swaps rows and columns. Nested loops can achieve this efficiently.

Example 5: Matrix Transpose

program matrix_transpose
integer :: i, j
real :: A(2,3), B(3,2)
A = reshape((/1,2,3,4,5,6/), (/2,3/))
B = 0.0
do i = 1, 2
    do j = 1, 3
        B(j,i) = A(i,j)
    end do
end do
print *, "Transpose of A:"
do i = 1, 3
    print *, B(i,1:2)
end do
end program matrix_transpose

Explanation:

  • Loops iterate over rows and columns
  • Assigns B(j,i) = A(i,j) to transpose
  • Nested loops allow precise element mapping

8. Scaling a Matrix

Scaling a matrix multiplies each element by a constant.

Example 6: Scaling Matrix Elements

program scale_matrix
integer :: i, j
real :: matrix(3,3)
real :: factor
matrix = reshape((/1,2,3,4,5,6,7,8,9/), (/3,3/))
factor = 0.5
do i = 1, 3
    do j = 1, 3
        matrix(i,j) = matrix(i,j) * factor
    end do
end do
print *, "Scaled matrix:"
do i = 1, 3
    print *, matrix(i,1:3)
end do
end program scale_matrix

Explanation:

  • Each element multiplied by factor
  • Demonstrates element-wise operations using nested loops

9. Element-wise Conditional Operations

You can combine nested loops with IF statements to perform conditional operations.

Example 7: Set Negative Elements to Zero

program conditional_matrix
integer :: i, j
real :: matrix(3,3)
matrix = reshape((/1,-2,3,-4,5,-6,7,-8,9/), (/3,3/))
do i = 1, 3
    do j = 1, 3
        if (matrix(i,j) < 0) then
            matrix(i,j) = 0
        end if
    end do
end do
print *, "Matrix after setting negatives to zero:"
do i = 1, 3
    print *, matrix(i,1:3)
end do
end program conditional_matrix

Explanation:

  • Applies conditional logic to each element
  • Useful for filtering, thresholding, or preprocessing data

10. Advanced Example: Sum of Rows and Columns

Nested loops can be used to calculate row sums and column sums.

program sum_rows_columns
integer :: i, j
real :: matrix(3,3), row_sum(3), col_sum(3)
matrix = reshape((/1,2,3,4,5,6,7,8,9/), (/3,3/))
row_sum = 0.0
col_sum = 0.0
do i = 1, 3
    do j = 1, 3
        row_sum(i) = row_sum(i) + matrix(i,j)
        col_sum(j) = col_sum(j) + matrix(i,j)
    end do
end do
print *, "Row sums:", row_sum
print *, "Column sums:", col_sum
end program sum_rows_columns

Explanation:

  • Outer loop iterates rows, inner loop iterates columns
  • Row sum accumulates values in each row
  • Column sum accumulates values in each column

11. Best Practices for Matrix Operations

  1. Use nested loops systematically: Outer loop for rows, inner for columns.
  2. Initialize arrays before use to prevent undefined behavior.
  3. Use reshape for multi-dimensional initialization.
  4. Apply element-wise operations carefully for clarity.
  5. Use descriptive variable names for matrices and loop indices.
  6. Consider array operations without loops when applicable (Fortran supports array arithmetic).
  7. Comment complex operations for readability.

12. Practical Applications

  • Linear algebra: Addition, subtraction, multiplication, transpose
  • Scientific simulations: Grids, tensors, spatial computations
  • Data analysis: Summing rows/columns, scaling, filtering
  • Computer graphics: Pixel arrays, transformations
  • Engineering computations: Stress/strain matrices, physical simulations

13. Comprehensive Example

program comprehensive_matrix_operations
integer :: i, j, n
real :: matrix(3,3), doubled(3,3), transpose(3,3), row_sum(3), col_sum(3)
! Initialize matrix
matrix = reshape((/1,2,3,4,5,6,7,8,9/), (/3,3/))
! Double each element
doubled = matrix
do i = 1, 3
    do j = 1, 3
        doubled(i,j) = doubled(i,j) * 2
    end do
end do
! Transpose the matrix
do i = 1, 3
    do j = 1, 3
        transpose(j,i) = doubled(i,j)
    end do
end do
! Calculate row and column sums
row_sum = 0.0
col_sum = 0.0
do i = 1, 3
    do j = 1, 3
        row_sum(i) = row_sum(i) + doubled(i,j)
        col_sum(j) = col_sum(j) + doubled(i,j)
    end do
end do
! Print results
print *, "Original matrix:"
do i = 1, 3
    print *, matrix(i,1:3)
end do
print *, "Doubled matrix:"
do i = 1, 3
    print *, doubled(i,1:3)
end do
print *, "Transpose of doubled matrix:"
do i = 1, 3
    print *, transpose(i,1:3)
end do
print *, "Row sums:", row_sum
print *, "Column sums:", col_sum
end program comprehensive_matrix_operations

Explanation:

  • Demonstrates nested loops, arithmetic, transpose, and sum operations
  • Combines multiple matrix operations in a single program
  • Illustrates structured and readable Fortran matrix processing

Comments

Leave a Reply

Your email address will not be published. Required fields are marked *