c=========================================================== c Implements matrix-matrix multiply c c c = a b c c where a, b and c are n x n (square) real*8 matrices. c=========================================================== subroutine dmmmult(a,b,c,n) implicit none integer n real*8 a(n,n), b(n,n), c(n,n) integer i, j, k do j = 1 , n do i = 1 , n c(i,j) = 0.0d0 do k = 1 , n c(i,j) = c(i,j) + a(i,k) * b(k,j) end do end do end do return end c=========================================================== c Writes a double precision matrix (two dimensional c array) to file 'fname'. If 'fname' is the c string '-', the matrix is written to standard input. c c This routine is modelled on 'dvto' previously c discussed in class: see ~phys410/f77/ex3/dvto.f c=========================================================== subroutine dmto(fname,a,d1,d2) c----------------------------------------------------------- c Arguments: c c fname: (I) File name c a: (I) Input matrix c d1: (I) First dimension of a c d2: (I) Second dimension of a c----------------------------------------------------------- implicit none integer indlnb, getu character*(*) fname integer d1, d2 real*8 a(d1,d2) integer ustdout parameter ( ustdout = 6 ) integer uto, rc c----------------------------------------------------------- c Parse fname: either "attach" 'uto' to stdout or c get a unit number using 'getu', and open the c file 'fname' for formatted I/O via 'uto' c----------------------------------------------------------- if( fname .eq. '-' ) then uto = ustdout else uto = getu() open(uto,file=fname(1:indlnb(fname)), & form='formatted',iostat=rc) if( rc .ne. 0 ) then write(0,*) 'dmto: Error opening ', & fname(1:indlnb(fname)) return end if end if c----------------------------------------------------------- c Write dimensions, then array elements c----------------------------------------------------------- write(uto,*,iostat=rc) d1, d2 if( rc .ne. 0 ) then write(0,*) 'dmto: Error writing dimensions' go to 500 end if write(uto,*,iostat=rc) a if( rc .ne. 0 ) then write(0,*) 'dmto: Error reading matrix' end if c----------------------------------------------------------- c Exit: Close file and return c----------------------------------------------------------- 500 continue if( uto .ne. ustdout ) then close(uto) end if return end c=========================================================== c Returns a double precision matrix (two dimensional c array) read from file 'fname'. If 'fname' is the c string '-', the matrix is read from standard input. c c The dimensions of the matrix must precede the matrix c elements themselves in the file. Specifically, the c file should have been created using the following c list-directed, formatted READ statement c (or equivalent): c c integer d1, d2 c real*8 a(d1,d2) c integer uout c c write(uout,*) d1, d2 c write(uout,*) a c c This routine is modelled on 'dvfrom' previously c discussed in class: see ~phys410/f77/ex3/dvfrom.f c c Note the use of helper routine 'dmfrom1' which c reads actual array values once bounds have been c extracted from file. c=========================================================== subroutine dmfrom(fname,a,d1,d2,asize) c----------------------------------------------------------- c Arguments: c c fname: (I) File name c a: (O) Return matrix c d1: (O) First dimension of a c d2: (O) Second dimension of a c asize: (I) Maximum size (d1 * d2) of a c----------------------------------------------------------- implicit none integer indlnb, getu character*(*) fname integer d1, d2, asize real*8 a(d1,d2) integer ustdin parameter ( ustdin = 5 ) integer ufrom, rc c----------------------------------------------------------- c Parse fname: either "attach" 'ufrom' to stdin or c get a unit number using 'getu', and open the c file 'fname' for formatted I/O via 'ufrom' c----------------------------------------------------------- if( fname .eq. '-' ) then ufrom = ustdin else ufrom = getu() open(ufrom,file=fname(1:indlnb(fname)), & form='formatted',iostat=rc,status='old') if( rc .ne. 0 ) then write(0,*) 'dmfrom: Error opening ', & fname(1:indlnb(fname)) return end if end if c----------------------------------------------------------- c Read dimensions and abort if there is insufficient c storage for the entire matrix. Note the 'go to' c to the 'exit block' since we've opened a file now c and should close it, even if there's an error. c Also, we set the dimensions to 0 for all error c conditions as a way of communicating failure to c the calling routine. c----------------------------------------------------------- read(ufrom,*,iostat=rc) d1, d2 if( rc .ne. 0 ) then write(0,*) 'dmfrom: Error reading dimensions' d1 = 0 d2 = 0 go to 500 end if if( (d1 * d2) .gt. asize ) then write(0,*) 'dmfrom: Insufficient storage' d1 = 0 d2 = 0 go to 500 end if c----------------------------------------------------------- c Now that dimensions have been determined call c helper routine to read values c----------------------------------------------------------- call dmfrom1(ufrom,a,d1,d2,rc) if( rc .ne. 0 ) then write(0,*) 'dmfrom: Error reading matrix' d1 = 0 d2 = 0 end if c----------------------------------------------------------- c Exit: Close file and return c----------------------------------------------------------- 500 continue if( ufrom .ne. ustdin ) then close(ufrom) end if return end c=========================================================== c Helper routine for dmfrom: Reads array values, returns c I/O status to calling routine via 'rc' c=========================================================== subroutine dmfrom1(ufrom,a,d1,d2,rc) implicit none integer d1, d2, ufrom, rc real*8 a(d1,d2) read(ufrom,*,iostat=rc) a return end