program fortran_introduction

c..   This program is to give an introduction to FORTRAN
c..   We will make some matrix operation, e.g., multiplication,
c..   addition, substraction, etc.

c..   We need to define arrays, type of variables

      implicit none

c..   define type

      integer nn,i,j,k

c..   give value to parameters

      parameter(nn=5)

c..   define 2-dimensional arrayies

      real a(nn,nn),b(nn,nn),c(nn,nn)

      write(6,*)'This is to introduce FORTRAN'

c..   open files for input data and output data

      open (1, file='input.array',status='unknown',form='formatted')
      open (2, file='output.array',status='unknown',form='formatted')

c..   Note: Don't use channel 5 and 6, which are screen input and output

c..   input array a and b

      read(1,*)
      do i=1,nn
         read(1,*)(a(i,j),j=1,nn)
      end do

      read(1,*)
      do i=1,nn
         read(1,*)(b(i,j),j=1,nn)
      end do

c..   c=a*b

      do i=1,nn
         do j=1,nn
            c(i,j)=0.0
         end do
      end do

      do i=1,nn
         do j=1,nn
            do k=1,nn
               c(i,j)=c(i,j) + a(i,k)*b(k,j)
            end do
         end do
      end do

c..   print the c array on screen

      write(6,*)'c=a*b'

      do i=1,nn
         write(6,*)(c(i,j),j=1,nn)
      end do

c..   print the c array to file "output.array"

      write(2,*)'c=a*b'

      do i=1,nn
         write(2,*)(c(i,j),j=1,nn)
      end do

c..   Next, we want to do addition c=a+b
c..   We will do it through calling a subroutine

      call add(a,b,c,nn)

c..   print the c array on screen

      write(6,*)'c=a+b'

      do i=1,nn
         write(6,*)(c(i,j),j=1,nn)
      end do

c..   print the c array to file "output.array"

      write(2,*)'c=a+b'

      do i=1,nn
         write(2,*)(c(i,j),j=1,nn)
      end do

c..   ending the program

      stop
      end


      subroutine add(a,b,c,n)

c..   this subroutine adds array a and b output the results in c
c..   n is the dimension of hte array

c..   INPUT: a,b,n
c..   OUTPUT: c

      implicit none

      integer n,i,j
      real a(n,n), b(n,n), c(n,n)

      do i = 1,n
         do j= 1,n

            c(i,j) = a(i,j) + b(i,j)

         end do
      end do

      return
      end