program gauss_main
       real a(3,3), b(3,1)
       integer n, np, m, mp, i
       parameter (n=3, np=3, m=1, mp=1)       
       open (unit=20, file='a.dat')
       do i=1,3
         read (20,*) a(i,1), a(i,2), a(i,3)
       enddo
       close(20)
       open (unit=20, file='b.dat')
       do i=1,3
         read (20,*) b(i,1)
       enddo
       close(20)
       
       call gaussj(a,n,np,b,m,mp)
      
       write(*,*) 'The inverse matrix is :' 
       do i=1,3
         write (*,*) a(i,1), a(i,2), a(i,3)
       enddo
       write (*,*)
       write (*,*) 'and the solution matrix is :'
       do i=1,3
         write (*,*) b(i,1)
       enddo

       end