?? dsrc2c.f90
字號:
! 90 IB1 = 1 IB2 = IB1+N IB3 = IB2+N IB4 = IB3+N IB5 = IB4+N IB6 = IB5+N IB7 = IB6+N IPARM(8) = 6*N+2*ITMAX IF (ISYM.NE.0) IPARM(8) = IPARM(8)+2*ITMAX IF (NW.GE.IPARM(8)) GO TO 110 IER = 42 IF (LEVEL.GE.0) WRITE (NOUT,100) NW,IPARM(8) 100 FORMAT ('0','*** F A T A L E R R O R ************'/'0', & & ' CALLED FROM ITPACK ROUTINE SSORCG '/' ', & & ' NOT ENOUGH WORKSPACE AT ',I10/' ',' SET IPARM(8) =',I10& & ,' (NW)') GO TO 390 ! ! ... PERMUTE TO RED-BLACK SYSTEM IF REQUESTED ! 110 NB = IPARM(9) IF (NB.LT.0) GO TO 170 N3 = 3*N CALL IVFILL (N3,IWKSP,0) CALL PRBNDX (N,NB,IA,JA,IWKSP,IWKSP(IB2),LEVEL,NOUT,IER) IF (IER.EQ.0) GO TO 130 IF (LEVEL.GE.0) WRITE (NOUT,120) IER,NB 120 FORMAT ('0','*** F A T A L E R R O R ************'/'0', & & ' CALLED FROM ITPACK ROUTINE SSORCG '/' ', & & ' ERROR DETECTED IN SUBROUTINE PRBNDX'/' ', & & ' WHICH COMPUTES THE RED-BLACK INDEXING'/' ',' IER = ',I5& & ,' IPARM(9) = ',I5,' (NB)') GO TO 390 ! ! ... PERMUTE MATRIX AND RHS! 130 IF (LEVEL.GE.2) WRITE (NOUT,140) NB 140 FORMAT (/10X,'ORDER OF BLACK SUBSYSTEM = ',I5,' (NB)') CALL PERMAT (N,IA,JA,A,IWKSP,IWKSP(IB3),ISYM,LEVEL,NOUT,IER) IF (IER.EQ.0) GO TO 160 IF (LEVEL.GE.0) WRITE (NOUT,150) IER 150 FORMAT ('0','*** F A T A L E R R O R ************'/'0', & & ' CALLED FROM ITPACK ROUTINE SSORCG '/' ', & & ' ERROR DETECTED IN SUBROUTINE PERMAT'/' ', & & ' WHICH DOES THE RED-BLACK PERMUTATION'/' ',' IER = ',I5) GO TO 390 160 CALL PERVEC (N,RHS,IWKSP) CALL PERVEC (N,U,IWKSP) ! ! ... SCALE LINEAR SYSTEM, U, AND RHS BY THE SQUARE ROOT OF THE ! ... DIAGONAL ELEMENTS. ! 170 CONTINUE CALL VFILL (IPARM(8),WKSP,0.0D0) CALL SCAL (N,IA,JA,A,RHS,U,WKSP,LEVEL,NOUT,IER) IF (IER.EQ.0) GO TO 190 IF (LEVEL.GE.0) WRITE (NOUT,180) IER 180 FORMAT ('0','*** F A T A L E R R O R ************'/'0', & & ' CALLED FROM ITPACK ROUTINE SSORCG '/' ', & & ' ERROR DETECTED IN SUBROUTINE SCAL '/' ', & & ' WHICH SCALES THE SYSTEM '/' ',' IER = ',I5) GO TO 390 190 IF (LEVEL.LE.2) GO TO 220 WRITE (NOUT,200) 200 FORMAT (///1X,'IN THE FOLLOWING, RHO AND GAMMA ARE',& & ' ACCELERATION PARAMETERS') WRITE (NOUT,210) 210 FORMAT (1X,'S-PRIME IS AN INITIAL ESTIMATE FOR NEW CME') 220 IF (IPARM(11).NE.0) GO TO 230 TIMI1 = TIMER(DUMMY) ! ! ... SPECIAL PROCEDURE FOR FULLY ADAPTIVE CASE.! 230 CONTINUE IF (.NOT.ADAPT) GO TO 250 IF (.NOT.BETADT) GO TO 240 CALL VFILL (N,WKSP(IB1),1.D0) BETNEW = PBETA(N,IA,JA,A,WKSP(IB1),WKSP(IB2),WKSP(IB3))/ & & DBLE(FLOAT(N)) BETAB = DMAX1(BETAB,.25D0,BETNEW) 240 CALL OMEG (0.D0,1) IS = 0 ! ! ... INITIALIZE FORWARD PSEUDO-RESIDUAL! 250 CALL DCOPY (N,RHS,1,WKSP(IB1),1) CALL DCOPY (N,U,1,WKSP(IB2),1) CALL PFSOR (N,IA,JA,A,WKSP(IB2),WKSP(IB1)) CALL VEVMW (N,WKSP(IB2),U) ! ! ... ITERATION SEQUENCE ! ITMAX1 = ITMAX+1 DO 270 LOOP = 1,ITMAX1 IN = LOOP-1 IF (MOD(IN,2).EQ.1) GO TO 260! ! ... CODE FOR THE EVEN ITERATIONS. ! ! U = U(IN) WKSP(IB2) = C(IN) ! WKSP(IB1) = U(IN-1) WKSP(IB3) = C(IN-1) ! CALL ITSRCG (N,IA,JA,A,RHS,U,WKSP(IB1),WKSP(IB2),WKSP(IB3),& & WKSP(IB4),WKSP(IB5),WKSP(IB6),WKSP(IB7)) ! IF (HALT) GO TO 300 GO TO 270! ! ... CODE FOR THE ODD ITERATIONS. ! ! U = U(IN-1) WKSP(IB2) = C(IN-1) ! WKSP(IB1) = U(IN) WKSP(IB3) =C(IN)! 260 CALL ITSRCG (N,IA,JA,A,RHS,WKSP(IB1),U,WKSP(IB3),WKSP(IB2),& & WKSP(IB4),WKSP(IB5),WKSP(IB6),WKSP(IB7)) ! IF (HALT) GO TO 300 270 CONTINUE ! ! ... ITMAX HAS BEEN REACHED! IF (IPARM(11).NE.0) GO TO 280 TIMI2 = TIMER(DUMMY) TIME1 = DBLE(TIMI2-TIMI1) 280 IF (LEVEL.GE.1) WRITE (NOUT,290) ITMAX 290 FORMAT ('0','*** W A R N I N G ************'/'0', & & ' IN ITPACK ROUTINE SSORCG'/' ',' FAILURE TO CONVERGE IN'& & ,I5,' ITERATIONS') IER = 43 IF (IPARM(3).EQ.0) RPARM(1) = STPTST GO TO 330 ! ! ... METHOD HAS CONVERGED ! 300 IF (IPARM(11).NE.0) GO TO 310 TIMI2 = TIMER(DUMMY) TIME1 = DBLE(TIMI2-TIMI1) 310 IF (LEVEL.GE.1) WRITE (NOUT,320) IN 320 FORMAT (/1X,'SSORCG HAS CONVERGED IN ',I5,' ITERATIONS') ! ! ... PUT SOLUTION INTO U IF NOT ALREADY THERE. ! 330 CONTINUE IF (MOD(IN,2).EQ.1) CALL DCOPY (N,WKSP(IB1),1,U,1) ! ! ... UNSCALE THE MATRIX, SOLUTION, AND RHS VECTORS. ! CALL UNSCAL (N,IA,JA,A,RHS,U,WKSP)! ! ... UN-PERMUTE MATRIX,RHS, AND SOLUTION ! IF (IPARM(9).LT.0) GO TO 360 CALL PERMAT (N,IA,JA,A,IWKSP(IB2),IWKSP(IB3),ISYM,LEVEL,NOUT, & & IERPER) IF (IERPER.EQ.0) GO TO 350 IF (LEVEL.GE.0) WRITE (NOUT,340) IERPER 340 FORMAT ('0','*** F A T A L E R R O R ************'/'0', & & ' CALLED FROM ITPACK ROUTINE SSORCG '/' ', & & ' ERROR DETECTED IN SUBROUTINE PERMAT'/' ', & & ' WHICH UNDOES THE RED-BLACK PERMUTATION '/' ', & & ' IER = ',I5) IF (IER.EQ.0) IER = IERPER GO TO 390 350 CALL PERVEC (N,RHS,IWKSP(IB2)) CALL PERVEC (N,U,IWKSP(IB2)) ! ! ... OPTIONAL ERROR ANALYSIS ! 360 IDGTS = IPARM(12) IF (IDGTS.LT.0) GO TO 370 IF (IPARM(2).LE.0) IDGTS = 0 CALL PERROR5 (N,IA,JA,A,RHS,U,WKSP,DIGIT1,DIGIT2,IDGTS)! ! ... SET RETURN PARAMETERS IN IPARM AND RPARM ! 370 IF (IPARM(11).NE.0) GO TO 380 TIMJ2 = TIMER(DUMMY) TIME2 = DBLE(TIMJ2-TIMJ1) 380 IPARM(8) = IPARM(8)-2*(ITMAX-IN) IF (ISYM.NE.0) IPARM(8) = IPARM(8)-2*(ITMAX-IN) IF (IPARM(3).NE.0) GO TO 390 IPARM(1) = IN IPARM(9) = NB RPARM(2) = CME RPARM(3) = SME RPARM(5) = OMEGA RPARM(6) = SPECR RPARM(7) = BETAB RPARM(9) = TIME1 RPARM(10) = TIME2 RPARM(11) = DIGIT1 RPARM(12) = DIGIT2 ! 390 CONTINUE IERR = IER IF (LEVEL.GE.3) CALL ECHALL (N,IA,JA,A,RHS,IPARM,RPARM,2) ! RETURN END SUBROUTINE SSORSI (NN,IA,JA,A,RHS,U,IWKSP,NW,WKSP,IPARM,RPARM,& & IERR) ! ! ITPACK 2C MAIN SUBROUTINE SSORSI (SYMMETRIC SUCCESSIVE RELAX- ! ATION SEMI-ITERATION) ! EACH OF THE MAIN SUBROUTINES: ! JCG, JSI, SOR, SSORCG, SSORSI, RSCG, RSSI ! CAN BE USED INDEPENDENTLY OF THE OTHERS ! ! ... FUNCTION: ! ! THIS SUBROUTINE, SSORSI, DRIVES THE SYMMETRIC SOR-SI ! ALGORITHM. ! ! ... PARAMETER LIST: ! ! N INPUT INTEGER. DIMENSION OF THE MATRIX. (= NN) ! IA,JA INPUT INTEGER VECTORS. THE TWO INTEGER ARRAYS OF ! THE SPARSE MATRIX REPRESENTATION. ! A INPUT D.P. VECTOR. THE D.P. ARRAY OF THE SPARSE ! MATRIX REPRESENTATION.! RHS INPUT D.P. VECTOR. CONTAINS THE RIGHT HAND SIDE ! OF THE MATRIX PROBLEM.! U INPUT/OUTPUT D.P. VECTOR. ON INPUT, U CONTAINS THE ! INITIAL GUESS TO THE SOLUTION. ON OUTPUT, IT CONTAINS ! THE LATEST ESTIMATE TO THE SOLUTION. ! IWKSP INTEGER VECTOR WORKSPACE OF LENGTH 3*N ! NW INPUT INTEGER. LENGTH OF AVAILABLE WKSP. ON OUTPUT, ! IPARM(8) IS AMOUNT USED. ! WKSP D.P. VECTOR USED FOR WORKING SPACE. SSORSI ! NEEDS THIS TO BE IN LENGTH AT LEAST 5*N! IPARM INTEGER VECTOR OF LENGTH 12. ALLOWS USER TO SPECIFY! SOME INTEGER PARAMETERS WHICH AFFECT THE METHOD. IF! RPARM D.P. VECTOR OF LENGTH 12. ALLOWS USER TO SPECIFY SOME ! D.P. PARAMETERS WHICH AFFECT THE METHOD.! IER OUTPUT INTEGER. ERROR FLAG. (= IERR) ! ! ... SSORSI SUBPROGRAM REFERENCES: ! ! FROM ITPACK BISRCH, CHEBY, CHGSI, DFAULT, ECHALL, ECHOUT, ! ITERM, TIMER, ITSRSI, IVFILL, OMEG, ! OMGSTR, PARSI, PBETA, PERMAT, PERROR5, ! PERVEC, PFSOR, PMULT, PRBNDX, PSSOR1, ! PSTOP, PVTBV, QSORT, SBELM, SCAL, DCOPY, ! DDOT, SUM3, TSTCHG, UNSCAL, VEVPW, VFILL, ! VOUT, WEVMW ! SYSTEM DABS, DLOG, DLOG10, DBLE(AMAX0), DMAX1, DBLE(F! MOD, DSQRT ! ! VERSION: ITPACK 2C (MARCH 1982)! ! CODE WRITTEN BY: DAVID KINCAID, ROGER GRIMES, JOHN RESPESS ! CENTER FOR NUMERICAL ANALYSIS ! UNIVERSITY OF TEXAS ! AUSTIN, TX 78712 ! (512) 471-1242! ! FOR ADDITIONAL DETAILS ON THE ! (A) SUBROUTINE SEE TOMS ARTICLE 1982 ! (B) ALGORITHM SEE CNA REPORT 150 ! ! BASED ON THEORY BY: DAVID YOUNG, DAVID KINCAID, LOU HAGEMAN ! ! REFERENCE THE BOOK: APPLIED ITERATIVE METHODS ! L. HAGEMAN, D. YOUNG ! ACADEMIC PRESS, 1981 ! ! ************************************************** ! * IMPORTANT NOTE * ! * * ! * WHEN INSTALLING ITPACK ROUTINES ON A * ! * DIFFERENT COMPUTER, RESET SOME OF THE VALUES * ! * IN SUBROUTNE DFAULT. MOST IMPORTANT ARE * ! * * ! * DRELPR MACHINE RELATIVE PRECISION * ! * RPARM(1) STOPPING CRITERION * ! * * ! * ALSO CHANGE SYSTEM-DEPENDENT ROUTINE * ! * SECOND USED IN TIMER * ! * * ! ************************************************** ! ! SPECIFICATIONS FOR ARGUMENTS ! INTEGER IA(1),JA(1),IWKSP(1),IPARM(12),NN,NW,IERR DOUBLE PRECISION A(1),RHS(NN),U(NN),WKSP(NW),RPARM(12)! ! SPECIFICATIONS FOR LOCAL VARIABLES! INTEGER IB1,IB2,IB3,IB4,IB5,IDGTS,IER,IERPER,ITMAX1,LOOP,N,NB,N3 DOUBLE PRECISION BETNEW,DIGIT1,DIGIT2,PBETA,TEMP,TIME1,TIME2,TOL! ! *** BEGIN: ITPACK COMMON ! INTEGER IN,IS,ISYM,ITMAX,LEVEL,NOUT COMMON /ITCOM1/ IN,IS,ISYM,ITMAX,LEVEL,NOUT ! LOGICAL ADAPT,BETADT,CASEII,HALT,PARTAD COMMON /ITCOM2/ ADAPT,BETADT,CASEII,HALT,PARTAD ! DOUBLE PRECISION BDELNM,BETAB,CME,DELNNM,DELSNM,FF,GAMMA,OMEGA,QA,& & QT,RHO,RRR,SIGE,SME,SPECR,SPR,DRELPR,STPTST,UDNM,ZETA COMMON /ITCOM3/ BDELNM,BETAB,CME,DELNNM,DELSNM,FF,GAMMA,OMEGA,QA, & & QT,RHO,RRR,SIGE,SME,SPECR,SPR,DRELPR,STPTST,UDNM,ZETA ! ! *** END : ITPACK COMMON ! ! ... VARIABLES IN COMMON BLOCK - ITCOM1! ! IN - ITERATION NUMBER ! ISYM - SYMMETRIC/NONSYMMETRIC STORAGE FORMAT SWITCH ! IS - ITERATION NUMBER WHEN PARAMETERS LAST CHANGED! ITMAX - MAXIMUM NUMBER OF ITERATIONS ALLOWED ! LEVEL - LEVEL OF OUTPUT CONTROL SWITCH ! NOUT - OUTPUT UNIT NUMBER ! ! ... VARIABLES IN COMMON BLOCK - ITCOM2! ! ADAPT - FULLY ADAPTIVE PROCEDURE SWITCH ! BETADT - SWITCH FOR ADAPTIVE DETERMINATION OF BETA ! CASEII - ADAPTIVE PROCEDURE CASE SWITCH ! HALT - STOPPING TEST SWITCH ! PARTAD - PARTIALLY ADAPTIVE PROCEDURE SWITCH! ! ... VARIABLES IN COMMON BLOCK - ITCOM3! ! BDELNM - TWO NORM OF B TIMES DELTA-SUPER-N! BETAB - ESTIMATE FOR THE SPECTRAL RADIUS OF LU MATRIX! CME - ESTIMATE OF LARGEST EIGENVALUE ! DELNNM - INNER PRODUCT OF PSEU
?? 快捷鍵說明
復(fù)制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -