C kyst2a.f C The authors of this software are Joseph B Kruskal, Forrest W Young, and C Judith B Seery. C Copyright (c) 1993 by AT&T. C Permission to use, copy, modify, and distribute this software for any C purpose without fee is hereby granted, provided that this entire C notice is included in all copies of any software which is or includes C a copy or modification of this software and in all copies of the C supporting documentation for such software. C THIS SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMP- C LIED WARRANTY. IN PARTICULAR, NEITHER THE AUTHORS NOR AT&T MAKE ANY C REPRESENTATION OR WARRANTY OF ANY KIND CONCERNING THE MERCHANTABILITY C OF THIS SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR PURPOSE. C This software comes from the SECOND MDS Package of AT&T Bell Laboratories C The manual of how to use this program is available in several different C formats: see mds/kyst2a_manual. C For explanation of the method this software implements, see C "Multidimensional Scaling" by Joseph B. Kruskal and Myron Wish, C a monograph published by Sage Publications, Beverly Hills, Calif. C in 1978 as series 07, item 011, in the Sage University Papers, and C "Multidimensional Scaling by Optimizing Goodness of Fit to a Nonmetric C Hypothesis" by Joseph B. Kruskal in Psychometrika in 1964, vol 29, C nos. 1 and 2, pages 1-27 and 115-129. C----+----@----+----@----+----@----+----@----+----@----+----@----+----@ CMAIN MAIN 1 C KYST-2A AUGUST, 1977 MAIN 2 C MAIN 3 C KRUSKAL-YOUNG-SHEPARD-TORGERSON MULTIDIMENSIONAL SCALING PROGRAM MAIN 4 C MAIN 5 C VERSION 1 JANUARY, 1973 MAIN 6 C VERSION 2 OCTOBER, 1976 MAIN 7 C VERSION 2A AUGUST, 1977 MAIN 8 C (THIS VERSION INCLUDES REVISED MAIN PROGRAM MENTIONED MAIN 9 C IN THE ERRATA SHEET OF JULY 1, 1977, AND INCLUDES THE MAIN 10 C CORRECTIONS FROM THE ERRATA SHEET OF JULY 31, 1977.) MAIN 11 C MAIN 12 C WRITTEN BY: MAIN 13 C MAIN 14 C DR. J. B. KRUSKAL MAIN 15 C BELL TELEPHONE LABORATORIES MAIN 16 C MURRAY HILL, N. J. MAIN 17 C MAIN 18 C DR. F. W. YOUNG MAIN 19 C PSYCHOMETRIC LABORATORY MAIN 20 C UNIVERSITY OF NORTH CAROLINA MAIN 21 C CHAPEL HILL, N. C. MAIN 22 C ASSISTED BY: MAIN 23 C JUDITH B. SEERY MAIN 24 C BELL TELEPHONE LABORATORIES MAIN 25 C MURRAY HILL, N. J. MAIN 26 C MAIN 27 C MAIN 28 C MAIN 29 C **********************************************************************MAIN 30 C MAIN 31 C MAIN 31 C GENERAL REMARKS. MAIN 32 C MAIN 33 C MAIN 34 C KYST INCLUDES THE FOLLOWING SUBROUTINES: MAIN 35 C MAIN 36 C BLDA EQSOLV PLOT SERCH MAIN 37 C BSEC1 FITM REGR SGEV MAIN 38 C CCACT FITMS RM1POW SORT MAIN 39 C CONFIG FITP RPOWER ST03 MAIN 40 C DATAPR INICON RROOT WTRAN MAIN 41 C DFLT NEWSTP RUNIFV XITEM MAIN 42 C DTRAN NRMLZ SBK NVIT1 MAIN 43 C MAIN 44 C ALL ROUTINES ARE WRITTEN IN FORTRAN IV, AND PUNCHED IN THE IBMEL 029 MAIN 45 C CHARACTER SET MAIN 46 C MAIN 47 C NO USE IS MADE OF SPECIAL OR NON-STANDARD SOFTWARE. MAIN 48 C MAIN 49 C ALL INPUT AND OUTPUT IS ONTO FORTRAN LOGICAL UNITS WITH NUMBERS MAIN 50 C CONTROLLED BY THESE VARIABLES: LREAD,LPRINT,LPUNCH,LSCRAT. MAIN 51 C UNIT NUMBERS 5,6,7,8 HAVE BEEN USED RESPECTIVELY. TO CHANGE THESEMAIN 52 C ASSIGNMENTS, CHANGE THE VALUES FOR THE FOUR VARIABLES SET IN THE MAIN 53 C BLOCK DATA SUBPROGRAM. MAIN 54 C MAIN 55 C THE SCRATCH UNIT IS USED IN A MINOR WAY BY CCACT ONLY. MANY MAIN 56 C INSTALLATIONS WILL HAVE ALTERNATE METHODS OF DOING THE SAME THING.MAIN 57 C MAIN 58 C KYST USES THREE MACHINE DEPENDENT CONSTANTS: PRECSN (RELATIVE MACHINE MAIN 59 C PRECISION), XMAG (LEAST POSITIVE MACHINE NUMBER), XMAXN (MINIMUM MAIN 60 C OF LARGEST POSITIVE MACHINE NUMBER AND MINUS LARGEST NEGATIVE MAIN 61 C MACHINE NUMBER), FLPRIN (FLOATING LINES PRINTED PER INCH OF PAPER)MAIN 62 C , CHPRIN (CHARACTERS PRINTED PER INCH ON A LINE). VALUES ARE MAIN 63 C 1.5E-8, 1.0E-38, 1.0E38, 6.0, AND 10.0 MAIN 64 C RESPECTIVELY. TO CHANGE THESE ASSIGNMENTS, CHANGE THE VALUES MAIN 65 C FOR THESE VARIABLES SET IN THE BLOCK DATA SUBPROGRAM. MAIN 66 C MAIN 67 C **********************************************************************MAIN 68 C MAIN 69 C MAIN 70 C MAIN PROGRAM FOR KYST JANUARY,1973 MAIN 71 C WRITTEN BY J.KRUSKAL MAIN 72 C MODIFIED FOR KYST BY J.KRUSKAL AND J.SEERY JANUARY,1973 MAIN 73 DIMENSION DATA(1800), IJ(1800), DIST(1800), DHAT(1800) MAIN 74 DIMENSION WW(1800),LDIST(1800) MAIN 75 DIMENSION X(100,6), GR(100,6), GL(100,6) MAIN 76 REAL PMAT(1800,2),RVEC(1800) MAIN 77 INTEGER GRNO(101), NOGRPS, LSPLIT, SPLITH MAIN 78 INTEGER GRSDIS(100), SDSWIT,SDSWT1 MAIN 79 INTEGER FMAT(80),PTID(100),XNUM(6),ONE,TWO,THREE,FOUR,FIVE,SIX MAIN 80 INTEGER AAA,BEE,CEE,DEE,FFF,AIE,EEE MAIN 81 INTEGER TITLE(80),CTITLE(80),LBLOCK(1800),ITEM(101) MAIN 82 REAL GRSTRS(100),PCOEFF(5),DUMMY(25) MAIN 83 DIMENSION STPL(10),DIMSV(10) MAIN 84 DIMENSION X2(100),Y2(100),PTID2(100),DUMEST(100) MAIN 85 INTEGER RTYPE MAIN 86 C MAIN 87 EXTERNAL WTRAN , DTRAN MAIN 88 C MAIN 89 COMMON /ACCUR/ PRECSN,XMAG,XMAXN,FLPRIN,CHPRIN MAIN 90 COMMON /KYST1/ DATA,WW,X,IJ MAIN 91 COMMON /KYST2/ GR,GL,RVEC,PCOEFF,DUMMY,PMAT MAIN 92 COMMON /MDSCL1/ LREAD, LPRINT, LPUNCH, LSCRAT MAIN 93 COMMON /MDSCL2/ MAIN 94 1 LDIMX, LDIMN, LDIMD, CUTOFF, STRMIN MAIN 95 2 , SFGRMN, COSAVW, ACSAVW, IFIRST, MATSWP MAIN 96 3 , SDSWIT, LCSWIT, LFITSW, RCOM, NOIT MAIN 97 4 , SRATST, LSCH, LPUNSW, LSPL, LRANDC MAIN 98 5 , LDAPRT, LDIPRT, LREG, LHIPRT, NUDASW MAIN 99 6 , LPOLSP, LCONSW, LNFIXZ, DCON1, DCON2 MAIN 100 7 , DCON3, DCON4, DCON5, WCON1, WCON2 MAIN 101 8 , WCON3, WCON4, WCON5, NOITIN, LPLCON MAIN 102 9 , LPLSCT, LCOORS, LDATAF,SPREAD MAIN 103 COMMON /METRIC/ RM1,RECR,R ,RTYPE MAIN 104 COMMON /PLTCHR/ PTID,ITEM MAIN 105 C MAIN 106 EQUIVALENCE (DIST(1),LDIST(1)) , (LBLOCK(1),RVEC(1)) MAIN 107 EQUIVALENCE (PMAT(1,1),DIST(1)), (PMAT(1,2),DHAT(1)) MAIN 108 EQUIVALENCE (XNUM(1),ONE),(XNUM(2),TWO),(XNUM(3),THREE), MAIN 109 . (XNUM(4),FOUR),(XNUM(5),FIVE),(XNUM(6),SIX), MAIN 110 1 (GR(1,1),X2(1)),(GR(1,2),Y2(1)),(GR(1,3),PTID2(1)) MAIN 111 C MAIN 112 DATA AAA,BEE,CEE,DEE,EEE,FFF /1HA,1HB,1HC,1HD,1HE,1HF/ MAIN 113 DATA ONE,TWO,THREE,FOUR,FIVE,SIX /1H1,1H2,1H3,1H4,1H5,1H6/ MAIN 114 DATA KPACK,IJPACK,NPREVZ /10000,100,0/ , AIE/1H-/ MAIN 115 DATA IVRSN /0/ MAIN 116 C MAIN 117 C INITIALIZE PARAMETERS MAIN 118 C MAIN 119 100 CONTINUE MAIN 120 C MAIN 121 C ALPHABETICAL ORDER MAIN 122 ACSAVW=0.66 MAIN 123 COSAVW=0.66 MAIN 124 CUTOFF=-1.23E+20 MAIN 125 DCON1=0.0 MAIN 126 DCON2=1.0 MAIN 127 DCON3=1.0 MAIN 128 DCON4=0.0 MAIN 129 DCON5=1.0 MAIN 130 GRNO(1) = 1 MAIN 131 IFIRST=1 MAIN 132 ISTC=0 MAIN 133 LCONSW=4 MAIN 134 LCOORS=2 MAIN 135 LCSWIT=1 MAIN 136 LDAPRT = 1 MAIN 137 LDATAF=LREAD MAIN 138 LDIMD=1 MAIN 139 LDIMN=2 MAIN 140 LDIMX=2 MAIN 141 LDIPRT = 1 MAIN 142 LFITSW=1 MAIN 143 LHIPRT = 2 MAIN 144 LNFIXZ = 0 MAIN 145 LPOLSP = 100 MAIN 146 LPLCON=2 MAIN 147 LPLSCT=3 MAIN 148 LPUNSW = 2 MAIN 149 LRANDC=-101 MAIN 150 LREG=0 MAIN 151 LSCH=1 MAIN 152 LSPL=302 MAIN 153 MATSWP=101 MAIN 154 NOIT=50 MAIN 155 NOITIN=1 MAIN 156 NUDASW=1 MAIN 157 R=2.0 MAIN 158 RCOM=-10.0 MAIN 159 SDSIN=1.0 MAIN 160 SDSWIT = 10 MAIN 161 SFGRMN=0.0 MAIN 162 SRATST=0.999 MAIN 163 SPREAD=10.0 MAIN 164 STRESS=1.0 MAIN 165 STRMIN=0.01 MAIN 166 WCON1=0.0 MAIN 167 WCON2=1.0 MAIN 168 WCON3=1.0 MAIN 169 WCON4=0.0 MAIN 170 WCON5=1.0 MAIN 171 WRITE(LPRINT, 17) MAIN 172 17 FORMAT(1H1) MAIN 173 IF(IVRSN.EQ.0) WRITE(LPRINT,9) MAIN 174 9 FORMAT(1H0,25X,43HOUTPUT FROM KYST2A--VERSION OF AUGUST, 1977) MAIN 175 IVRSN=1 MAIN 176 C MAIN 177 C CONTROL CARD READ MAIN 178 C MAIN 179 1000 LCSWIT=1 MAIN 180 CALL CCACT MAIN 181 GO TO (1000, 1100, 1165, 1182, 1200, 2000, 9000, 1190 ), LCSWIT MAIN 182 C MAIN 183 C DATA READ MAIN 184 C MAIN 185 1100 CONTINUE MAIN 186 NREPL1=1 MAIN 187 NREPL3 = 1 MAIN 188 LSPLIT=LSPL/100 MAIN 189 SPLITH = MOD(LSPL,100) MAIN 190 MATSW = MOD(MATSWP,100) MAIN 191 LBLKDS = MATSWP/100 MAIN 192 LPOSW = MOD(LPOLSP,100) MAIN 193 LFITRM = LPOLSP/100 MAIN 194 IF(LREG.EQ.0) GO TO 1104 MAIN 195 SDSWIT=LREG MAIN 196 IF(LREG.LT.10) SDSWIT = SDSWIT+LPOSW MAIN 197 1104 CONTINUE MAIN 198 SDSWIT = SDSWIT+100*LFITRM MAIN 199 C MAIN 200 IF(NUDASW.NE.1) MLASTD=MM MAIN 201 IF(NUDASW.NE.1) GO TO 1106 MAIN 202 NUDASW=2 MAIN 203 MLASTD=0 MAIN 204 MMIN=0 MAIN 205 NOGRPS=0 MAIN 206 1106 M=MLASTD MAIN 207 MA = M+1 MAIN 208 MSERA=MA MAIN 209 C MAIN 210 READ(LREAD,10) TITLE MAIN 211 WRITE(LPRINT,11) TITLE MAIN 212 IF(MATSW.GE.4) GO TO 1102 MAIN 213 1101 READ (LREAD,12) NPART,NREPL2,NREPL4 MAIN 214 WRITE (LPRINT,13) NPART,NREPL2,NREPL4 MAIN 215 NROWS = NPART MAIN 216 NCOLS = NPART MAIN 217 GO TO 1103 MAIN 218 1102 READ (LREAD,12) NROWS,NCOLS,NREPL2,NREPL4 MAIN 219 WRITE (LPRINT,13) NROWS,NCOLS,NREPL2,NREPL4 MAIN 220 NPART = NROWS+NCOLS MAIN 221 C MAIN 222 1103 N=NPART MAIN 223 IF(NREPL4.NE.0) NREPL3 = NREPL4 MAIN 224 IF(NREPL2.NE.0) NREPL1 = NREPL2 MAIN 225 IF(LBLKDS.EQ.2) N = NPART*NREPL3 MAIN 226 IF(NREPL2.EQ.0 .OR. NREPL4.EQ.0 ) WRITE (LPRINT,16) MAIN 227 1105 READ(LREAD,10) FMAT MAIN 228 WRITE(LPRINT,11) FMAT MAIN 229 C MAIN 230 DO 1162 NREPL=1,NREPL3 MAIN 231 IJBLKD=0 MAIN 232 IF(LBLKDS.EQ.2) IJBLKD =(NREPL-1)*NPART MAIN 233 MLASTG=M MAIN 234 IA = 1 MAIN 235 IB = NROWS MAIN 236 IF(MATSW.EQ.2) IA = IFIRST MAIN 237 IF(MATSW.EQ.3) IB = NROWS-(IFIRST-1) MAIN 238 C MAIN 239 DO 1160 I = IA,IB MAIN 240 MLASTR=M MAIN 241 LTEMP = NCOLS MAIN 242 IF(MATSW.EQ.2) LTEMP = I-IFIRST+1 MAIN 243 IF(MATSW.EQ.3) LTEMP = NCOLS-I-IFIRST+2 MAIN 244 MB = MA + LTEMP * NREPL1 - 1 MAIN 245 MSERB=MSERA+LTEMP*NREPL1-1 MAIN 246 ITRUE = I+IJBLKD MAIN 247 IF(MATSW.EQ.4) ITRUE = ITRUE+NCOLS MAIN 248 C MAIN 249 1130 READ (LDATAF,FMAT) (DATA(MP),MP=MA,MB) MAIN 250 C MAIN 251 DO 1150 MP=MA,MB MAIN 252 IF( DATA(MP)-CUTOFF ) 1150, 1150, 1140 MAIN 253 1140 CONTINUE MAIN 254 M = M + 1 MAIN 255 DATA(M)=DATA(MP) MAIN 256 WW(M) = 1.0 MAIN 257 J=((MP-MA)/NREPL1)+1 MAIN 258 IF(MATSW.EQ.3) J = J + (I-1) + (IFIRST-1) MAIN 259 J = J+IJBLKD MAIN 260 IF(MATSW.EQ.5) J = J + NROWS MAIN 261 IJ(M)=IJPACK*(ITRUE-1)+J-1 MAIN 262 LDIST(M) = MSERA+MP-MA MAIN 263 1150 CONTINUE MAIN 264 C MAIN 265 IF(LSPLIT.EQ.1 .AND. M.GT.MLASTR) GO TO 1155 MAIN 266 GO TO 1158 MAIN 267 1155 NOGRPS=NOGRPS+1 MAIN 268 GRNO(NOGRPS+1) = M+1 MAIN 269 GRSDIS(NOGRPS)=SDSWIT MAIN 270 1158 MA = M + 1 MAIN 271 1160 MSERA=MSERB+1 MAIN 272 C MAIN 273 IF(LSPLIT.EQ.2 .AND. M.GT.MLASTG) GO TO 1161 MAIN 274 GO TO 1162 MAIN 275 1161 NOGRPS=NOGRPS+1 MAIN 276 GRNO(NOGRPS+1) = M+1 MAIN 277 GRSDIS(NOGRPS)=SDSWIT MAIN 278 1162 CONTINUE MAIN 279 C MAIN 280 IF(LSPLIT.EQ.3 .AND. M.GT.MLASTD) GO TO 1163 MAIN 281 IF(LSPLIT.EQ.4 .AND. SPLITH.EQ.2 .AND. M.GT.MLASTD) GO TO 1163 MAIN 282 GO TO 1164 MAIN 283 1163 NOGRPS=NOGRPS+1 MAIN 284 GRSDIS(NOGRPS)=SDSWIT MAIN 285 C MAIN 286 1164 CONTINUE MAIN 287 GRNO(NOGRPS+1) = M+1 MAIN 288 IF(LSPLIT.LE.3) SPLITH=2 MAIN 289 IF(LSPLIT.EQ.4) SPLITH=1 MAIN 290 IF((NOGRPS.EQ.1).OR.(MMIN.EQ.0)) MMIN=M MAIN 291 C MAIN 292 LSPL = 100 * LSPLIT + SPLITH MAIN 293 C MAIN 294 MM = M MAIN 295 LDATAF=LREAD MAIN 296 GO TO 1000 MAIN 297 C MAIN 298 C WEIGHTS READ MAIN 299 C MAIN 300 1165 CONTINUE MAIN 301 M=MLASTD MAIN 302 MA=M+1 MAIN 303 MSERA=MA MAIN 304 DO 1181 NREPL=1,NREPL3 MAIN 305 IA = 1 MAIN 306 IB = NROWS MAIN 307 IF(MATSW.EQ.2) IA = IFIRST MAIN 308 IF(MATSW.EQ.3) IB = NROWS-(IFIRST-1) MAIN 309 DO 1180 I = IA,IB MAIN 310 LTEMP = NCOLS MAIN 311 IF(MATSW.EQ.2) LTEMP = I-IFIRST+1 MAIN 312 IF(MATSW.EQ.3) LTEMP = NCOLS-I-IFIRST+2 MAIN 313 MB = MA + LTEMP*NREPL1-1 MAIN 314 MSERB=MSERA+LTEMP*NREPL1-1 MAIN 315 1172 READ (LDATAF,FMAT) (WW(MP),MP=MA,MB) MAIN 316 DO 1177 MP=MA,MB MAIN 317 IF(LDIST(M+1).GT.(MSERA+MP-MA)) GO TO 1177 MAIN 318 M = M + 1 MAIN 319 WW(M)=WW(MP) MAIN 320 1177 CONTINUE MAIN 321 MA = M + 1 MAIN 322 1180 MSERA=MSERB+1 MAIN 323 1181 CONTINUE MAIN 324 LDATAF = LREAD MAIN 325 GO TO 1000 MAIN 326 C MAIN 327 C WEIGHT FORMATION BY WFUNCTION MAIN 328 C MAIN 329 1182 CONTINUE MAIN 330 WRITE(LPRINT,18) WCON1,WCON2,WCON3,WCON4,WCON5 MAIN 331 18 FORMAT(1H0,10X,69H** WEIGHTS FORMED ACCORDING TO RULE, WW(I,J)=S+TMAIN 332 .*((A+B*DATA(I,J))**P)/ 15X, 9H WHERE A=,F15.7,5X,2HB=,F15.7,5X, MAIN 333 .2HP=,F15.7,5X,2HS=,F15.7,5X,2HT=,F15.7) MAIN 334 MA=MLASTD+1 MAIN 335 DO 1185 M=MA,MM MAIN 336 TEMP1=DATA(M) MAIN 337 WW(M) = WTRAN(TEMP1) MAIN 338 1185 CONTINUE MAIN 339 GO TO 1000 MAIN 340 C MAIN 341 C DATA TRANSFORMATION BY DFUNCTION MAIN 342 C MAIN 343 1190 CONTINUE MAIN 344 WRITE(LPRINT,19) DCON1,DCON2,DCON3,DCON4,DCON5 MAIN 345 19 FORMAT(1H0,10X,73H** DATA TRANSFORMED ACCORDING TO RULE, DATA(I,J)MAIN 346 .=S+T*((A=B*DATA(I,J))**P) / 15X, 9H WHERE A=,F15.7,5X,2HB=, MAIN 347 .F15.7,5X,2HP=,F15.7,5X,2HS=,F15.7,5X,2HT=,F15.7) MAIN 348 MA=MLASTD+1 MAIN 349 DO 1195 M=MA,MM MAIN 350 DATA(M)=DTRAN(DATA(M)) MAIN 351 1195 CONTINUE MAIN 352 GO TO 1000 MAIN 353 C MAIN 354 C CONFIGURATION READ MAIN 355 C MAIN 356 1200 LCONSW=2 MAIN 357 READ(LREAD,10) CTITLE MAIN 358 WRITE(LPRINT,11)CTITLE MAIN 359 READ (LREAD, 12) NCON, LDIMCO MAIN 360 WRITE (LPRINT, 13) NCON, LDIMCO MAIN 361 READ(LREAD,10) FMAT MAIN 362 WRITE(LPRINT,11) FMAT MAIN 363 DO 1210 I=1,NCON MAIN 364 READ (LDATAF,FMAT) (X(I,L),L=1,LDIMCO) MAIN 365 WRITE(LPRINT,4008) I,(X(I,L),L=1,LDIMCO) MAIN 366 4008 FORMAT(1X,I2,10F7.3) MAIN 367 1210 CONTINUE MAIN 368 LDATAF = LREAD MAIN 369 GO TO 1000 MAIN 370 C MAIN 371 C SOME INPUT FORMATS MAIN 372 C MAIN 373 10 FORMAT (80A1) MAIN 374 11 FORMAT (1H0,80A1) MAIN 375 12 FORMAT(24I3) MAIN 376 13 FORMAT(1X,24I3) MAIN 377 15 FORMAT(20F4.0) MAIN 378 16 FORMAT(118H THE NUMBER OF INTERNAL REPLICATES, OR THE NUMBER OF EMAIN 379 1XTERNAL REPLICATES HAS NOT BEEN SPECIFIED. IT IS TAKEN TO BE 1 . )MAIN 380 C MAIN 381 C COMPUTATION *************************************MAIN 382 2000 CONTINUE MAIN 383 C MAIN 384 FN=FLOAT (N) MAIN 385 SQRTN=SQRT (FN) MAIN 386 FNGRPS = FLOAT(NOGRPS) MAIN 387 LDIM=LDIMX MAIN 388 IF(RCOM.NE.(-10.0)) R=RCOM MAIN 389 RR=R MAIN 390 RTYPE=3 MAIN 391 IF(R.EQ.1.0) RTYPE=1 MAIN 392 IF(R .EQ.2.0) RTYPE=2 MAIN 393 RM1=R-1.0 MAIN 394 RECR=1.0/R MAIN 395 IF(LDAPRT.EQ.2) MAIN 396 1 CALL DATAPR(GRNO,MM,NOGRPS,1) MAIN 397 WRITE (LPRINT,23) MAIN 398 C MAIN 399 C FINISH STARTING CONFIGURATION. LCONSW = 2,3,4. MAIN 400 C 4=USE INICON, 3=SAVED, 2=WAS READ IN. MAIN 401 C TO FILL IN ADDITIONAL POINTS, IF NEEDED: MAIN 402 C LRANDC<0, USE ARBITRARY, ELSE USE RANDOM. MAIN 403 C MAIN 404 2100 CONTINUE MAIN 405 WRITE(LPRINT,11) TITLE MAIN 406 IF((LCONSW.EQ.4) .AND. (LRANDC.EQ.(-101)))GO TO 2140 MAIN 407 ISTCON = 1 MAIN 408 IF( LCONSW .EQ. 2 ) ISTCON = NCON + 1 MAIN 409 IF( LCONSW .EQ. 3 ) ISTCON = NPREVZ + 1 MAIN 410 IF( ISTCON .GT. N ) GO TO 2200 MAIN 411 IF(LRANDC.LE.0) GO TO 2110 MAIN 412 DO 2105 K=1,LRANDC MAIN 413 2105 TEMP=RUNIFV(1.0) MAIN 414 2110 TEMP1=0.01 MAIN 415 DO 2130 I=ISTCON,N MAIN 416 DO 2120 L=1,LDIMX MAIN 417 2120 X(I,L)=0.0 MAIN 418 K= MOD (I-1,LDIM) +1 MAIN 419 X(I,K)=TEMP1 MAIN 420 IF(LRANDC.LT.0) GO TO 2130 MAIN 421 DO 2125 L=1,LDIMX MAIN 422 TEMP = RUNIFV(1.0) MAIN 423 2125 X(I,L) = ALOG( TEMP/(1.0-TEMP)) MAIN 424 2130 TEMP1=TEMP1+0.01 MAIN 425 GO TO 2200 MAIN 426 2140 SDSWIT=MOD(GRSDIS(1),100) MAIN 427 IF(SDSWIT.EQ.11) SDSIN=-1.0 MAIN 428 IF(N.GT.60) NOITIN=0 MAIN 429 CALL INICON(N,LDIM,MMIN,NOITIN,LHIPRT,SDSIN,LSCH,SPREAD) MAIN 430 C MAIN 431 C SORT DATA AND IJ AND WW. MAIN 432 C ALSO RECORD BLOCKS OF EQUAL DATA VALUES. MAIN 433 C MAIN 434 2200 CONTINUE MAIN 435 NPREVZ = N MAIN 436 DO 2250 NG = 1,NOGRPS MAIN 437 MX = GRNO(NG) MAIN 438 MY = GRNO(NG+1) - 1 MAIN 439 MZ = MY - MX + 1 MAIN 440 SDSWIT = GRSDIS(NG) MAIN 441 SDSWIT = MOD(SDSWIT,100) MAIN 442 IF(SDSWIT.EQ.11)SDSWIT = -10 MAIN 443 CALL SORT( DATA(MX),MZ,IJ(MX),WW(MX),DUMMY,2,SDSWIT) MAIN 444 C MAIN 445 C 'SORT' WILL SORT THE MM ELEMNTS OF 'DATA' IN ALGEBRAIC MAIN 446 C ORDER, ASCENDING OR DESCENDING ACCORDING TO WHETHER MAIN 447 C SDSWIT IS + OR -. MAIN 448 C AT THE SAME TIME, THE ELEMENTS IN 'IJ' AND IN 'WW' WILL BMAIN 449 C REARRANGED IN EXACTLY THE SAME ORDER. THUS THE MAIN 450 C CORRESPONDENCE BETWEEN THE ELEMENTS OF 'DATA' AND 'IJ' MAIN 451 C AND 'WW' IS PRESERVED. MAIN 452 C MAIN 453 DO 2240 MB = MX,MY MAIN 454 IF((DATA(MB+1).EQ.DATA(MB)).AND.(MB.NE.MY)) GO TO 2240 MAIN 455 IJ(MB)=MOD(IJ(MB),KPACK) MAIN 456 IJ(MB)=IJ(MB)+KPACK MAIN 457 2240 CONTINUE MAIN 458 2250 CONTINUE MAIN 459 C MAIN 460 C START COMPUTATION IN CURRENT DIMENSION MAIN 461 C MAIN 462 2300 FLDIM=FLOAT (LDIM) MAIN 463 ITNO=0 MAIN 464 COSAV=0.0 MAIN 465 SRATAV=0.8 MAIN 466 ACSAV=0.0 MAIN 467 STEP = 0.0 MAIN 468 NBAKUP = 0 MAIN 469 C MAIN 470 C INITIALIZE GRADIENT MAIN 471 C MAIN 472 2400 TEMP1=SQRT (1.0/FLDIM) MAIN 473 DO 2410 I=1,N MAIN 474 DO 2410 L=1,LDIM MAIN 475 2410 GR(I,L)=TEMP1 MAIN 476 SFGR=SQRTN MAIN 477 C MAIN 478 C PRINT HEADING FOR INFORMATION ABOUT SCALING IN CURRENT DIMENSION MAIN 479 C MAIN 480 2500 WRITE (LPRINT, 20) N, MM, NOGRPS, LDIM MAIN 481 IF(LHIPRT.EQ.2) WRITE (LPRINT,21) MAIN 482 WRITE (LPRINT, 22) MAIN 483 20 FORMAT(28H0HISTORY OF COMPUTATION. N= , I4, MAIN 484 1 15H. THERE ARE , I6, MAIN 485 2 26H DATA VALUES, SPLIT INTO , I4, MAIN 486 3 27H LISTS. DIMENSION = , I4 ) MAIN 487 21 FORMAT(52H0ITERATION STRESS SRAT SRATAV CAGRGL COSAV ACSAV, MAIN 488 1 16H SFGR STEP ) MAIN 489 22 FORMAT(1X) MAIN 490 23 FORMAT(1H1) MAIN 491 C MAIN 492 C START CURRENT ITERATION *************************************MAIN 493 C MAIN 494 C NORMALIZE CONFIGURATION. MOVE AND CLEAR GRADIENT. MAIN 495 C MAIN 496 3000 TEMP1=1.0 MAIN 497 IF(LCOORS.EQ.0) GO TO 3035 MAIN 498 TEMP1=0.0 MAIN 499 DO 3030 L=1,LDIM MAIN 500 TEMP2=0.0 MAIN 501 DO 3010 I=1,N MAIN 502 3010 TEMP2=TEMP2+X(I,L) MAIN 503 TEMP2=TEMP2/FN MAIN 504 DO 3020 I=1,N MAIN 505 X(I,L)=X(I,L)-TEMP2 MAIN 506 3020 TEMP1=TEMP1+X(I,L)**2 MAIN 507 3030 CONTINUE MAIN 508 TEMP1=SQRT (FN/TEMP1) MAIN 509 3035 DO 3040 L=1,LDIM MAIN 510 DO 3040 I=1,N MAIN 511 X(I,L)=TEMP1*X(I,L) MAIN 512 GL(I,L)=TEMP1*GR(I,L) MAIN 513 3040 GR(I,L)=0.0 MAIN 514 SFGL=TEMP1*SFGR MAIN 515 C MAIN 516 STBAMU = TEMP1 MAIN 517 C LOOP THROUGH THE DATA GROUPS MAIN 518 C MAIN 519 STRLST = STRESS MAIN 520 STRESS = 0.0 MAIN 521 DO 3340 NG = 1,NOGRPS MAIN 522 MX = GRNO(NG) MAIN 523 MY = GRNO(NG+1) - 1 MAIN 524 MZ = MY - MX + 1 MAIN 525 SDSWIT = GRSDIS(NG) MAIN 526 LFITRM = SDSWIT/100 MAIN 527 SDSWIT = MOD(SDSWIT,100) MAIN 528 IF(SDSWIT.EQ.11)SDSWIT = -10 MAIN 529 C MAIN 530 C COMPUTE DISTANCES AND FIND BEST-FITTING MONOTONE PSEUDO-DISTANCESMAIN 531 C MAIN 532 SUMW = 0.0 MAIN 533 DBAR = 0.0 MAIN 534 DO 3120 M=MX,MY MAIN 535 LTEMP1=MOD(IJ(M),KPACK) MAIN 536 I=LTEMP1/IJPACK+1 MAIN 537 J=MOD(LTEMP1,IJPACK)+1 MAIN 538 TEMP1=0.0 MAIN 539 DO 3110 L=1,LDIM MAIN 540 3110 TEMP1=TEMP1+RPOWER (X(I,L)-X(J,L)) MAIN 541 DIST(M)=RROOT (TEMP1) MAIN 542 DBAR=DBAR+DIST(M)*WW(M) MAIN 543 SUMW = SUMW + WW(M) MAIN 544 3120 CONTINUE MAIN 545 DBAR=DBAR/SUMW MAIN 546 C DBAS IS EITHER 0 OR DBAR ACCORDING TO WHETHER MAIN 547 C STRESS FORMULA 1 OR 2 IS BEING USED. MAIN 548 DBAS = 0.0 MAIN 549 IF(LSCH.EQ.2) DBAS = DBAR MAIN 550 IF(IABS(SDSWIT).GE.10) MAIN 551 1 CALL FITM(MX,MY,LFITSW) MAIN 552 IF(IABS(SDSWIT).LT.10) MAIN 553 1 CALL FITP(MX,MY,SDSWIT,LFITRM) MAIN 554 C MAIN 555 C CALCULATE U, T, AND STRESS MAIN 556 C MAIN 557 3200 U=0.0 MAIN 558 T=0.0 MAIN 559 DO 3210 M=MX,MY MAIN 560 U=U+(DIST(M)-DHAT(M))**2*WW(M) MAIN 561 3210 T=T+(DIST(M)-DBAS)**2*WW(M) MAIN 562 3215 U=SQRT (U) MAIN 563 TEMP1=T MAIN 564 T=SQRT (T) MAIN 565 GRSTRS(NG) = U/T MAIN 566 STRESS = STRESS + GRSTRS(NG)**2 MAIN 567 IF(U.EQ.0.0) GO TO 3340 MAIN 568 3220 RUT=1.0/(U*T) MAIN 569 UOT3=U/(TEMP1*T) MAIN 570 C MAIN 571 C CALCULATE THE (NEGATIVE) GRADIENT MAIN 572 C MAIN 573 3300 DO 3330 M = MX,MY MAIN 574 IF(DIST(M).EQ.0.0) GO TO 3330 MAIN 575 LTEMP1=MOD(IJ(M),KPACK) MAIN 576 I=LTEMP1/IJPACK+1 MAIN 577 J=MOD(LTEMP1,IJPACK)+1 MAIN 578 FACTOR=UOT3*(DIST(M)-DBAS)-RUT*(DIST(M)-DHAT(M)) MAIN 579 FACTOR = (FACTOR/RM1POW(DIST(M)) ) * WW(M) MAIN 580 FACTOR = FACTOR * GRSTRS(NG) MAIN 581 DO 3320 L=1,LDIM MAIN 582 TEMP1 = FACTOR * RM1POW(X(I,L)-X(J,L)) MAIN 583 GR(I,L)=GR(I,L)+TEMP1 MAIN 584 3320 GR(J,L)=GR(J,L)-TEMP1 MAIN 585 3330 CONTINUE MAIN 586 IF(SDSWIT.GE.10) GO TO 3340 MAIN 587 C SAVE POLYNOMIAL COEFFICIENTS MAIN 588 M=MX MAIN 589 DO 3335 M1=1,SDSWIT MAIN 590 DHAT(M)=PCOEFF(M1) MAIN 591 3335 M=M+1 MAIN 592 IF(LFITRM.EQ.2) DHAT(MX)=0.0 MAIN 593 3340 CONTINUE MAIN 594 IF(STRESS .EQ. 0.0 ) GO TO 3700 MAIN 595 STRESS = SQRT( STRESS / FNGRPS ) MAIN 596 C MAIN 597 C AVOID MOVING FIXED POINTS MAIN 598 C MAIN 599 IF( LNFIXZ .EQ. 0) GO TO 3400 MAIN 600 DO 3360 I=1,LNFIXZ MAIN 601 DO 3360 L=1,LDIM MAIN 602 3360 GR(I,L) = 0.0 MAIN 603 C MAIN 604 C FIND SCALE FACTOR OF GRADIENT, ANGLE-COSINE BETWEEN GRADIENT MAIN 605 C AND PREVIOUS GRADIENT. MAIN 606 C MAIN 607 3400 SFGR=0.0 MAIN 608 CAGRGL=0.0 MAIN 609 DO 3410 I=1,N MAIN 610 DO 3410 L=1,LDIM MAIN 611 SFGR=SFGR+GR(I,L)**2 MAIN 612 3410 CAGRGL=CAGRGL+GR(I,L)*GL(I,L) MAIN 613 SFGR=SQRT (SFGR/FN) MAIN 614 C IF GRADIENT = 0.0, SKIP AHEAD. MAIN 615 IF(SFGR) 3420,3700,3420 MAIN 616 3420 TEMP1=SFGR*SFGL*FN MAIN 617 CAGRGL=CAGRGL/TEMP1 MAIN 618 C MAIN 619 IF(ITNO.EQ.0 .OR. NBAKUP.GE.4) GO TO 3500 MAIN 620 IF(CAGRGL.GT.(-0.95) .AND. STRESS/STRLST.LT. 1.2001 ) GOTO 3500MAIN 621 C MAIN 622 C BACK UP ALONG LAST GRADIENT MAIN 623 C MAIN 624 NBAKUP = NBAKUP + 1 MAIN 625 IF(NBAKUP.EQ.1) STBASC=1.0 MAIN 626 STBASC=STBASC*STBAMU MAIN 627 TEMP1 = STEP MAIN 628 STEP = STEP / 10.0 MAIN 629 WRITE (LPRINT,38) STRESS, CAGRGL, STEP MAIN 630 38 FORMAT(10X, F7.3, 14X, F7.3, 22X, F8.4 ) MAIN 631 TEMP1 = (TEMP1 - STEP) / SFGL MAIN 632 TEMP1=STBASC*TEMP1 MAIN 633 DO 3430 I = 1,N MAIN 634 DO 3430 L = 1,LDIM MAIN 635 X(I,L) = X(I,L) - TEMP1*GL(I,L) MAIN 636 3430 GR(I,L) = GL(I,L) MAIN 637 SFGR = SFGL MAIN 638 STRESS = STRLST MAIN 639 GO TO 3000 MAIN 640 C MAIN 641 C STEP SIZE CALCULATIONS MAIN 642 C MAIN 643 3500 IF(ITNO) 9999, 3510, 3520 MAIN 644 3510 SRAT=0.8 MAIN 645 GO TO 3530 MAIN 646 3520 SRAT=STRESS/STRLST MAIN 647 3530 CALL NEWSTP( STEP, ITNO, SFGR, STRESS, MAIN 648 1 CAGRGL, COSAV, ACSAV, COSAVW, ACSAVW, SRAT, SRATAV ) MAIN 649 NBAKUP = 0 MAIN 650 C MAIN 651 C PRINT CURRENT STATUS OF COMPUTATION MAIN 652 C MAIN 653 3700 IF(LHIPRT.EQ.2) WRITE (LPRINT,30) ITNO,STRESS,SRAT,SRATAV, MAIN 654 1 CAGRGL,COSAV,ACSAV,SFGR,STEP MAIN 655 30 FORMAT(I10,F7.3,F7.3,F7.3,F7.3,F7.3,F7.3,F8.4, F8.4) MAIN 656 C MAIN 657 C DECIDE WHETHER TO CONTINUE ITERATING MAIN 658 C MAIN 659 3800 IF(STRESS) 9999, 3840, 3810 MAIN 660 3810 IF(SFGR-SFGRMN) 3850, 3850, 3815 MAIN 661 3815 TEMP1 = 0.5 * (1.0+SRATST) MAIN 662 TEMP2 = 1.0 - TEMP1 MAIN 663 IF( ABS (SRAT-TEMP1) - TEMP2 ) 3816, 3816, 3820 MAIN 664 3816 IF( ABS (SRATAV-TEMP1) - TEMP2 ) 3850, 3850, 3820 MAIN 665 3820 IF(STRESS-STRMIN) 3860, 3860, 3830 MAIN 666 3830 IF(ITNO-NOIT) 3900, 3870, 9999 MAIN 667 3840 CONTINUE MAIN 668 WRITE (LPRINT, 31) MAIN 669 31 FORMAT(24H0ZERO STRESS WAS REACHED ) MAIN 670 GO TO 4000 MAIN 671 3850 CONTINUE MAIN 672 WRITE (LPRINT, 32) MAIN 673 32 FORMAT(21H0MINIMUM WAS ACHIEVED ) MAIN 674 GO TO 4000 MAIN 675 3860 CONTINUE MAIN 676 WRITE (LPRINT, 33) MAIN 677 33 FORMAT(32H0SATISFACTORY STRESS WAS REACHED ) MAIN 678 GO TO 4000 MAIN 679 3870 CONTINUE MAIN 680 WRITE (LPRINT, 34) MAIN 681 34 FORMAT(39H0MAXIMUM NUMBER OF ITERATIONS WERE USED ) MAIN 682 GO TO 4000 MAIN 683 C MAIN 684 C CONTINUE ITERATING MAIN 685 C MAIN 686 3900 ITNO=ITNO+1 MAIN 687 TEMP1=STEP/SFGR MAIN 688 DO 3910 I=1,N MAIN 689 DO 3910 L=1,LDIM MAIN 690 3910 X(I,L)=X(I,L)+TEMP1*GR(I,L) MAIN 691 GO TO 3000 MAIN 692 C MAIN 693 C STOP ITERATING *************************************MAIN 694 C MAIN 695 C ROTATE FINAL CONFIGURATION TO PRINCIPAL COMPONENTS MAIN 696 C MAIN 697 4000 IF(LCOORS.LT.2) GO TO 4002 MAIN 698 C COMPUTE X TRANSPOSE TIMES X AND STORE UPPER HALF IN RVEC MAIN 699 KK=0 MAIN 700 DO 4410 K=1,LDIM MAIN 701 DO 4410 J=1,K MAIN 702 SUM=0.0 MAIN 703 DO 4405 I=1,N MAIN 704 4405 SUM=SUM+X(I,J)*X(I,K) MAIN 705 KK=KK+1 MAIN 706 4410 RVEC(KK)=SUM MAIN 707 C COMPUTE MATRIX OF EIGENVECTORS, GL MAIN 708 CALL SGEV(LDIM,LDIM,GL) MAIN 709 C COMPUTE X TIMES GL MAIN 710 DO 4420 K=1,LDIM MAIN 711 DO 4420 J=1,N MAIN 712 SUM=0.0 MAIN 713 DO 4415 I=1,LDIM MAIN 714 4415 SUM=SUM+X(J,I)*GL(I,K) MAIN 715 4420 GR(J,K)=SUM MAIN 716 DO 4425 I=1,N MAIN 717 DO 4425 J=1,LDIM MAIN 718 4425 X(I,J)=GR(I,J) MAIN 719 WRITE(LPRINT,152) MAIN 720 152 FORMAT(66H0THE FINAL CONFIGURATION HAS BEEN ROTATED TO PRINCIPAL CMAIN 721 1OMPONENTS.) MAIN 722 GO TO 4004 MAIN 723 C MAIN 724 4002 IF(LCOORS.EQ.1) WRITE(LPRINT,153) MAIN 725 153 FORMAT(110H0THE CONFIGURATION HAS BEEN NORMALIZED DURING THE ITERAMAIN 726 1TIONS BUT HAS NOT BEEN ROTATED TO PRINCIPAL COMPONENTS.) MAIN 727 IF(LCOORS.EQ.0) WRITE(LPRINT,154) MAIN 728 154 FORMAT(70H0NO NORMALIZATION WAS DONE TO THE CONFIGURATION DURING TMAIN 729 1HE ITERATIONS.) MAIN 730 4004 WRITE (LPRINT, 40)N,LDIM,STRESS,LSCH,(L,L=1,LDIM) MAIN 731 40 FORMAT(27H0THE FINAL CONFIGURATION OF,I4, MAIN 732 1 10H POINTS IN,I3, 22H DIMENSIONS HAS STRESS,F7.3,8H FORMULA ,I2 MAIN 733 2 /1X/30H LABEL FOR CONFIGURATION PLOTS,10X,19HFINAL CONFIGURATION MAIN 734 3 / 39X,10I7) MAIN 735 IF(LPUNSW.EQ.2) GO TO 4005 MAIN 736 WRITE (LPUNCH,41) (TITLE(I),I=1,80),N,LDIM MAIN 737 41 FORMAT (14H CONFIGURATION/80A1/2I3) MAIN 738 WRITE (LPUNCH,52) MAIN 739 52 FORMAT (12H (2X,10F7.3)) MAIN 740 4005 DO 4010 I=1,N MAIN 741 WRITE(LPRINT,42) PTID(I),I,(X(I,L),L=1,LDIM) MAIN 742 IF(LPUNSW.EQ.2) GO TO 4010 MAIN 743 WRITE (LPUNCH, 43)I,(X(I,L),L=1,LDIM) MAIN 744 4010 CONTINUE MAIN 745 42 FORMAT(20X,A1,18X,I2,10F7.3) MAIN 746 43 FORMAT(I2,10F7.3) MAIN 747 WRITE (LPRINT,46) MAIN 748 DO 4020 NG=1,NOGRPS MAIN 749 MX=GRNO(NG) MAIN 750 MZ=GRNO(NG+1)-MX MAIN 751 SDSWIT=GRSDIS(NG) MAIN 752 LFITRM=SDSWIT/100 MAIN 753 SDSWIT=MOD(SDSWIT,100) MAIN 754 IF(SDSWIT-10) 4016,4013,4014 MAIN 755 4013 WRITE(LPRINT,60) NG,MZ,GRSTRS(NG) MAIN 756 60 FORMAT(1X,I5,2X,I5,F7.3,11H ASCENDING) MAIN 757 GO TO 4020 MAIN 758 4014 WRITE(LPRINT,62) NG,MZ,GRSTRS(NG) MAIN 759 62 FORMAT(1X,I5,2X,I5,F7.3,11H DESCENDING) MAIN 760 GO TO 4020 MAIN 761 4016 M1=SDSWIT+MX-1 MAIN 762 M=MX MAIN 763 DO 4017 I=1,SDSWIT MAIN 764 PCOEFF(I)=DHAT(M) MAIN 765 4017 M=M+1 MAIN 766 WRITE(LPRINT,61) NG,MZ,GRSTRS(NG),(PCOEFF(I),I=1,SDSWIT) MAIN 767 61 FORMAT(1X,I5,2X,I5,F7.3,12H POLYNOMIAL,5F15.7) MAIN 768 C RESTORE DHAT VALUES MAIN 769 DO 4019 M=MX,M1 MAIN 770 DA=DATA(M) MAIN 771 TEMP=0.0 MAIN 772 DO 4018 I=LFITRM,SDSWIT MAIN 773 4018 TEMP=TEMP+PCOEFF(I)*REGR(DA,I) MAIN 774 4019 DHAT(M)=TEMP MAIN 775 4020 CONTINUE MAIN 776 46 FORMAT(14H0DATA GROUP(S) /73H0SERIAL COUNT STRESS REGRESSION COEMAIN 777 .FFICIENTS (FROM DEGREE 0 TO MAX OF 4)) MAIN 778 4030 CONTINUE MAIN 779 IF(LDIPRT.EQ.2) MAIN 780 1 CALL DATAPR(GRNO,MM,NOGRPS,2) MAIN 781 C MAIN 782 C PLOTTING SECTION MAIN 783 C MAIN 784 IF(LPLSCT.EQ.0) GO TO 5100 MAIN 785 C MAIN 786 C PLOT SCATTER DIAGRAMS OF DIST AND DHAT VS. DATA MAIN 787 PTID(1)=DEE MAIN 788 PTID(2)=AIE MAIN 789 C MAIN 790 C PLOT UNSPLIT SCATTER DIAGRAMS MAIN 791 C IF NOGRPS=1 THERE IS ONLY ONE SPLIT MAIN 792 C IF LPLSCT=3 PRODUCE A JOINT PLOT OF ALL SPLITS MAIN 793 C MAIN 794 IF(NOGRPS.GT.1) GO TO 5091 MAIN 795 WRITE(LPRINT,5093) LDIM,LSCH,STRESS,TITLE MAIN 796 5093 FORMAT(52H1DIST(D) AND DHAT(-) (Y-AXIS) VS. DATA (X-AXIS), FOR,I3,MAIN 797 . 27H DIMENSIONS. STRESS,FORMULA,I2,2H,=,F8.4/ 30X,80A1) MAIN 798 LNYSCT=2 MAIN 799 CALL PLOT(DATA,PMAT,0.,0.,0.,0.,MM,LNYSCT,-1,1800,1) MAIN 800 GO TO 5094 MAIN 801 C MAIN 802 5091 IF(LPLSCT.LT.3) GO TO 5092 MAIN 803 WRITE(LPRINT,44) LDIM,LSCH,STRESS,TITLE MAIN 804 44 FORMAT(45H1PLOT OF DIST (Y-AXIS) VS. DATA (X-AXIS), FOR,I3, MAIN 805 . 27H DIMENSIONS. STRESS,FORMULA,I2,2H,=,F8.4/ 30X,80A1) MAIN 806 LNYSCT=-2 MAIN 807 CALL PLOT(DATA,PMAT,0.,0.,0.,0.,MM,LNYSCT,-1,1800,1) MAIN 808 GO TO 5094 MAIN 809 C MAIN 810 C PLOT SPLIT SCATTER DIAGRAMS MAIN 811 C IF LPLSCT=1 PLOT NOT MORE THAN THE FIRST FIVE SPLITS MAIN 812 C IF LPLSCT=2 PLOT ALL SPLITS MAIN 813 C MAIN 814 5092 KPLT=NOGRPS MAIN 815 IF(LPLSCT.EQ.1) KPLT=MIN0(KPLT,5) MAIN 816 DO 5060 NG=1,KPLT MAIN 817 MX=GRNO(NG) MAIN 818 MY=GRNO(NG+1)-1 MAIN 819 MZ=MY-MX+1 MAIN 820 WRITE(LPRINT,4062) NG,LSCH,GRSTRS(NG),LDIM,TITLE MAIN 821 4062 FORMAT(43H1DIST(D) AND DHAT(-) VS. DATA FOR GROUP NO.,I3,17H. STRMAIN 822 .ESS,FORMULA,I2,2H,=,F8.4,12H DIMENSION=,I3 / 30X,80A1) MAIN 823 LNYSCT=2 MAIN 824 5060 CALL PLOT(DATA,PMAT,0.,0.,0.,0.,MZ,LNYSCT,-1,1800,MX) MAIN 825 5094 PTID(1)=AAA MAIN 826 PTID(2)=BEE MAIN 827 C MAIN 828 C PLOT CONFIGURATION MAIN 829 C MAIN 830 5100 IF(LPLCON.EQ.0)GO TO 5220 MAIN 831 IF(LDIM.EQ.1)GO TO 5190 MAIN 832 DO 5138 J=1,N MAIN 833 5138 PTID2(J)=PTID(J) MAIN 834 5125 IF(LPLCON.EQ.1) GO TO 5160 MAIN 835 C MAIN 836 C PLOT ALL PAIRS OF DIMENSIONS MAIN 837 C MAIN 838 DO 5145 J=2,LDIM MAIN 839 JMIN1=J-1 MAIN 840 DO 5145 I=1,JMIN1 MAIN 841 DO 5150 K=1,N MAIN 842 X2(K)=X(K,I) MAIN 843 Y2(K)=X(K,J) MAIN 844 5150 PTID(K)=PTID2(K) MAIN 845 WRITE(LPRINT,5015) J,I,TITLE MAIN 846 5015 FORMAT(31H1CONFIGURATION PLOT. DIMENSION,I2,23H (Y-AXIS) VS. DIMEMAIN 847 .NSION,I3,9H (X-AXIS) / 30X,80A1) MAIN 848 CALL SORT(Y2,N,PTID,X2,DUMEST,2,1) MAIN 849 LNYSCT=1 MAIN 850 5145 CALL PLOT(X2,Y2,0.,0.,0.,0.,N,LNYSCT,2,100,1) MAIN 851 GO TO 5200 MAIN 852 C MAIN 853 C PLOT PAIRS OF DIMENSIONS SO THAT EACH DIMEN PLOTTED ONLY ONCE MAIN 854 C MAIN 855 5160 LNYSCT=1 MAIN 856 DO 5185 J=1,LDIM,2 MAIN 857 JJ1=J+1 MAIN 858 JJ2=J MAIN 859 IF((MOD(LDIM,2).EQ.0).OR.(J.EQ.1)) GO TO 5175 MAIN 860 JJ1=J MAIN 861 JJ2=J-1 MAIN 862 IF(J.EQ.3) JJ2=1 MAIN 863 5175 DO 5180 K=1,N MAIN 864 X2(K)=X(K,JJ2) MAIN 865 Y2(K)=X(K,JJ1) MAIN 866 5180 PTID(K)=PTID2(K) MAIN 867 WRITE(LPRINT,5015) JJ1,JJ2,TITLE MAIN 868 CALL SORT(Y2,N,PTID,X2,DUMEST,2,1) MAIN 869 5185 CALL PLOT(X2,Y2,0.,0.,0.,0.,N,LNYSCT,2,100,1) MAIN 870 GO TO 5200 MAIN 871 C MAIN 872 C ONE DIMENSIONAL CONFIGURATION PLOT MAIN 873 C MAIN 874 5190 WRITE(LPRINT,5195) TITLE MAIN 875 5195 FORMAT(1H1,20X,34HPLOT OF CONFIGURATION, DIMENSION 1/30X,80A1) MAIN 876 LNYSCT=1 MAIN 877 CALL PLOT(X(1,1),X(1,1),0.,0.,0.,0.,N,LNYSCT,3,100,1) MAIN 878 GO TO 5220 MAIN 879 C MAIN 880 C CHANGE DIMENSION MAIN 881 C MAIN 882 5200 DO 5210 I=1,N MAIN 883 5210 PTID(I)=PTID2(I) MAIN 884 5220 ISTC=ISTC+1 MAIN 885 STPL(ISTC)= STRESS MAIN 886 DIMSV(ISTC)=LDIM MAIN 887 4100 LDIM=LDIM-LDIMD MAIN 888 IF(LDIM-LDIMN)4110,4101,4101 MAIN 889 4101 WRITE(LPRINT,4102) MAIN 890 4102 FORMAT(1H1) MAIN 891 GO TO 2300 MAIN 892 C MAIN 893 C PLOT STRESS VERSUS DIMENSION IF MORE THAN TWO POINTS MAIN 894 C MAIN 895 4110 CONTINUE MAIN 896 IF(ISTC.LT.3) GO TO 100 MAIN 897 WRITE(LPRINT, 50) TITLE MAIN 898 50 FORMAT (1H1,10X,31HPLOT OF STRESS VERSUS DIMENSION / 30X,80A1) MAIN 899 C INVERT TO IMPROVE READING OF PLOT MAIN 900 NVT=ISTC/2 MAIN 901 DO 8050 I=1,NVT MAIN 902 IC=ISTC+1-I MAIN 903 TEMP = DIMSV(I) MAIN 904 SEMP = STPL(I) MAIN 905 DIMSV(I) = DIMSV(IC) MAIN 906 STPL(I) =STPL(IC) MAIN 907 DIMSV(IC) = TEMP MAIN 908 STPL(IC) = SEMP MAIN 909 8050 CONTINUE MAIN 910 C FIND MAX STRESS FOR Y-AXIS , AND FILL IN PLOTTING CHARACTERS MAIN 911 YMA=0.0 MAIN 912 DO 8231 I=1,ISTC MAIN 913 IDIM=DIMSV(I) MAIN 914 PTID(I)=XNUM(IDIM) MAIN 915 8231 IF(YMA .LE. STPL(I)) YMA=STPL(I) MAIN 916 CALL SORT(STPL,ISTC,PTID,DIMSV,DUMEST,2,1) MAIN 917 LNYSCT=1 MAIN 918 CALL PLOT(DIMSV,STPL,10.,YMA,0.,0.,ISTC,LNYSCT,-2,10,1) MAIN 919 PTID(1)=AAA MAIN 920 PTID(2)=BEE MAIN 921 PTID(3)=CEE MAIN 922 PTID(4)=DEE MAIN 923 PTID(5)=EEE MAIN 924 PTID(6)=FFF MAIN 925 C MAIN 926 C REINITIALIZE, AND RETURN FOR MORE INPUT. MAIN 927 GO TO 100 MAIN 928 C MAIN 929 C NORMAL TERMINATION, AFTER READING 'STOP' ONCONTROL CARD MAIN 930 9000 STOP MAIN 931 C MAIN 932 C TROUBLE EXIT MAIN 933 C MAIN 934 9999 WRITE (LPRINT, 99) MAIN 935 99 FORMAT(52H0KRUSKAL. IMPOSSIBLE BRANCH TAKEN FROM IF STATEMENT. ) MAIN 936 STOP MAIN 937 C MAIN 938 END MAIN 939 CBLOCK DATA BLOCK 1 BLOCK DATA BLOCK 2 C BLDA FOR KYST JANUARY,1973 BLOCK 3 C WRITTEN BY J.KRUSKAL BLOCK 4 C MODIFIED FOR KYST BY J.KRUSKAL AND J.SEERY JANUARY,1973 BLOCK 5 C KYST VERSION 2 OCTOBER, 1976 BLOCK 6 C BLOCK 7 C THIS PROGRAM HOLDS THE TABLE OF WORDS WHICH CCACT CONSULTS BLOCK 8 C LOGICALLY IT CONSISTS OF SEVERAL TABLES BLOCK 9 C THE VALUES IN MTAB INDICATE WHICH ROWS START NEW TABLES BLOCK 10 C EACH ENTRY IN THE TABLE CONTAINS 4 ITEMS BLOCK 11 C THE FIRST IS THE CODED ALPHABETIC WORD BLOCK 12 C (ONE WORD MAY BE ENTERED SEVERAL TIMES ) BLOCK 13 C THE SECOND INDICATES THE NATURE OF THIS ENTRY BLOCK 14 C 1 INDICATES INTEGER PARAMETERS WHICH MUST BE EXPLIC9TLY BLOCK 15 C READ IN BLOCK 16 C 2 INDICATES REAL PARAMETERS WHICH MUST BE EXPLICITLY BLOCK 17 C READ IN BLOCK 18 C 3 INDICATES AN INTEGER PARAMETER WHICH IS IMPLICITLY BLOCK 19 C DEFINED BLOCK 20 C 4 INDICATES A REAL PARAMETER WHICH IS IMPLICITLY BLOCK 21 C DEFINED BLOCK 22 C 5 INDICATES A PARAMETER WHICH BELONGS TO CCACT ONLY BLOCK 23 C 6 INDICATES AN IMPLICITLY DEFINED INTEGER PARAMETER BLOCK 24 C EFFECTING PROGRAM EXECUTION UPON RETURN FROM CCACT BLOCK 25 C THE THIRD INDICATES WHICH PARAMTER IS INVOLVED BLOCK 26 C SEPARATE NUMBERING FOR MAIN PROGRAM PARAMETERS AND CCACT BLOCK 27 C THE FOURTH GIVES THE VALUE FOR ANY INTERNALLY DEFINED PARAMETERS BLOCK 28 C IF THE PARAMETER IS OF TYPE 3,4 OR 6.IF OF TYPE 1 OR 2, THIS BLOCK 29 C ENTRY GIVES THE NUMBER OF ITEMS TO BE READ IN. BLOCK 30 C BLOCK 31 INTEGER LPAR( 44 ) BLOCK 32 REAL PAR( 44 ) BLOCK 33 C BLOCK 34 INTEGER MTAB(13) BLOCK 35 INTEGER TAB1(4,18),TAB2(4,17),TAB3(4,16),TAB4(4,19),TAB5(4,19) BLOCK 36 INTEGER PTID(100),ITEM(101) BLOCK 37 C BLOCK 38 COMMON /ACCUR/ PRECSN,XMAG,XMAXN,FLPRIN,CHPRIN BLOCK 39 COMMON /MDSCL1/ LREAD, LPRINT, LPUNCH, LSCRAT BLOCK 40 COMMON /MDSCL2/ LPAR BLOCK 41 COMMON /MDSCL3/ MTAB,TAB1,TAB2,TAB3,TAB4,TAB5 BLOCK 42 COMMON /PLTCHR/ PTID,ITEM BLOCK 43 EQUIVALENCE (LPAR(1),PAR(1)) BLOCK 44 C BLOCK 45 DATA PRECSN, XMAG, XMAXN, FLPRIN, CHPRIN BLOCK 46 . /1.5E-8, 1.0E-38, 1.0E38, 6.0, 10.0/, BLOCK 47 . LREAD,LPRINT,LPUNCH,LSCRAT /5,6,7,8/ BLOCK 48 DATA MTAB(1),MTAB(2),MTAB(3),MTAB(4),MTAB(5),MTAB(6),MTAB(7), BLOCK 49 . MTAB(8),MTAB(9),MTAB(10),MTAB(11),MTAB(12),MTAB(13) BLOCK 50 . /1,52,54,59,65,67,76,78,80,83,87,90,90/ BLOCK 51 C BLOCK 52 C TAB1 ENTRIES PERTAIN SUCCESSIVELY TO THESE CONTROL WORDS: BLOCK 53 C DIMMAX, DIMMIN, DIMDIF, CUTOFF, STRMIN, BLOCK 54 C SFGRMN, COSAVW, ACSAVW, DIAGON, MATRIX, BLOCK 55 C HALFMA, LOWERH, UPPERH, LOWERC, ARBITR, BLOCK 56 C TORSCA, DATAFI, SPREAD BLOCK 57 DATA TAB1(1,1),TAB1(2,1),TAB1(3,1),TAB1(4,1),TAB1(1,2),TAB1(2,2) BLOCK 58 .,TAB1(3,2),TAB1(4,2),TAB1(1,3),TAB1(2,3),TAB1(3,3),TAB1(4,3) BLOCK 59 .,TAB1(1,4),TAB1(2,4),TAB1(3,4),TAB1(4,4),TAB1(1,5),TAB1(2,5) BLOCK 60 .,TAB1(3,5),TAB1(4,5),TAB1(1,6),TAB1(2,6),TAB1(3,6),TAB1(4,6) BLOCK 61 .,TAB1(1,7),TAB1(2,7),TAB1(3,7),TAB1(4,7),TAB1(1,8),TAB1(2,8) BLOCK 62 .,TAB1(3,8),TAB1(4,8),TAB1(1,9),TAB1(2,9),TAB1(3,9),TAB1(4,9) BLOCK 63 .,TAB1(1,10),TAB1(2,10),TAB1(3,10),TAB1(4,10),TAB1(1,11) BLOCK 64 .,TAB1(2,11),TAB1(3,11),TAB1(4,11)/ BLOCK 65 1 6989, 1, 1, 1, BLOCK 66 2 5375, 1, 2, 1, BLOCK 67 3 6923, 1, 3, 1, BLOCK 68 4 13520, 2, 4, 1, BLOCK 69 5 390, 2, 5, 1, BLOCK 70 6 2553, 2, 6, 1, BLOCK 71 7 7303, 2, 7, 1, BLOCK 72 8 11226, 2, 8, 1, BLOCK 73 9 11136, 5, 1, 2, BLOCK 74 A 9277, 3, 10, 1, BLOCK 75 B 3527, 3, 10, 2/ BLOCK 76 DATA TAB1(1,12),TAB1(2,12),TAB1(3,12),TAB1(4,12),TAB1(1,13), BLOCK 77 . TAB1(2,13),TAB1(3,13),TAB1(4,13),TAB1(1,14),TAB1(2,14), BLOCK 78 . TAB1(3,14),TAB1(4,14),TAB1(1,15),TAB1(2,15),TAB1(3,15), BLOCK 79 . TAB1(4,15),TAB1(1,16),TAB1(2,16),TAB1(3,16),TAB1(4,16), BLOCK 80 . TAB1(1,17),TAB1(2,17),TAB1(3,17),TAB1(4,17),TAB1(1,18), BLOCK 81 . TAB1(2,18),TAB1(3,18),TAB1(4,18)/ BLOCK 82 C 6069, 3, 10, 2, BLOCK 83 D 8938, 3, 10, 3, BLOCK 84 E 1754, 3, 10, 4, BLOCK 85 F 10499, 3, 20, -99, BLOCK 86 G 2048, 3, 27, 4, BLOCK 87 H 3323, 1, 43, 1, BLOCK 88 I 6910, 2, 44, 1/ BLOCK 89 C BLOCK 90 C TAB2 ENTRIES PERTAIN SUCCESSIVELY TO THESE CONTROL WORDS: BLOCK 91 C UPPERC, BLOCKD, DISSIM, DISSIM, SIMILA, BLOCK 92 C SIMILA, CONFIG, COMPUT, PRIMAR, SECOND, BLOCK 93 C R , ITERAT, FIX , SRATST, STOP , BLOCK 94 C PRE-IT, COORDI BLOCK 95 DATA TAB2(1,1),TAB2(2,1),TAB2(3,1),TAB2(4,1),TAB2(1,2),TAB2(2,2) BLOCK 96 .,TAB2(3,2),TAB2(4,2),TAB2(1,3),TAB2(2,3),TAB2(3,3),TAB2(4,3) BLOCK 97 .,TAB2(1,4),TAB2(2,4),TAB2(3,4),TAB2(4,4),TAB2(1,5),TAB2(2,5) BLOCK 98 .,TAB2(3,5),TAB2(4,5),TAB2(1,6),TAB2(2,6),TAB2(3,6),TAB2(4,6) BLOCK 99 .,TAB2(1,7),TAB2(2,7),TAB2(3,7),TAB2(4,7),TAB2(1,8),TAB2(2,8) BLOCK100 .,TAB2(3,8),TAB2(4,8),TAB2(1,9),TAB2(2,9),TAB2(3,9),TAB2(4,9) BLOCK101 .,TAB2(1,10),TAB2(2,10),TAB2(3,10),TAB2(4,10),TAB2(1,11) BLOCK102 .,TAB2(2,11),TAB2(3,11),TAB2(4,11)/ BLOCK103 1 4623, 3, 10, 5, BLOCK104 2 8683, 5, 1, 5, BLOCK105 3 15096, 6, 12, 2, BLOCK106 4 15096, 3, 11, 10, BLOCK107 5 6868, 6, 12, 2, BLOCK108 6 6868, 3, 11, 11, BLOCK109 7 14846, 6, 12, 5, BLOCK110 8 11754, 6, 12, 6, BLOCK111 9 765, 3, 13, 1, BLOCK112 A 4119, 3, 13, 2, BLOCK113 B 16326, 2, 14, 1/ BLOCK114 DATA TAB2(1,12),TAB2(2,12),TAB2(3,12),TAB2(4,12),TAB2(1,13), BLOCK115 . TAB2(2,13),TAB2(3,13),TAB2(4,13),TAB2(1,14),TAB2(2,14), BLOCK116 . TAB2(3,14),TAB2(4,14),TAB2(1,15),TAB2(2,15),TAB2(3,15), BLOCK117 . TAB2(4,15),TAB2(1,16),TAB2(2,16),TAB2(3,16),TAB2(4,16), BLOCK118 . TAB2(1,17),TAB2(2,17),TAB2(3,17),TAB2(4,17)/ BLOCK119 C 15170, 1, 15, 1, BLOCK120 D 1855, 1, 28, 1, BLOCK121 E 3720, 2, 16, 1, BLOCK122 F 13171, 6, 12, 7, BLOCK123 G 10738, 1, 39, 1, BLOCK124 H 7261, 5, 1, 11/ BLOCK125 C BLOCK126 C TAB3 ENTRIES PERTAIN SUCCESSIVELY TO THESE CONTROL WORDS: BLOCK127 C WEIGHT, WFUNCT, SFORM1, SFORM2, CARDS, BLOCK128 C NOCARD, SPLIT , DATA , SAVE , RANDOM, BLOCK129 C PRINT , DFUNCT, REGRES, DCONST, WCONST, BLOCK130 C PLOT BLOCK131 DATA TAB3(1,1),TAB3(2,1),TAB3(3,1),TAB3(4,1),TAB3(1,2),TAB3(2,2) BLOCK132 .,TAB3(3,2),TAB3(4,2),TAB3(1,3),TAB3(2,3),TAB3(3,3),TAB3(4,3) BLOCK133 .,TAB3(1,4),TAB3(2,4),TAB3(3,4),TAB3(4,4),TAB3(1,5),TAB3(2,5) BLOCK134 .,TAB3(3,5),TAB3(4,5),TAB3(1,6),TAB3(2,6),TAB3(3,6),TAB3(4,6) BLOCK135 .,TAB3(1,7),TAB3(2,7),TAB3(3,7),TAB3(4,7),TAB3(1,8),TAB3(2,8) BLOCK136 .,TAB3(3,8),TAB3(4,8),TAB3(1,9),TAB3(2,9),TAB3(3,9),TAB3(4,9) BLOCK137 .,TAB3(1,10),TAB3(2,10),TAB3(3,10),TAB3(4,10),TAB3(1,11) BLOCK138 .,TAB3(2,11),TAB3(3,11),TAB3(4,11)/ BLOCK139 1 14543, 6, 12, 3, BLOCK140 2 11427, 6, 12, 4, BLOCK141 3 5318, 3, 17, 1, BLOCK142 4 6181, 3, 17, 2, BLOCK143 5 6927, 3, 18, 1, BLOCK144 6 16009, 3, 18, 2, BLOCK145 7 1966, 5, 1, 3, BLOCK146 8 6675, 6, 12, 2, BLOCK147 9 9189, 5, 1, 7, BLOCK148 A 8330, 1, 20, 1, BLOCK149 B 2775, 5, 1, 4/ BLOCK150 DATA TAB3(1,12),TAB3(2,12),TAB3(3,12),TAB3(4,12),TAB3(1,13), BLOCK151 . TAB3(2,13),TAB3(3,13),TAB3(4,13),TAB3(1,14),TAB3(2,14), BLOCK152 . TAB3(3,14),TAB3(4,14),TAB3(1,15),TAB3(2,15),TAB3(3,15), BLOCK153 . TAB3(4,15),TAB3(1,16),TAB3(2,16),TAB3(3,16),TAB3(4,16)/ BLOCK154 C 10575, 6, 12, 8, BLOCK155 D 14439, 5, 1, 6, BLOCK156 E 267, 2, 29, 5, BLOCK157 F 1119, 2, 34, 5, BLOCK158 G 6878, 5, 1, 8/ BLOCK159 C BLOCK160 C TAB4 ENTRIES PERTAIN SUCCESSIVELY TO THESE CONTROL WORDS: BLOCK161 C ABSENT, PRESEN, BYROWS, BYGROU, BYDECK, BLOCK162 C NOMORE, NOMORE, DATA , DISTAN, HISTOR, BLOCK163 C NODATA, NODIST, NOHIST, NO , YES , BLOCK164 C MONOTO, ASCEND, DESCEN, POLYNO BLOCK165 DATA TAB4(1,1),TAB4(2,1),TAB4(3,1),TAB4(4,1),TAB4(1,2),TAB4(2,2) BLOCK166 .,TAB4(3,2),TAB4(4,2),TAB4(1,3),TAB4(2,3),TAB4(3,3),TAB4(4,3) BLOCK167 .,TAB4(1,4),TAB4(2,4),TAB4(3,4),TAB4(4,4),TAB4(1,5),TAB4(2,5) BLOCK168 .,TAB4(3,5),TAB4(4,5),TAB4(1,6),TAB4(2,6),TAB4(3,6),TAB4(4,6) BLOCK169 .,TAB4(1,7),TAB4(2,7),TAB4(3,7),TAB4(4,7),TAB4(1,8),TAB4(2,8) BLOCK170 .,TAB4(3,8),TAB4(4,8),TAB4(1,9),TAB4(2,9),TAB4(3,9),TAB4(4,9) BLOCK171 .,TAB4(1,10),TAB4(2,10),TAB4(3,10),TAB4(4,10),TAB4(1,11) BLOCK172 .,TAB4(2,11),TAB4(3,11),TAB4(4,11)/ BLOCK173 1 4258, 3, 9, 2, BLOCK174 2 2575, 3, 9, 1, BLOCK175 3 7761, 3, 19, 100, BLOCK176 4 11782, 3, 19, 200, BLOCK177 5 11288, 3, 19, 300, BLOCK178 6 5274, 3, 19, 400, BLOCK179 7 5274, 3, 19, 2, BLOCK180 8 6675, 3, 21, 2, BLOCK181 9 9824, 3, 22, 2, BLOCK182 A 12801, 3, 24, 2, BLOCK183 B 16057, 3, 21, 1/ BLOCK184 DATA TAB4(1,12),TAB4(2,12),TAB4(3,12),TAB4(4,12),TAB4(1,13), BLOCK185 . TAB4(2,13),TAB4(3,13),TAB4(4,13),TAB4(1,14),TAB4(2,14), BLOCK186 . TAB4(3,14),TAB4(4,14),TAB4(1,15),TAB4(2,15),TAB4(3,15), BLOCK187 . TAB4(4,15),TAB4(1,16),TAB4(2,16),TAB4(3,16),TAB4(4,16), BLOCK188 . TAB4(1,17),TAB4(2,17),TAB4(3,17),TAB4(4,17),TAB4(1,18), BLOCK189 . TAB4(2,18),TAB4(3,18),TAB4(4,18),TAB4(1,19),TAB4(2,19), BLOCK190 . TAB4(3,19),TAB4(4,19)/ BLOCK191 C 5863, 3, 22, 1, BLOCK192 D 9395, 3, 24, 1, BLOCK193 E 9622, 3, 10, 100, BLOCK194 F 11125, 3, 10, 200, BLOCK195 G 15634, 3, 23, 0, BLOCK196 H 7782, 3, 23, 10, BLOCK197 I 11188, 3, 23, 11, BLOCK198 J 3756, 3, 26, 1/ BLOCK199 C BLOCK200 C TAB5 ENTRIES PERTAIN SUCCESSIVELY TO THESE CONTROL WORDS: BLOCK201 C POLYNO, CONSTA, NOCONS, MULTIV, MULTIV, BLOCK202 C DATA , CONFIG, CONFIG, SCATTE, ALL , BLOCK203 C SOME , NONE , ALL , SOME , NONE , BLOCK204 C JOINT , ROTATE, STANDA, AS-IS BLOCK205 DATA TAB5(1,1),TAB5(2,1),TAB5(3,1),TAB5(4,1),TAB5(1,2),TAB5(2,2) BLOCK206 .,TAB5(3,2),TAB5(4,2),TAB5(1,3),TAB5(2,3),TAB5(3,3),TAB5(4,3) BLOCK207 .,TAB5(1,4),TAB5(2,4),TAB5(3,4),TAB5(4,4),TAB5(1,5),TAB5(2,5) BLOCK208 .,TAB5(3,5),TAB5(4,5),TAB5(1,6),TAB5(2,6),TAB5(3,6),TAB5(4,6) BLOCK209 .,TAB5(1,7),TAB5(2,7),TAB5(3,7),TAB5(4,7),TAB5(1,8),TAB5(2,8) BLOCK210 .,TAB5(3,8),TAB5(4,8),TAB5(1,9),TAB5(2,9),TAB5(3,9),TAB5(4,9) BLOCK211 .,TAB5(1,10),TAB5(2,10),TAB5(3,10),TAB5(4,10),TAB5(1,11) BLOCK212 .,TAB5(2,11),TAB5(3,11),TAB5(4,11)/ BLOCK213 1 3756, 1, 23, 1, BLOCK214 2 14387, 3, 26, 100, BLOCK215 3 5018, 3, 26, 200, BLOCK216 4 3608, 3, 26, 0, BLOCK217 5 3608, 1, 23, 1, BLOCK218 6 6675, 3, 25, 2, BLOCK219 7 14846, 3, 27, 3, BLOCK220 8 14846, 5, 1, 9, BLOCK221 9 11109, 5, 1, 10, BLOCK222 A 5766, 3, 40, 2, BLOCK223 B 13660, 3, 40, 1/ BLOCK224 DATA TAB5(1,12),TAB5(2,12),TAB5(3,12),TAB5(4,12),TAB5(1,13), BLOCK225 . TAB5(2,13),TAB5(3,13),TAB5(4,13),TAB5(1,14),TAB5(2,14), BLOCK226 . TAB5(3,14),TAB5(4,14),TAB5(1,15),TAB5(2,15),TAB5(3,15), BLOCK227 . TAB5(4,15),TAB5(1,16),TAB5(2,16),TAB5(3,16),TAB5(4,16), BLOCK228 . TAB5(1,17),TAB5(2,17),TAB5(3,17),TAB5(4,17),TAB5(1,18), BLOCK229 . TAB5(2,18),TAB5(3,18),TAB5(4,18),TAB5(1,19),TAB5(2,19), BLOCK230 . TAB5(3,19),TAB5(4,19)/ BLOCK231 C 10008, 3, 40, 0, BLOCK232 D 5766, 3, 41, 2, BLOCK233 E 13660, 3, 41, 1, BLOCK234 F 10008, 3, 41, 0, BLOCK235 G 11053, 3, 41, 3, BLOCK236 H 4503, 3, 42, 2, BLOCK237 I 3418, 3, 42, 1, BLOCK238 J 9499, 3, 42, 0/ BLOCK239 C BLOCK240 DATA PTID(1),PTID(2),PTID(3),PTID(4),PTID(5),PTID(6),PTID(7), BLOCK241 . PTID(8),PTID(9),PTID(10),PTID(11),PTID(12),PTID(13),PTID(14), BLOCK242 . PTID(15),PTID(16),PTID(17),PTID(18),PTID(19),PTID(20),PTID(21), BLOCK243 . PTID(22),PTID(23),PTID(24),PTID(25),PTID(26),PTID(27),PTID(28), BLOCK244 . PTID(29),PTID(30),PTID(31),PTID(32),PTID(33),PTID(34),PTID(35), BLOCK245 . PTID(36),PTID(37),PTID(38),PTID(39),PTID(40),PTID(41),PTID(42), BLOCK246 . PTID(43),PTID(44),PTID(45),PTID(46),PTID(47),PTID(48),PTID(49), BLOCK247 . PTID(50)/ BLOCK248 . 1HA,1HB,1HC,1HD,1HE,1HF,1HG,1HH,1HI,1HJ,1HK,1HL,1HM, BLOCK249 1 1HN,1HO,1HP,1HQ,1HR,1HS,1HT,1HU,1HV,1HW,1HX,1HY, BLOCK250 2 1HZ,1H1,1H2,1H3,1H4,1H5,1H6,1H7,1H8,1H9,1HA,1HB,1HC, BLOCK251 3 1HD,1HE,1HF,1HG,1HH,1HI,1HJ,1HK,1HL,1HM,1HN,1HO/ BLOCK252 DATA PTID(51),PTID(52),PTID(53),PTID(54),PTID(55),PTID(56), BLOCK253 . PTID(57),PTID(58),PTID(59),PTID(60),PTID(61),PTID(62),PTID(63), BLOCK254 . PTID(64),PTID(65),PTID(66),PTID(67),PTID(68),PTID(69),PTID(70), BLOCK255 . PTID(71),PTID(72),PTID(73),PTID(74),PTID(75),PTID(76),PTID(77), BLOCK256 . PTID(78),PTID(79),PTID(80),PTID(81),PTID(82),PTID(83),PTID(84), BLOCK257 . PTID(85),PTID(86),PTID(87),PTID(88),PTID(89),PTID(90),PTID(91), BLOCK258 . PTID(92),PTID(93),PTID(94),PTID(95),PTID(96),PTID(97), BLOCK259 . PTID(98),PTID(99),PTID(100)/ BLOCK260 4 1HP,1HQ,1HR,1HS,1HT,1HU,1HV,1HW,1HX,1HY,1HZ,1H1,1H2, BLOCK261 5 1H3,1H4,1H5,1H6,1H7,1H8,1H9,1HA,1HB,1HC,1HD,1HE, BLOCK262 6 1HF,1HG,1HH,1HI,1HJ,1HK,1HL,1HM,1HN,1HO,1HP,1HQ, BLOCK263 7 1HR,1HS,1HT,1HU,1HV,1HW,1HX,1HY,1HZ,1H1,1H2,1H3,1H4/ BLOCK264 DATA ITEM(1),ITEM(2),ITEM(3),ITEM(4),ITEM(5),ITEM(6),ITEM(7), BLOCK265 . ITEM(8),ITEM(9),ITEM(10),ITEM(11),ITEM(12),ITEM(13),ITEM(14), BLOCK266 . ITEM(15),ITEM(16),ITEM(17),ITEM(18),ITEM(19),ITEM(20),ITEM(21), BLOCK267 . ITEM(22),ITEM(23),ITEM(24),ITEM(25),ITEM(26),ITEM(27),ITEM(28), BLOCK268 . ITEM(29),ITEM(30),ITEM(31),ITEM(32),ITEM(33),ITEM(34),ITEM(35), BLOCK269 . ITEM(36),ITEM(37),ITEM(38),ITEM(39),ITEM(40),ITEM(41),ITEM(42), BLOCK270 . ITEM(43),ITEM(44),ITEM(45),ITEM(46),ITEM(47),ITEM(48),ITEM(49)/ BLOCK271 1 1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H , BLOCK272 2 1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H , BLOCK273 3 1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H , BLOCK274 4 1H ,1H ,1H ,1H / BLOCK275 DATA BLOCK276 . ITEM(50),ITEM(51),ITEM(52),ITEM(53),ITEM(54),ITEM(55),ITEM(56), BLOCK277 . ITEM(57),ITEM(58),ITEM(59),ITEM(60),ITEM(61),ITEM(62),ITEM(63), BLOCK278 . ITEM(64),ITEM(65),ITEM(66),ITEM(67),ITEM(68),ITEM(69),ITEM(70), BLOCK279 . ITEM(71),ITEM(72),ITEM(73),ITEM(74),ITEM(75),ITEM(76),ITEM(77), BLOCK280 . ITEM(78),ITEM(79),ITEM(80),ITEM(81),ITEM(82),ITEM(83),ITEM(84), BLOCK281 . ITEM(85),ITEM(86),ITEM(87),ITEM(88),ITEM(89),ITEM(90),ITEM(91), BLOCK282 . ITEM(92),ITEM(93),ITEM(94),ITEM(95),ITEM(96),ITEM(97),ITEM(98), BLOCK283 . ITEM(99),ITEM(100),ITEM(101)/ BLOCK284 5 1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H , BLOCK285 6 1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H , BLOCK286 7 1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H ,1H , BLOCK287 8 1H ,1H ,1H ,1H ,1H ,1H ,1H / BLOCK288 END BLOCK289 CBSEC1 BSEC1 1 SUBROUTINE BSEC1(IB,N,ETA,A,B,W,E,S) BSEC1 2 C BSEC1 FOR KYST JANUARY,1973 BSEC1 3 C WRITTEN AND TECHNIQUES ADVOCATED IN BSEC1 4 C REFERENCE IMPLEMENTED FOR SGEV BSEC1 5 C BY P.BUSINGER JANUARY,1972 BSEC1 6 C KYST VERSION 2 OCTOBER, 1976 BSEC1 7 DIMENSION A(1),B(1),W(1) BSEC1 8 REAL ETA,E,S BSEC1 9 C BSEC1 10 C W.KAHAN, ACCURATE EIGENVALUES OF A SYM- BSEC1 11 C METRIC TRIDIAGONAL MATRIX. TECH.REPORT BSEC1 12 C NO.CS41, JULY 22,1966, COMP.SC.DEPT., BSEC1 13 C STANFORD UNIVERSITY. BSEC1 14 C BSEC1 15 C BSEC1 IS CALLED BY SGEV, AND BSEC1 16 C USES IB BISECTION STEPS TO FIND THE ALGEBRAICALLY BSEC1 17 C GREATEST EIGENVALUE E OF THE SYMMETRIC BSEC1 18 C TRIDIAGONAL MATRIX WHOSE DIAGONAL ELEMENTS ARE BSEC1 19 C A(1) THROUGH A(N) BSEC1 20 C AND WHOSE SUB-AND-SUPERDIAGONAL ELEMENTS ARE BSEC1 21 C B(1) THROUGH B(N-1) BSEC1 22 C FURTHERMORE, BSEC1 23 C W(1) THROUGH W(N) BSEC1 24 C ARE SET EQUAL TO THE DIAGONAL ELEMENTS OF THE BSEC1 25 C MATRIX U IN THE LU-DECOMPOSITION REQUIRED BSEC1 26 C FOR SUBSEQUENT INVERSE ITERATION. A AND B BSEC1 27 C ARE RESCALED, THE CHANGE IF SCALE BEING BSEC1 28 C REFLECTED IN S. ETA IS LEAST POSITIVE MA- BSEC1 29 C CHINE NUMBER. BSEC1 30 C BSEC1 31 REAL BB,D,T,X BSEC1 32 IF(N.GT.1)GOTO 5 BSEC1 33 E=A(1) BSEC1 34 W(1)=ETA BSEC1 35 GOTO 70 BSEC1 36 5 B(N)=0.E0 BSEC1 37 T=0.E0 BSEC1 38 DO 10 I=1,N BSEC1 39 T=AMAX1(T,ABS(A(I))) BSEC1 40 10 T=AMAX1(T,ABS(B(I))) BSEC1 41 IF(T.NE.0.E0)GOTO 20 BSEC1 42 E=0.E0 BSEC1 43 DO 19 I=1,N BSEC1 44 19 W(I)=ETA BSEC1 45 GOTO 70 BSEC1 46 20 S=S*T BSEC1 47 DO 30 I=1,N BSEC1 48 A(I)=A(I)/T BSEC1 49 30 B(I)=B(I)/T BSEC1 50 E=3.E0 BSEC1 51 D=-E BSEC1 52 DO 50 K=1,IB BSEC1 53 X=(D+E)/2.E0 BSEC1 54 U=1.E0 BSEC1 55 NU=0 BSEC1 56 BB=0.E0 BSEC1 57 DO 40 I=1,N BSEC1 58 U=(A(I)-BB/U)-X BSEC1 59 IF(U.GT.ETA)GOTO 35 BSEC1 60 U=AMIN1(U,-ETA) BSEC1 61 NU=NU+1 BSEC1 62 35 BB=B(I)**2 BSEC1 63 40 W(I)=U BSEC1 64 IF(NU.LT.N)D=X BSEC1 65 IF(NU.EQ.N)E=X BSEC1 66 50 CONTINUE BSEC1 67 IF(NU.EQ.N)GOTO 70 BSEC1 68 U=1.E0 BSEC1 69 BB=0.E0 BSEC1 70 DO 60 I=1,N BSEC1 71 U=(A(I)-BB/U)-E BSEC1 72 IF(U.GT.ETA)GOTO 55 BSEC1 73 U=AMIN1(U,-ETA) BSEC1 74 55 BB=B(I)**2 BSEC1 75 60 W(I)=U BSEC1 76 70 RETURN BSEC1 77 END BSEC1 78 CCACT CACT 1 SUBROUTINE CCACT CACT 2 C CCACT FOR KYST JANUARY,1973 CACT 3 C WRITTEN BY J.KRUSKAL CACT 4 C MODIFIED FOR KYST BY J.KRUSKAL AND J.SEERY JANUARY,1973 CACT 5 C KYST VERSION 2 OCTOBER, 1976 CACT 6 C CCACT--CONTROL CARD ACTIVITY. READS AND INTERPRETS CONTROL CARDS.CACT 7 C CACT 8 REAL DUMMY(27),ATAB(4,89) CACT 9 INTEGER TAB(4,89),IFAC(6) CACT 10 INTEGER MTAB(13) CACT 11 C CACT 12 INTEGER LPAR( 44 ) CACT 13 REAL PAR( 44 ) CACT 14 C CACT 15 INTEGER IPARAM(5) CACT 16 C CACT 17 INTEGER DFLTSW CACT 18 INTEGER C(81),BLANK,DOT,EQUALS,COMMA,DOLLAR,MAP(40),X CACT 19 INTEGER CHTAB(13),WORD(18) CACT 20 INTEGER BLSW,NUMSW,DECSW,TYPE,XTYPE,T,TA, PARNO, TABNO CACT 21 C CACT 22 COMMON /MDSCL1/ LREAD, LPRINT, LPUNCH, LSCRAT CACT 23 COMMON /MDSCL2/ LPAR CACT 24 COMMON /MDSCL3/ MTAB, TAB CACT 25 EQUIVALENCE (LPAR(1),PAR(1)) CACT 26 EQUIVALENCE (MAP(1),DUMMY(1)), (MAP(28),CHTAB(1)) CACT 27 EQUIVALENCE (ATAB(1),TAB(1)) CACT 28 C CACT 29 DATA BLANK, DOT, EQUALS, COMMA, C(81) /1H ,1H.,1H=,1H,,1H / CACT 30 DATA DOLLAR/1H$/ CACT 31 DATA MODP,IFAC(1),IFAC(2),IFAC(3),IFAC(4),IFAC(5),IFAC(6) CACT 32 1 / 16381, 907, 887, 883, 881, 877, 863/ CACT 33 C CACT 34 DATA MAP(1),MAP(2),MAP(3),MAP(4),MAP(5),MAP(6),MAP(7),MAP(8), CACT 35 . MAP(9),MAP(10),MAP(11),MAP(12),MAP(13),MAP(14),MAP(15),MAP(16), CACT 36 . MAP(17),MAP(18),MAP(19),MAP(20),MAP(21),MAP(22),MAP(23),MAP(24), CACT 37 . MAP(25),MAP(26),MAP(27),MAP(28),MAP(29),MAP(30),MAP(31),MAP(32), CACT 38 . MAP(33),MAP(34),MAP(35),MAP(36),MAP(37),MAP(38),MAP(39),MAP(40)/ CACT 39 . 1H ,1HA,1HB,1HC,1HD,1HE,1HF,1HG,1HH,1HI,1HJ,1HK,1HL, CACT 40 . 1HM,1HN,1HO,1HP,1HQ,1HR,1HS,1HT,1HU,1HV,1HW,1HX,1HY, CACT 41 . 1HZ,1H0,1H1,1H2,1H3,1H4,1H5,1H6,1H7,1H8,1H9,1H-,1H+,1H./CACT 42 C CACT 43 C CACT 44 C CACT 45 1 FORMAT(80A1) CACT 46 11 FORMAT(1H0,80A1) CACT 47 2 FORMAT(1X,I1) CACT 48 3 FORMAT(1X,80A1) CACT 49 4 FORMAT(1X,I18) CACT 50 5 FORMAT(1X,F18.3) CACT 51 9 FORMAT(98H CCACT DIAGNOSTIC. A CONTROL CARD WAS EXPECTED ABOVE, BCACT 52 .UT THIS CARD IS NOT A PROPER CONTROL CARD.) CACT 53 99 FORMAT( 12H ITEM NUMBER, I5, 20H HAS EXPECTED TYPE , I5, CACT 54 1 18H AND ACTUAL TYPE , I5, 16H . THIS ITEM IS , 18A1) CACT 55 98 FORMAT(12H ITEM NUMBER,I5,39H HAS ILLEGAL CHARACTER. THIS ITEM ISCACT 56 . ,18A1) CACT 57 C CACT 58 C READ AND PRINT CONTROL CARD CACT 59 C CACT 60 100 READ (LREAD,1) (C(I),I=1,80) CACT 61 WRITE(LPRINT,11) (C(I),I=1,80) CACT 62 C ANALYZE CONTROL CARD INTO TOKENS CACT 63 C CACT 64 C EACH 'TOKEN' IS DELIMITED BY BLANKS CACT 65 C THERE ARE FOUR TYPES OF TOKENS. CACT 66 C ALPHABETIC, INTEGER, DECIMAL, AND DEFAULT CACT 67 C ALPHABETIC UNLESS FIRST CHARACTER IS DIGIT OR DEC POINT CACT 68 C OR PLUS OR MINUS OR DOLLAR SIGN CACT 69 C DECIMALS DISTINGUISHED BY DECIMAL POINT CACT 70 C DEFAULT = $ CACT 71 C CACT 72 REWIND LSCRAT CACT 73 BLSW=1 CACT 74 NUMSW=1 CACT 75 DECSW=1 CACT 76 DFLTSW=1 CACT 77 K=0 CACT 78 C CACT 79 DO 300 I=1,81 CACT 80 C CACT 81 X=C(I) CACT 82 GO TO (110,120), BLSW CACT 83 C CACT 84 110 IF(X.EQ.BLANK .OR. X.EQ.EQUALS .OR. X.EQ. COMMA ) GO TO 300CACT 85 BLSW=2 CACT 86 JA=I CACT 87 DO 115 KX=1,13 CACT 88 115 IF(X.EQ.CHTAB(KX))NUMSW=2 CACT 89 IF(X.EQ.DOLLAR) DFLTSW=2 CACT 90 IF(X.EQ.DOT) DECSW=2 CACT 91 GO TO 300 CACT 92 C CACT 93 120 IF(X.EQ.BLANK .OR. X.EQ.EQUALS .OR. X.EQ. COMMA ) GO TO 130CACT 94 IF(X.EQ.DOT) DECSW=2 CACT 95 GO TO 300 CACT 96 C CACT 97 130 K=K+1 CACT 98 JB = MIN0 (I-1,JA+16) CACT 99 JC=18-(JB-JA+1) CACT 100 TYPE=1 CACT 101 IF(NUMSW.EQ.2) TYPE=NUMSW+DECSW-1 CACT 102 IF(DFLTSW.EQ.2) TYPE=4 CACT 103 C CACT 104 WRITE (LSCRAT,2) TYPE CACT 105 IF(TYPE.EQ.4) GO TO 160 CACT 106 GO TO (140,150),NUMSW CACT 107 140 WRITE (LSCRAT,3) (C(J),J=JA,JB), (BLANK,J=1,JC) CACT 108 GO TO 160 CACT 109 150 WRITE (LSCRAT,3) (BLANK,J=1,JC), (C(J),J=JA,JB) CACT 110 160 BLSW=1 CACT 111 NUMSW=1 CACT 112 DECSW=1 CACT 113 DFLTSW=1 CACT 114 GO TO 300 CACT 115 C CACT 116 300 CONTINUE CACT 117 C CACT 118 C ANALYZE TOKENS AND SET PARAMETER VALUES ACCORDINGLY CACT 119 C CACT 120 KB=K CACT 121 IF(KB.EQ.0) RETURN CACT 122 REWIND LSCRAT CACT 123 XTYPE=1 CACT 124 IPARAM(1) = 1 CACT 125 NT=0 CACT 126 NTOKNS=0 CACT 127 NOPCC=0 CACT 128 C CACT 129 DO 1000 K=1,KB CACT 130 C CACT 131 READ (LSCRAT,2) TYPE CACT 132 IF( (XTYPE.EQ.1) .AND. (TYPE.NE.1) ) GO TO 995 CACT 133 IF( (XTYPE.NE.TYPE). AND. (TYPE.NE.4) ) GO TO 995 CACT 134 350 GO TO (400,410,420,430), TYPE CACT 135 C CACT 136 400 READ (LSCRAT,3) (WORD(L),L=1,18) CACT 137 C CONVERT FIRST SIX LETTERS OF ALPHABETIC WORD INTO CODE NUMBER. CACT 138 C LARGEST POSSIBLE CODE NUMBER IS LESS THAN 2**15. CACT 139 ICODE=0 CACT 140 DO 408 M1=1,6 CACT 141 X=WORD(M1) CACT 142 DO 402 M2=1,38 CACT 143 IF( X.EQ.MAP(M2) ) GO TO 404 CACT 144 402 CONTINUE CACT 145 C ILLEGAL CHARACTER CACT 146 WRITE (LPRINT,9) CACT 147 WRITE (LPRINT,98) K,(WORD(L),L=1,18) CACT 148 GO TO 999 CACT 149 404 M3=(M2-1)*IFAC(M1) CACT 150 IF( M3.GE.MODP ) M3=M3-MODP CACT 151 ICODE=ICODE+M3 CACT 152 IF( ICODE.GE.MODP ) ICODE=ICODE-MODP CACT 153 408 CONTINUE CACT 154 GO TO 510 CACT 155 C CACT 156 410 READ (LSCRAT,4) INTPAR CACT 157 ISUB=PARNO+NT CACT 158 LPAR(ISUB)=INTPAR CACT 159 GO TO 430 CACT 160 C CACT 161 420 READ (LSCRAT,5) DECPAR CACT 162 ISUB=PARNO+NT CACT 163 PAR(ISUB)=DECPAR CACT 164 C CACT 165 430 NT=NT+1 CACT 166 IF(NT.EQ.NTOKNS) XTYPE=1 CACT 167 GO TO 1000 CACT 168 C CACT 169 510 TABNO = IPARAM(1) CACT 170 MA = MTAB(TABNO) CACT 171 MB = MTAB(TABNO+1) - 1 CACT 172 C CACT 173 LMISSW = 0 CACT 174 DO 700 M=MA,MB CACT 175 IF( ICODE.NE.TAB(1,M) ) GO TO 700 CACT 176 LMISSW = 1 CACT 177 600 XTYPE=1 CACT 178 NT=0 CACT 179 IPARAM(1) = 1 CACT 180 PARNO=TAB(3,M) CACT 181 LTEMP=TAB(2,M) CACT 182 NTOKNS=TAB(4,M) CACT 183 IF(LTEMP.GT.2) NTOKNS=0 CACT 184 GO TO(610,620,630,640,650,660), LTEMP CACT 185 C CACT 186 C NAME OF INTEGER PARAMETER CACT 187 610 XTYPE=2 CACT 188 GO TO 700 CACT 189 C CACT 190 C NAME OF DECIMAL PARAMETER CACT 191 620 XTYPE=3 CACT 192 GO TO 700 CACT 193 C CACT 194 C IMPLICITLY SPECIFIED INTEGER PARAMETER CACT 195 C A SINGLE IMPLICITLY SPECIFIED PARAMTER CAN IN PRINCIPLE HOLD CACT 196 C AS MANY AS THREE PARAMTERS WHICH THE MAIN PROGRAM CONSIDERS CACT 197 C CONCEPTUALLY DISTINCT. ONE IS IN THE UNITS POSITION, ANOTHER IN CACT 198 C IN THE HUNDREDS POSITION, AND ANOTHER IN THE TEN-THOUSANDS CACT 199 C POSITION. CACT 200 630 LP = 1 CACT 201 LT=TAB(4,M) CACT 202 IF(LT.GE.100) LP = 100 CACT 203 IF(LT.GE.10000) LP = 10000 CACT 204 LQ = 100*LP CACT 205 LA = LPAR(PARNO) CACT 206 LPAR(PARNO)= LA-(MOD(LA,LQ)/LP)*LP +TAB(4,M) CACT 207 GO TO 700 CACT 208 C CACT 209 C IMPLICITLY SPECIFIED DECIMAL PARAMETER CACT 210 640 PAR(PARNO)=ATAB(4,M) CACT 211 GO TO 700 CACT 212 C CACT 213 C INTERNAL PARAMETER OF CCACT PROGRAM CACT 214 650 IPARAM(PARNO)=TAB(4,M) CACT 215 GO TO 700 CACT 216 C CACT 217 C IMPLICITLY SPECIFIED INTEGER PARAMETER CACT 218 C DETERMINES PROGRAM FLOW UPON RETURN FROM CCACT--SO ONLY ONE OF CACT 219 C THIS TYPE PER CONTROL CARD ALLOWED. CACT 220 660 IF(NOPCC)9999,665,996 CACT 221 665 NOPCC=1 CACT 222 LPAR(PARNO)=TAB(4,M) CACT 223 GO TO 700 CACT 224 C CACT 225 700 CONTINUE CACT 226 IF(LMISSW.EQ.0) GO TO 994 CACT 227 1000 CONTINUE CACT 228 IF(NT.LT.NTOKNS) GO TO 997 CACT 229 C CACT 230 1001 RETURN CACT 231 C CACT 232 994 WRITE(LPRINT,9) CACT 233 WRITE(LPRINT,998) (WORD(L),L=1,18) CACT 234 998 FORMAT(27H UNDEFINED CONTROL PHRASE ,18A1) CACT 235 GO TO 999 CACT 236 995 WRITE (LPRINT,9) CACT 237 IF(TYPE .LT.4) READ(LSCRAT,3) (WORD(L),L=1,18) CACT 238 IF (TYPE.NE.4) GO TO 9000 CACT 239 WORD(1)=DOLLAR CACT 240 DO 900 L=2,18 CACT 241 WORD(L)=BLANK CACT 242 900 CONTINUE CACT 243 9000 WRITE (LPRINT,99) K,XTYPE,TYPE, (WORD(L),L=1,18) CACT 244 WRITE(LPRINT,96) CACT 245 96 FORMAT(34H IS THIS ITEM SPELLED INCORRECTLY ) CACT 246 GO TO 999 CACT 247 996 WRITE(LPRINT,9) CACT 248 WRITE(LPRINT,6) CACT 249 6 FORMAT(69H0NO MORE THAN ONE CONTROL PHRASE OF TYPE 6 ALLOWED ON A CACT 250 .CONTROL CARD.) CACT 251 GO TO 999 CACT 252 997 WRITE(LPRINT,9) CACT 253 WRITE(LPRINT,97) CACT 254 97 FORMAT(47H CONTROL PHRASE NOT COMPLETED ON A SINGLE CARD.) CACT 255 999 STOP CACT 256 C CACT 257 9999 WRITE(LPRINT,7) CACT 258 7 FORMAT(23H0ERROR EXIT FROM CCACT.) CACT 259 GO TO 999 CACT 260 END CACT 261 CCONFIG CONFIG 1 SUBROUTINE CONFIG(N,ND) CONFIG 2 C CONFIG FOR KYST JANUARY,1973 CONFIG 3 C WRITTEN BY F.YOUNG CONFIG 4 C MODIFIED FOR KYST BY J.KRUSKAL AND J.SEERY JANUARY,1973 CONFIG 5 C KYST VERSION 2 OCTOBER, 1976 CONFIG 6 C ROUTINE TO OBTAIN A METRIC MULTIDIMENSIONAL SCALING SOLUTION CONFIG 7 C CONFIG 8 DIMENSION DATA(1800),WW(1800),IJ(1800),X(100,6) CONFIG 9 DIMENSION STORE(494),RWMEAN(100),ROOTS(6),XBEST(100,6), CONFIG10 . DATAIN(5430) CONFIG11 REAL MEAN CONFIG12 COMMON /KYST1/ DATA,WW,X,IJ CONFIG13 COMMON /KYST2/ STORE,RWMEAN,ROOTS,XBEST,DATAIN CONFIG14 C CONFIG15 C PARAMETERS CONFIG16 C X - CONFIGURATION AT OUTPUT CONFIG17 C N - NUMBER OF POINTS IN SCALING CONFIG18 C ND - NUMBER OF DIMENSIONS OF SOLUTION CONFIG19 C DATAIN - ON INPUT, DATA FROM INICON. DESTROYED BY CONFIG CONFIG20 C RWMEAN - ROW AVERAGES CONFIG21 C ROOTS - ON RETURN FROM SGEV, THE FIRST ND EIGENVALUES CONFIG22 C CONFIG23 ISUB(I1,I2)=((I2-1)*(I2-2))/2+I1 CONFIG24 MEAN=0.0 CONFIG25 FN=N CONFIG26 C CONFIG27 C COMPUTE ROW,COLUMN, AND OVERALL AVERAGES OF DATIN VALUES SQUARED CONFIG28 C CONFIG29 DO 130 J=1,N