Penadés Martínez, José

SUBROUTINE CHEB(M, N, MDIM, NDIM, A, B, TOL, RELERR, X, RANK, CHE 10
* RESMAX, ITER, OCODE)
C THIS SUBROUTINE USES A MODIFICATION OF THE SIMPLEX METHOD
C OF LINEAR PROGRAMMING TO CALCULATE A CHEBYSHEV SOLUTION TO
C AN OVER-DETERMINED SYSTEM OF LINEAR EQUATIONS.
C DESCRIPTION OF PARAMETERS.
C M NUMBER OF EQUATIONS.
C N NUMBER OF UNKNOWNS (N MUST NOT EXCEED M).
C MDIMTHE NUMBER OF COLUMNS OF A, AT LEAST M+1.
C NDIM THE NUMBER OF ROWS OF A, AT LEAST N+3.
C A TWO DIMENSIONAL REAL ARRAY OF SIZE (NDIM,MDIM).
C ON ENTRY,THE TRANSPOSE OF THE MATRIX OF
C COEFFICIENTS OF THE OVER-DETERMINED SYSTEM MUST
C BE STORED IN THE FIRST M COLUMNS AND N ROWS OF A.
C THESE VALUES ARE DESTROYED BY THE SUBROUTINE.
C B ONE DIMENSIONALREAL ARRAY OF SIZE MDIM. ON ENTRY,
C B MUST CONTAIN THE RIGHT-HAND SIDES OF THE
C EQUATIONS IN ITS FIRST M LOCATIONS. ON EXIT, B
C CONTAINS THE RESIDUALS FOR THE EQUATIONS IN ITS
C FIRST M LOCATIONS (SEE DESCRIPTION).
C TOL A SMALL POSITIVE TOLERANCE. EMPIRICAL EVIDENCE
C SUGGESTS TOL=10**(-D+1) WHERE D REPRESENTS THE
C NUMBER OF DECIMAL DIGITSOF ACCURACY AVAILABLE
C (SEE DESCRIPTION).
C RELERR A REAL VARIABLE WHICH ON ENTRY MUST HAVE THE VALUE
C 0.0 IF A CHEBYSHEV SOLUTION IS REQUIRED. IF RELERR
C IS POSITIVE, THE SUBROUTINE CALCULATES AN
C APPROXIMATE SOLUTION WITH RELERR AS AN UPPER BOUND
C ON THE RELATIVE ERROR OF ITS LARGEST RESIDUAL (SEE
C DESCRIPTION-INEQUALITY (2)). ON EXIT, THEVALUE OF
C RELERR GIVES A SMALLER UPPER BOUND FOR THIS
C RELATIVE ERROR.
C X ONE DIMENSIONAL REAL ARRAY OF SIZE NDIM. ON EXIT,
C THIS ARRAY CONTAINS A SOLUTION TO THE PROBLEM IN
C ITS FIRST N LOCATIONS.
C RANK AN INTEGER WHICH GIVES ON EXIT THE RANK OF THE
C MATRIX OF COEFFICIENTS.
C RESMAX THE LARGEST RESIDUAL IN MAGNITUDE.
C ITER THE NUMBER OFSIMPLEX ITERATIONS PERFORMED.
C OCODE AN EXIT CODE WITH VALUES..
C 0 - OPTIMAL SOLUTION WHICH IS PROBABLY
C NON-UNIQUE (SEE DESCRIPTION).
C 1 - UNIQUE OPTIMAL SOLUTION.
C 2 - CALCULATIONS TERMINATED PREMATURELY
C DUE TO ROUNDING ERRORS.
C IF YOUR FORTRAN COMPILER PERMITS A SINGLE COLUMN OF A TWO
C DIMENSIONAL ARRAY TOBE PASSED TO A ONE DIMENSIONAL ARRAY
C THROUGH A SUBROUTINE CALL, CONSIDERABLE SAVINGS IN
C EXECUTION TIME MAY BE ACHIEVED THROUGH THE USE OF THE
C FOLLOWING SUBROUTINE WHICH OPERATES ON COLUMN VECTORS.
C SUBROUTINE COL (V1,V2,MLT,NOTROW,I1,NP2)
C THIS SUBROUTINE SUBTRACTS FROM THE VECTOR V1 A MULTIPLE OF
C THE VECTOR V2 STARTING AT THE I1*TH ELEMENT UP TO THE
C NP2*TH ELEMENT, EXCEPTFOR THE NOTROW*TH ELEMENT.
C REAL V1(NP2),V2(NP2),MLT
C DO 1 I=I1,NP2
C IF(I.EQ.NOTROW) GO TO 1
C V1(I)=V1(I)-MLT*V2(I)
C 1 CONTINUE
C RETURN
C END
C SEE COMMENTS FOLLOWING STATEMENT NUMBER 340 FOR
C INSTRUCTIONS ON THE IMPLEMENTATION OF THIS MODIFICATION.
DIMENSION A(NDIM,MDIM), B(MDIM), X(NDIM)
INTEGER PROW, PCOL, RANK, RANKP1, OCODE
C BIGMUST BE SET EQUAL TO ANY VERY LARGE REAL CONSTANT.
C ITS VALUE HERE IS APPROPRIATE FOR THE IBM/370/145.
DATA BIG /1.E+75/
C INITIALIZATION.
MP1 = M + 1
NP1 = N + 1
NP2 = N + 2
NP3 = N + 3
NP1MR = 1
RANK = N
RELTMP = RELERR
RELERR = 0.
DO 10 J=1,M
A(NP1,J) = 1.
A(NP2,J) = -B(J)
A(NP3,J) = N + J
10CONTINUE
A(NP1,MP1) = 0.
ITER = 0
OCODE = 1
DO 20 I=1,N
X(I) = 0.
A(I,MP1) = I
20 CONTINUE
C LEVEL 1.
LEV = 1
K = 0
30 K = K + 1
KP1 = K + 1
NP1MK = NP1 - K
MODE = 0
DO 40 J=K,M
B(J) = 1.
40 CONTINUE
C DETERMINE THE VECTOR TO ENTER THE BASIS.
50 D = -BIG
DO 60 J=K,M...