Help understanding code!

[mfm17]

Hi, I am new with Fortran and have some understanding of what this code does but would anyone be able to describe what is going on inside each loop? Thanks!

SUBROUTINE gaussj(a,n,np,b,m,mp)
INTEGER m,mp,n,np,NMAX
REAL a(np,np),b(np,mp)
PARAMETER (NMAX=50)
INTEGER i,icol,irow,j,k,l,ll,indxc(NMAX),indxr(NMAX),ipiv(NMAX)
REAL big,dum,pivinv
do 11 j=1,n
ipiv(j)=0
11 continue
do 22 i=1,n
big=0.
do 13 j=1,n
if(ipiv(j).ne.1)then
do 12 k=1,n
if (ipiv(k).eq.0) then
if (abs(a(j,k)).ge.big)then
big=abs(a(j,k))
irow=j
icol=k
endif
else if (ipiv(k).gt.1) then
pause 'singular matrix in gaussj'
endif
12 continue
endif
13 continue
ipiv(icol)=ipiv(icol)+1
if (irow.ne.icol) then
do 14 l=1,n
dum=a(irow,l)
a(irow,l)=a(icol,l)
a(icol,l)=dum
14 continue
do 15 l=1,m
dum=b(irow,l)
b(irow,l)=b(icol,l)
b(icol,l)=dum
15 continue
endif
indxr(i)=irow
indxc(i)=icol
if (a(icol,icol).eq.0.) pause 'singular matrix in gaussj'
pivinv=1./a(icol,icol)
a(icol,icol)=1.
do 16 l=1,n
a(icol,l)=a(icol,l)*pivinv
16 continue
do 17 l=1,m
b(icol,l)=b(icol,l)*pivinv
17 continue
do 21 ll=1,n
if(ll.ne.icol)then
dum=a(ll,icol)
a(ll,icol)=0.
do 18 l=1,n
a(ll,l)=a(ll,l)-a(icol,l)*dum
18 continue
do 19 l=1,m
b(ll,l)=b(ll,l)-b(icol,l)*dum
19 continue
endif
21 continue
22 continue
do 24 l=n,1,-1
if(indxr(l).ne.indxc(l))then
do 23 k=1,n
dum=a(k,indxr(l))
a(k,indxr(l))=a(k,indxc(l))
a(k,indxc(l))=dum
23 continue
endif
24 continue
return
END

 

[xwb]

It is very old (possibly Fortran IV) code for possibly solving an equation using gauss jordan (the subroutine name is gaussj). Try this

1) find a simultaneous equation where you know the answer
2) create a program that calls this subroutine and prints the results
3) stick in the figures
4) check if the result is what you expect it to be
5) re-run it in a debugger and see what it does by yourself

That will help you understand the program

 

[gsal]

Oh, please, do not exagerate; I don't think this is Fortran IV...not with "then - endif". Fortran IV only had logical and numerical 'if'