1 subroutine lubksb(aa,n,np,indx,bb) 2c--------------------------------------------------------------------- 3c 4c LU back substitution routine. 5c 6c--------------------------------------------------------------------- 7 include "common.h" 8 9 dimension indx(n),aa(np,np),bb(n) 10 11 ii=0 12 do i=1,n 13 ll=indx(i) 14 sum=bb(ll) 15 bb(ll)=bb(i) 16 if (ii.ne.0)then 17 do j=ii,i-1 18 sum=sum-aa(i,j)*bb(j) 19 enddo 20 else if (sum.ne.0.0) then 21 ii=i 22 endif 23 bb(i)=sum 24 enddo 25 do i=n,1,-1 26 sum=bb(i) 27 do j=i+1,n 28 sum=sum-aa(i,j)*bb(j) 29 enddo 30 bb(i)=sum/aa(i,i) 31 enddo 32 33 return 34 end 35 36 37 38