C #LNEQB V1A 17-APR-72. C LAST UPDATE: C SUBROUTINE LNEQB(A,B,N,M,L,DET) DIMENSION A(L,L), B(L,M) C C *PURPOSE. C TO SOLVE THE SET OF LINEAR EQUATIONS: C A * X = B C WHERE "A" IS A SQUARE MATRIX AND "B" AND "X" RECTANGULAR. C AS A SPECIAL CASE "X" AND "B" MAY BE VECTORS. C C *PARAMETERS: C C A - THE REAL ARRAY OF DIMENSION "L" X "L" CONTAINING THE C MATRIX OF THE L.H.S. OF THE EQUATIONS IN THE TOP LEFT C "N" X "N" SQUARE. THIS MATRIX IS DESTROYED ON RETURN. C B - THE REAL ARRAY OF DIMENSION "L" X "M" CONTAINING THE C "M" VECTORS OF THE R.H.S. OF THE EQUATIONS IN ITS TOP C LEFT "N" X "M" RECTANGLE. ON RETURN THE SOLUTION VECTORS C ARE IN PLACE OF THEIR CORRESPONDING R.H.S VECTORS. C N - THE LOGICAL SIZE OF THE MATRIX OF THE R.H.S. THE NUMBER C OF INDEPENDENT UNKNOWNS. C M - THE NUMBER OF RIGHT HAND SIDE VECTORS BEING SOLVED C SIMULTANEOUSLY. M=1 CORRESPONDS TO SOLVING A SINGLE C SET OF LINEAR EQUATIONS. C L - THE PHYSICAL SIZE OF THE ARRAY CONTAINING THE MATRICES. C DET - A REAL VALUE WHICH RETURNS THE DETERMINANT OF THE "N" X "N" C MATRIX "A". C C *METHOD. C THE METHOD IS THE SAME AS USED IN INVTB FOR INVERSION BUT IS C ADAPTED HERE TO THE DIRECT SOLUTION OF EQUATIONS AND IS MORE EFFICIENT C WHEN THE INVERSE OF "A" IS NOT REQUIRED EXPLICITLY. C C *ACCURACY. C SEE INVTB. C C *RESTRICTIONS. C C *ERROR CONDITIONS. C DET IS RETURNED ZERO IF THE MATRIX A WAS SINGULAR IN WHICH CASE C THE RESULTS ARE MEANINGLESS. C C *NON STANDARD ROUTINES CALLED. C C *TYPICAL TIMES. C THE TIME TAKEN IS ROUGHLY PROPORTIONAL TO (N**3.+M*N**2). C C *SIZE. 366 WORDS. C C *LANGUAGE. FORTRAN 4. C C *ORIGIN. M.R.MANNING. C C *COMMENTS. C IF ONE ACTUALLY WANTS TO SOLVE A SET OF MATRIX EQUATIONS THIS C ROUTINE SHOULD BE USED IN PREFERENCE TO INVTB AS IT TAKES LESS SPACE C AND IS FASTER. HOWEVER WHERE ONE REALLY REQUIRES THE INVERSE C MATRIX THE INVERSION ROUTINE SHOULD BE USED. THE SIMPLE CASE OF M=1 C WILL BE THE MOST FREQUENT ONE AND IN THIS CASE IT IS NOT NECESSARY TO C DECLARE B IN THE CALLING PROGRAM AS A TWO DIMENSIONAL ARRAY! C C C #END. C DET=1. DO 15 K=1,N KP=K+1 FAC=0. DO 3 J=K,N IF(ABS(A(J,K)).LE.ABS(FAC)) GO TO 3 JST=J FAC=A(J,K) 3 CONTINUE IF(FAC) 5,4,5 4 DET=0. RETURN 5 DET=DET*FAC IF(JST.EQ.K) GO TO 8 DO 6 I=1,M TEM=B(K,I) B(K,I)=B(JST,I) 6 B(JST,I)=TEM DO 7 I=K,N TEM=A(K,I) A(K,I)=A(JST,I) 7 A(JST,I)=TEM 8 DO 9 I=1,M 9 B(K,I)=B(K,I)/FAC A(K,K)=1. DO 10 I=KP,N 10 A(K,I)=A(K,I)/FAC DO 15 I=1,N IF(I.EQ.K) GO TO 15 FAC=A(I,K) A(I,K)=0. DO 11 J=KP,N 11 A(I,J)=A(I,J)-FAC*A(K,J) DO 12 J=1,M 12 B(I,J)=B(I,J)-FAC*B(K,J) 15 CONTINUE RETURN END