Google
 

Trailing-Edge - PDP-10 Archives - decus_20tap2_198111 - decus/20-0026/dapch.ssp
There are 2 other files named dapch.ssp in the archive. Click here to see a list.
C                                                                       DAPC  10
C     ..................................................................DAPC  20
C                                                                       DAPC  30
C        SUBROUTINE DAPCH                                               DAPC  40
C                                                                       DAPC  50
C        PURPOSE                                                        DAPC  60
C           SET UP NORMAL EQUATIONS OF LEAST SQUARES FIT IN TERMS OF    DAPC  70
C           CHEBYSHEV POLYNOMIALS FOR A GIVEN DISCRETE FUNCTION         DAPC  80
C                                                                       DAPC  90
C        USAGE                                                          DAPC 100
C           CALL DAPCH(DATI,N,IP,XD,X0,WORK,IER)                        DAPC 110
C                                                                       DAPC 120
C        DESCRIPTION OF PARAMETERS                                      DAPC 130
C           DATI  - VECTOR OF DIMENSION 3*N (OR DIMENSION 2*N+1)        DAPC 140
C                   CONTAINING THE GIVEN ARGUMENTS, FOLLOWED BY THE     DAPC 150
C                   FUNCTION VALUES AND N (RESPECTIVELY 1) WEIGHT       DAPC 160
C                   VALUES. THE CONTENT OF VECTOR DATI REMAINS          DAPC 170
C                   UNCHANGED.                                          DAPC 180
C                   DATI MUST BE OF DOUBLE PRECISION                    DAPC 190
C           N     - NUMBER OF GIVEN POINTS                              DAPC 200
C           IP    - DIMENSION OF LEAST SQUARES FIT, I.E. NUMBER OF      DAPC 210
C                   CHEBYSHEV POLYNOMIALS USED AS FUNDAMENTAL FUNCTIONS DAPC 220
C                   IP SHOULD NOT EXCEED N                              DAPC 230
C           XD    - RESULTANT MULTIPLICATIVE CONSTANT FOR LINEAR        DAPC 240
C                   TRANSFORMATION OF ARGUMENT RANGE                    DAPC 250
C                   XD MUST BE DOUBLE PRECISION                         DAPC 260
C           X0    - RESULTANT ADDITIVE CONSTANT FOR LINEAR              DAPC 270
C                   TRANSFORMATION OF ARGUMENT RANGE                    DAPC 280
C                   X0 MUST BE DOUBLE PRECISION                         DAPC 290
C           WORK  - WORKING STORAGE OF DIMENSION (IP+1)*(IP+2)/2        DAPC 300
C                   ON RETURN WORK CONTAINS THE SYMMETRIC COEFFICIENT   DAPC 310
C                   MATRIX OF THE NORMAL EQUATIONS IN COMPRESSED FORM   DAPC 320
C                   FOLLOWED IMMEDIATELY BY RIGHT HAND SIDE             DAPC 330
C                   AND SQUARE SUM OF FUNCTION VALUES                   DAPC 340
C                   WORK MUST BE OF DOUBLE PRECISION                    DAPC 350
C           IER   - RESULTING ERROR PARAMETER                           DAPC 360
C                   IER =-1 MEANS FORMAL ERRORS IN DIMENSION            DAPC 370
C                   IER = 0 MEANS NO ERRORS                             DAPC 380
C                   IER = 1 MEANS COINCIDING ARGUMENTS                  DAPC 390
C                                                                       DAPC 400
C        REMARKS                                                        DAPC 410
C           NO WEIGHTS ARE USED IF THE VALUE OF DATI(2*N+1) IS          DAPC 420
C           NOT POSITIVE.                                               DAPC 430
C           EXECUTION OF SUBROUTINE DAPCH IS A PREPARATORY STEP FOR     DAPC 440
C           CALCULATION OF LEAST SQUARES FITS IN CHEBYSHEV POLYNOMIALS  DAPC 450
C           IT SHOULD BE FOLLOWED BY EXECUTION OF SUBROUTINE DAPFS      DAPC 460
C                                                                       DAPC 470
C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED                  DAPC 480
C           NONE                                                        DAPC 490
C                                                                       DAPC 500
C        METHOD                                                         DAPC 510
C           THE LEAST SQUARE FIT IS DETERMINED USING CHEBYSHEV          DAPC 520
C           POLYNOMIALS AS FUNDAMENTAL FUNCTION SYSTEM.                 DAPC 530
C           THE METHOD IS DISCUSSED IN THE ARTICLE                      DAPC 540
C           A.T.BERZTISS, LEAST SQUARES FITTING TO IRREGULARLY SPACED   DAPC 550
C           DATA, SIAM REVIEW, VOL.6, ISS.3, 1964, PP. 203-227.         DAPC 560
C                                                                       DAPC 570
C     ..................................................................DAPC 580
C                                                                       DAPC 590
      SUBROUTINE DAPCH(DATI,N,IP,XD,X0,WORK,IER)                        DAPC 600
C                                                                       DAPC 610
C                                                                       DAPC 620
C       DIMENSIONED DUMMY VARIABLES                                     DAPC 630
      DIMENSION DATI(1),WORK(1)                                         DAPC 640
      DOUBLE PRECISION DATI,WORK,XD,X0,XA,XE,XM,DF,T,SUM                DAPC 650
C                                                                       DAPC 660
C        CHECK FOR FORMAL ERRORS IN SPECIFIED DIMENSIONS                DAPC 670
      IF(N-1)19,20,1                                                    DAPC 680
    1 IF(IP)19,19,2                                                     DAPC 690
C                                                                       DAPC 700
C        SEARCH SMALLEST AND LARGEST ARGUMENT                           DAPC 710
    2 IF(IP-N)3,3,19                                                    DAPC 720
    3 XA=DATI(1)                                                        DAPC 730
      X0=XA                                                             DAPC 740
      XE=0.D0                                                           DAPC 750
      DO 7 I=1,N                                                        DAPC 760
      XM=DATI(I)                                                        DAPC 770
      IF(XA-XM)5,5,4                                                    DAPC 780
    4 XA=XM                                                             DAPC 790
    5 IF(X0-XM)6,7,7                                                    DAPC 800
    6 X0=XM                                                             DAPC 810
    7 CONTINUE                                                          DAPC 820
C                                                                       DAPC 830
C        INITIALIZE CALCULATION OF NORMAL EQUATIONS                     DAPC 840
      XD=X0-XA                                                          DAPC 850
      M=(IP*(IP+1))/2                                                   DAPC 860
      IEND=M+IP+1                                                       DAPC 870
      MT2=IP+IP                                                         DAPC 880
      MT2M=MT2-1                                                        DAPC 890
C                                                                       DAPC 900
C        SET WORKING STORAGE AND RIGHT HAND SIDE TO ZERO                DAPC 910
      DO 8 I=1,IP                                                       DAPC 920
      J=MT2-I                                                           DAPC 930
      WORK(J)=0.D0                                                      DAPC 940
      WORK(I)=0.D0                                                      DAPC 950
      K=M+I                                                             DAPC 960
    8 WORK(K)=0.D0                                                      DAPC 970
C                                                                       DAPC 980
C        CHECK FOR DEGENERATE ARGUMENT RANGE                            DAPC 990
      IF(XD)20,20,9                                                     DAPC1000
C                                                                       DAPC1010
C        CALCULATE CONSTANTS FOR REDUCTION OF ARGUMENTS                 DAPC1020
    9 X0=-(X0+XA)/XD                                                    DAPC1030
      XD=2.D0/XD                                                        DAPC1040
      SUM=0.D0                                                          DAPC1050
C                                                                       DAPC1060
C        START GREAT LOOP OVER ALL GIVEN POINTS                         DAPC1070
      DO 15 I=1,N                                                       DAPC1080
      T=DATI(I)*XD+X0                                                   DAPC1090
      J=I+N                                                             DAPC1100
      DF=DATI(J)                                                        DAPC1110
C                                                                       DAPC1120
C        CALCULATE AND STORE VALUES OF CHEBYSHEV POLYNOMIALS            DAPC1130
C        FOR ARGUMENT T                                                 DAPC1140
      XA=1.D0                                                           DAPC1150
      XM=T                                                              DAPC1160
      IF(DATI(2*N+1))11,11,10                                           DAPC1170
   10 J=J+N                                                             DAPC1180
      XA=DATI(J)                                                        DAPC1190
      XM=T*XA                                                           DAPC1200
   11 T=T+T                                                             DAPC1210
      SUM=SUM+DF*DF*XA                                                  DAPC1220
      DF=DF+DF                                                          DAPC1230
      J=1                                                               DAPC1240
   12 K=M+J                                                             DAPC1250
      WORK(K)=WORK(K)+DF*XA                                             DAPC1260
   13 WORK(J)=WORK(J)+XA                                                DAPC1270
      IF(J-MT2M)14,15,15                                                DAPC1280
   14 J=J+1                                                             DAPC1290
      XE=T*XM-XA                                                        DAPC1300
      XA=XM                                                             DAPC1310
      XM=XE                                                             DAPC1320
      IF(J-IP)12,12,13                                                  DAPC1330
   15 CONTINUE                                                          DAPC1340
      WORK(IEND)=SUM+SUM                                                DAPC1350
C                                                                       DAPC1360
C        CALCULATE MATRIX OF NORMAL EQUATIONS                           DAPC1370
      LL=M                                                              DAPC1380
      KK=MT2M                                                           DAPC1390
      JJ=1                                                              DAPC1400
      K=KK                                                              DAPC1410
      DO 18 J=1,M                                                       DAPC1420
      WORK(LL)=WORK(K)+WORK(JJ)                                         DAPC1430
      LL=LL-1                                                           DAPC1440
      IF(K-JJ)16,16,17                                                  DAPC1450
   16 KK=KK-2                                                           DAPC1460
      K=KK                                                              DAPC1470
      JJ=1                                                              DAPC1480
      GOTO 18                                                           DAPC1490
   17 JJ=JJ+1                                                           DAPC1500
      K=K-1                                                             DAPC1510
   18 CONTINUE                                                          DAPC1520
      IER=0                                                             DAPC1530
      RETURN                                                            DAPC1540
C                                                                       DAPC1550
C        ERROR RETURN IN CASE OF FORMAL ERRORS                          DAPC1560
   19 IER=-1                                                            DAPC1570
      RETURN                                                            DAPC1580
C                                                                       DAPC1590
C        ERROR RETURN IN CASE OF COINCIDING ARGUMENTS                   DAPC1600
   20 IER=1                                                             DAPC1610
      RETURN                                                            DAPC1620
      END                                                               DAPC1630