Silverfrost Forums

Welcome to our forums

How to loop over different matrix

10 Aug 2017 7:55 #19984

I want to use fortran90 to solve the following problem:

I have two different types of matrix:

A1(n1,1), A2(n2,1)......Ak(nk,1)

B1(n1,1), B2(n2,1)......Bk(nk,1)

I want to loop over i and calculate:

c'i' = A'i'+B'i' (i = 1....K)

What is the best way to do so?

Thanks

11 Aug 2017 1:06 #19988

There are two main answers to your question;

  1. is to 'loop' over your calculation with a DO loop,
  2. is to use 'Fortran 90' array syntax

The following (incomplete) Fortran code gives examples of both approaches:

! snippets of fortran code

  integer :: n1,n2,i,k
  real, allocatable :: a(:,:), b(:,:), c(:)

  n1 = ??
  n2 = ??
  allocate ( a(n1,n2), b(n1,n2), c(n1) )
!
! define values of a(:,1) and b(:,1)
!
! DO loop approach
  DO i = 1,k
    c(i) = a(i,1) + b(i,1)
  END DO
!
! Fortran 90 array syntax
  c(1:k) = a(1:k,1) + b(1:k,1)
!
! or if dimension 1 is 1:k, ie k=n1      
  c(:) = a(:,1) + b(:,1)

An explicit DO loop provides more control. There is no assumption on the order of array syntax, so if the calculation was c(i) =c(i-1) + a(i,1) + b(i,1), the array syntax approach may not be appropriate.

11 Aug 2017 2:57 #19997

Thank! This works well for arrays. However, if I have matix, I feel like this method is difficult to apply. I am wondering if there is any way in fortran that would allow me to loop over the index of the matrix (from matrix A1 to matrix Ak) without actually combining these matix into a huge one.

Thanks again!

12 Aug 2017 1:18 #20002

I am not sure I understand your question. 'matrix' can commonly refer to an array of rank 2. Array operations can also be applied to rank 2 arrays, which can give more options for array syntax, for example:

! snippets of fortran code 

  integer :: n1,n2,i,j,k 
  real, allocatable :: a(:,:), b(:,:), c(:,:) 

  n1 = ?? 
  n2 = ?? 
  allocate ( a(n1,n2), b(n1,n2), c(n1,n2) ) 
! 
! define values of a(:,1) and b(:,1) 
! 
! DO loop approach 
  DO j = 1,k 
   DO i = 1,k
    c(i,j) = a(i,j) + b(i,j) 
   END DO 
  END DO 
! 
! Fortran 90 array syntax 
  c(1:k,1:k) = a(1:k,1:k) + b(1:k,1:k) 
! 
! or if dimension 1&2 are 1:k, ie k=n1 and k=n2     
  c(:,:) = a(:,:) + b(:,:) 
!
! or simply
  c = a + b 

I write matrix multiplication using a mixture of array syntax and DO loops: subroutine stream_matmul_dp (a,b,c,nra,nca,ncb) use precision ! matrix multiplication : single thread, array syntax

       integer*4 nra,nca,ncb,  j,k
       real(dp) :: a(nra,nca), b(nca,ncb)
       real(dp) :: c(nra,ncb)
!
      do j = 1,ncb
        c(:,j) = 0
        do k = 1,nca
          c(1:nra,j) = c(1:nra,j) + a(1:nra,k) * b(k,j)
        end do
      end do
!
   end subroutine stream_matmul_dp

If you are concerned about the size of huge matrices, you should review recent posts that have referred to sparse matrices. However, there are many different ways of referencing non-zero terms in a matrix that utilise the non-zero patterns or groupings of the particular matrix structure.

Please login to reply.