CTITLEMFSD -- SQUARE ROOT MATRIX OF POS DEF SYM MATRIX 00000010 C*********************************************************************** C COPYRIGHT ATLANTIC RICHFIELD COMPANY 1991 * C*********************************************************************** CA AUTHOR IBM SSP 00000020 CA DESIGNER IBM SSP 00000030 CA LANGUAGE FORTRAN 77 00000040 CA SYSTEM IBM & CRAY 00000041 CA WRITTEN IBM 00000050 C REVISED MO-DA-YR 00000060 C REVISED 05-20-85 TWH. ADAPT TO RUN ON IBM & 00000061 C CRAY SYSTEMS. 00000062 C 00000070 CA 00000080 CA 00000090 CA CALL MFSD (A, N, EPS, IER) 00000100 CA 00000110 CA 00000120 CA FOR PARAMETER DEFINITION, SEE IBM SSP MANUAL 00000130 CA 00000140 CA 00000150 CA THIS ROUTINE FACTORS A POSITIVE DEFINITE SYMMETRIC 00000160 CA MATRIX INTO ITS UPPER TRIANGULAR SQUARE ROOT MATRIX. 00000170 CA 00000180 CA 00000190 C 00000200 C SUBROUTINES CALLED: ABS 00000210 C DBLE 00000220 C DSQRT 00000230 C 00000240 C EJECT 00000250 C 00000260 SUBROUTINE MFSD (A, N, EPS, IER) 00000270 CAEND 00000280 C 00000290 DIMENSION A(1) 00000300 DOUBLE PRECISION DPIV,DSUM 00000310 IF(N-1) 120 , 10 , 10 00000320 C 00000330 10 IER=0 00000340 KPIV=0 00000350 C 00000360 DO 110 00000370 * K=1,N 00000380 KPIV=KPIV+K 00000390 IND=KPIV 00000400 LEND=K-1 00000410 TOL=ABS(EPS*A(KPIV)) 00000420 C 00000430 DO 110 00000440 * I=K,N 00000450 DSUM=0.D0 00000460 IF(LEND) 20 , 40 , 20 00000470 C 00000480 20 DO 30 00000490 * L=1,LEND 00000500 LANF=KPIV-L 00000510 LIND=IND-L 00000520 C 00000530 30 DSUM=DSUM+DBLE(A(LANF)*A(LIND)) 00000540 C 00000550 40 DSUM=DBLE(A(IND))-DSUM 00000560 IF(I-K) 100 , 50 , 100 00000570 C 00000580 50 IF(SNGL(DSUM)-TOL) 60 , 60 , 90 00000590 C 00000600 60 IF(DSUM) 120 , 120 , 70 00000610 C 00000620 70 IF(IER) 80 , 80 , 90 00000630 C 00000640 80 IER=K-1 00000650 C 00000660 90 DPIV=DSQRT(DSUM) 00000670 A(KPIV)=DPIV 00000680 DPIV=1.D0/DPIV 00000690 GO TO 110 00000700 C 00000710 100 A(IND)=DSUM*DPIV 00000720 C 00000730 110 IND=IND+I 00000740 RETURN 00000750 C 00000760 120 IER=-1 00000770 RETURN 00000780 END 00000790