1 : C CGRES SOURCE CHAT 05/01/12 21:53:35 5004 2 : SUBROUTINE CGRES(KSTO, 3 : 1 NL,ILG, 4 : 2 IMAT,IA,JA,KA, 5 : 3 IZB,IZP,ICOLD, 6 : 4 NPT,NPITE,NEFF,ICONV,EPI,IPOU,VPOU, 7 : 5 NIMPR,IPAT) 8 : C 9 : CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 10 : C C 11 : C RESOLUTION ITERATIVE D'UN SYSTEME C 12 : C SYMETRIQUE CREUX PAR LA METHODE DU C 13 : C GRADIENT CONJUGUE SANS PRECONDITIONNEMENT C 14 : C C 15 : C***** INFOS SUR LA METHODE (LIGNE 0) C 16 : C C 17 : C DEUX MODES DE STOCKAGE SONT PROPOSES C 18 : C MORSE OU COMPRESSE. C 19 : C C 20 : C MORSE (KSTO=0) C 21 : C COMPRESSE (KSTO=1) C 22 : C C 23 : C POUR LE CAS SANS PRECONDITIONNEMENT ON RESOUD C 24 : C C 25 : C A*xi = b C 26 : C C 27 : C***** INFOS SUR LES DIMENSIONS (LIGNE 1) C 28 : C C 29 : C NL (L'ORDRE DU SYSTEME) C 30 : C C 31 : C ILG (LONGUEUR LA POUR LE MORSE) C 32 : C (NOMBRE MAXI DE VALEURS PAR LIGNE DE LA C 33 : C MATRICE POUR LE COMPRESSE) C 34 : C C 35 : C***** INFOS SUR LA MATRICE (LIGNE 2) C 36 : C C 37 : C IMAT (LA MATRICE DANS UN LISTREEL) C 38 : C C 39 : C IA (TABLEAU DE POINTEURS DU STOCKAGE MORSE) C 40 : C C 41 : C IL DOIT ETRE DIMENSIONNE A NL+1 AU MINIMUM C 42 : C C 43 : C JA (TABLEAU DE CONNECTIVITES DE LA MATRICE) C 44 : C C 45 : C IL EST STOCKE EN MORSE ET DE LONGUEUR LA C 46 : C C 47 : C KA (TABLEAU DE CONNECTIVITES EN MODE COMPRESSE) C 48 : C C 49 : C IL EST DE LA FORME KA(NL,NNZ) AVEC NNZ C 50 : C LE NOMBRE MAXI DE TERMES NON NULS POUR UNE C 51 : C LIGNE DE LA MATRICE C 52 : C C 53 : C NB : EN MORSE, KA N'EST PAS UTILISE C 54 : C EN COMPRESSE, IA ET JA NE SONT PAS UTILISES C 55 : C C 56 : C IA,JA ET KA SONT DES LISTENTI. C 57 : C C 58 : C***** INFOS SUR LE SECOND MEMBRE, LA SOLUTION INITIALE C 59 : C***** ET LA SOLUTION FINALE. (LIGNE 3) C 60 : C C 61 : C IZB SECOND MEMBRE EN ENTREE C 62 : C PAS MODIFIE. C 63 : C C 64 : C IZP ESTIMATION DE LA SOLUTION EN ENTREE C 65 : C LA SOLUTION EN SORTIE, IL FAUT LA C 66 : C CONSERVER POUR LA RESOLUTION SUIVANTE C 67 : C C 68 : C IZB ET IZP SONT DES LISTREEL DIMENSIONNES A NL C 69 : C C 70 : C ICOLD =1 (DEMARRAGE FROID) C 71 : C ICOLD!=1 (DEMARRAGE CHAUD A PARTIR DE LA C 72 : C SOLUTION INITIALE FOURNIE DANS IZP) C 73 : C C 74 : C***** INFOS SUR LA CONVERGENCE (LIGNE 4) C 75 : C C 76 : C NPT NOMBRE MAXI D'ITERATIONS EN ENTREE C 77 : C SI IL EST < 10, IL EST MIS A 10 C 78 : C NPITE FREQUENCE DES TESTS C 79 : C NEFF NOMBRE REEL EN SORTIE C 80 : C ICONV 0 NON CONVERGE, 1 CONVERGE C 81 : C EPI PRECISION DEMANDEE C 82 : C IPOU=1 LA PRECISION VA ETRE POUSSE D'UN FACTEUR C 83 : C VPOU. CECI PEUT S'AVERER UTILE LORS D'UN C 84 : C CALCUL ITERATIF OU ON VEUT DES LE DEBUT C 85 : C UNE SOLUTION PRECISE POUR EVITER DE PARTIR C 86 : C SUR UNE MAUVAISE PISTE. VPOU A UNE VALEUR C 87 : C PAR DEFAUT DE 100. C 88 : C C 89 : C***** INFOS SUR LE CONTROLE DES IMPRESSIONS (LIGNE 5) C 90 : C C 91 : C NIMPR =1 (INFOS SUR LA CONVERGENCE) C 92 : C NIMPR!=1 (SILENCIEUX SAUF EN CAS DE DIVERGENCE) C 93 : C C 94 : C IPAT = (FACULTATIF, PAR EXEMPLE UN PAS DE TEMPS)C 95 : C C 96 : C RSETD(B,A) B <- A C 97 : C ADIVEC(A,B,C,alpha) C <- A + alpha*B C 98 : C PIMDV(A,B,C) C <- A*B C 99 : C C 100 : CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 101 : C 102 : IMPLICIT INTEGER(I-N) 103 : IMPLICIT REAL*8 (A-H,O-Z) 104 : -INC SMLENTI 105 : POINTEUR IA.MLENTI 106 : POINTEUR JA.MLENTI 107 : POINTEUR KA.MLENTI 108 : -INC SMLREEL 109 : POINTEUR IMAT.MLREEL 110 : POINTEUR IZB.MLREEL 111 : POINTEUR IZBB.MLREEL 112 : POINTEUR IZP.MLREEL 113 : POINTEUR IZPS.MLREEL 114 : POINTEUR IZZPR.MLREEL 115 : POINTEUR IZPMV.MLREEL 116 : POINTEUR IZTMP.MLREEL 117 : POINTEUR IZDIR.MLREEL 118 : 119 : DATA VPOU0/100./ 120 : 121 : IF(KSTO.EQ.0) LA=ILG 122 : IF(KSTO.EQ.1) NNZ=ILG 123 : IF(IPOU.EQ.1) THEN 124 : VPOUE=VPOU0 125 : IF(VPOU.NE.VPOU0) VPOUE=VPOU 126 : ENDIF 127 : 127 : CALL OOOMRU(1) 128 : 129 : C On sauve le second membre en le dupliquant dans IZBB 130 : C qui lui sera ecrase. 131 : SEGACT IZB 132 : JG=IZB.PROG(/1) 133 : SEGINI IZBB 133 : CALL RSETD(IZBB.PROG,IZB.PROG,JG) 134 : SEGDES IZB 135 : 136 : SEGACT IZP 137 : 138 : JG=NL 139 : SEGINI IZPS 140 : SEGINI IZZPR 141 : SEGINI IZPMV 142 : SEGINI IZTMP 143 : SEGINI IZDIR 144 : 145 : SEGACT IMAT 146 : IF(KSTO.EQ.0) SEGACT IA 147 : IF(KSTO.EQ.0) SEGACT JA 148 : IF(KSTO.EQ.1) SEGACT KA 149 : 150 : IF(NPT.LT.10) NPT=10 151 : IF(NPT.GT.NL) NPT=NL 152 : 152 : IF(ICOLD.EQ.1) CALL INITD (IZP.PROG,NL,0.D0) 152 : CALL RSETD(IZPS.PROG,IZP.PROG,NL) 152 : CALL RSETD(IZTMP.PROG,IZP.PROG,NL) 153 : IF(KSTO.EQ.0) THEN 153 : C CALL PMVM(A,IA,JA,X,Y,NL,LA) 153 : CALL PMVM(IMAT.PROG,IA.LECT,JA.LECT, 154 : 1 IZTMP.PROG,IZPMV.PROG,NL,LA) 155 : ELSE 155 : C CALL PMVC(A,KA,X,Y,NL,NNZ) 155 : CALL PMVC(IMAT.PROG,KA.LECT, 156 : 1 IZTMP.PROG,IZPMV.PROG,NL,NNZ) 157 : ENDIF 157 : CALL ADIVEC (IZBB.PROG,IZPMV.PROG,IZBB.PROG,-1.d0,NL) 158 : ANORB=SQRT(PISCAL(IZBB.PROG,IZBB.PROG,NL)) 158 : CALL RSETD(IZZPR.PROG,IZBB.PROG,NL) 158 : CALL RSETD(IZDIR.PROG,IZZPR.PROG,NL) 159 : SCAR1=PISCAL(IZBB.PROG,IZDIR.PROG,NL) 160 : C 161 : C===========DEBUT DE LA BOUCLE ITERATIVE=====600======== 162 : C 163 : NEFF=0 164 : NCOMP=0 165 : ICONV=0 166 : DO 600 I=1,NPT 167 : NEFF=NEFF+1 168 : NCOMP=NCOMP+1 168 : CALL RSETD(IZTMP.PROG,IZDIR.PROG,NL) 169 : IF(KSTO.EQ.0) THEN 169 : C CALL PMVM(A,IA,JA,X,Y,NL,LA) 169 : CALL PMVM(IMAT.PROG,IA.LECT,JA.LECT, 170 : 1 IZTMP.PROG,IZPMV.PROG,NL,LA) 171 : ELSE 171 : C CALL PMVC(A,KA,X,Y,NL,NNZ) 171 : CALL PMVC(IMAT.PROG,KA.LECT, 172 : 1 IZTMP.PROG,IZPMV.PROG,NL,NNZ) 173 : ENDIF 174 : SCADIR=PISCAL(IZDIR.PROG,IZPMV.PROG,NL) 175 : ALPHA=SCAR1/SCADIR 175 : CALL ADIVEC(IZP.PROG,IZDIR.PROG,IZP.PROG,ALPHA,NL) 176 : ALP=-1.D0*ALPHA 176 : CALL ADIVEC (IZBB.PROG,IZPMV.PROG,IZBB.PROG,ALP,NL) 176 : CALL RSETD(IZZPR.PROG,IZBB.PROG,NL) 177 : SCAR2=PISCAL(IZBB.PROG,IZZPR.PROG,NL) 178 : BETA=SCAR2/SCAR1 178 : CALL ADIVEC (IZZPR.PROG,IZDIR.PROG,IZDIR.PROG,BETA,NL) 179 : IF (NCOMP.EQ.NPITE.OR.I.EQ.NPT) THEN 180 : INDIC=0 181 : IF (I.EQ.NPT) INDIC=1 182 : C On pousse la precision si necessaire 183 : ZEPS=EPI 184 : IF(IPOU.EQ.1) ZEPS=EPI/VPOUE 184 : CALL COTEST (IZP.PROG,IZPS.PROG, 185 : 1 ZEPS,DELTA,DELTAP, 186 : 2 INDIC,ICONV,NL) 187 : SEPA=DELTA 188 : RSEPA=DELTA/DELTAP 188 : CALL RSETD(IZPS.PROG,IZP.PROG,NL) 189 : NCOMP=0 190 : ENDIF 191 : IF (ICONV.EQ.1) GOTO 601 192 : SCAR1=SCAR2 193 : 600 CONTINUE 194 : 601 CONTINUE 194 : CALL MTABL(IZBB.PROG,RES,IMA,NL) 195 : RRES=RES/ANORB 196 : IF (ICONV.NE.1) THEN 197 : WRITE(6,1111)IPAT,NPT,RSEPA,DELTAP,RES,IZBB.PROG(IMA) 198 : ELSE 199 : IF(NIMPR.EQ.1) WRITE(6,2222)IPAT,EPI, 200 : 1 RES,SEPA,I, 201 : 2 RRES,RSEPA,IMA 202 : ENDIF 203 : 204 : SEGDES IZP 205 : 206 : SEGSUP IZBB 207 : SEGSUP IZPS 208 : SEGSUP IZZPR 209 : SEGSUP IZDIR 210 : SEGSUP IZTMP 211 : SEGSUP IZPMV 212 : 213 : SEGDES IMAT*(NOMOD,MRU) 214 : IF(KSTO.EQ.0) SEGDES IA*(NOMOD,MRU) 215 : IF(KSTO.EQ.0) SEGDES JA*(NOMOD,MRU) 216 : IF(KSTO.EQ.1) SEGDES KA*(NOMOD,MRU) 217 : 217 : CALL OOOMRU(0) 218 : 219 : C 220 : 1111 FORMAT(1X,'* CGRES* ',I6,' NON CONV EN ',I5, 221 : 1' ITERA, RSEPA = ',D8.2,1X,'DELTAP = ',D8.2,1X, 222 : 2' RES = ',D8.2,' P = ',D8.2) 223 : 2222 FORMAT(1X,'* CGRES* ',I6,' *PREC ',D8.2, 224 : 1 ' * RESIDU-SEPARA @ NITE ',D8.2,1X,D8.2,1X,I5, 225 : 2 ' * RRES-RSEP ',D8.2,1X,D8.2, 226 : 3 ' * ELEMENT ',I7) 227 : C 228 : C FIN RESOLUTION 229 : C 230 : RETURN 231 : END 232 : 233 :
© Cast3M 2003 - All rights reserved.
Disclaimer