1 : C NTAPE6 SOURCE CHAT 05/01/13 02:02:22 5004 2 : SUBROUTINE NTAPE6(MCP,MCQ,IVMINU,IVMINL,IVMAXU,IVMAXL,IVLAMB, 3 : * M,N,NVD,IVFP,IVFQ,MVDU,MVDL,IVB,IVD,IVN,II,KK,IVDR,IDVD, 4 : * NDR,TERMIN,IVLL,IVUL,IPBASE) 5 : IMPLICIT INTEGER(I-N) 6 : IMPLICIT REAL*8(A-H,O-Z) 7 : LOGICAL TERMIN 8 : -INC TMXMAT 9 : -INC SMLENTI 10 : -INC CCOPTIO 11 : -INC SMLREEL 12 : POINTEUR MLREE4.MLREEL,MLREE5.MLREEL,MLREE6.MLREEL,MLREE7.MLREEL 13 : IVXL1=0 14 : IVXU1=0 15 : IVU1=0 16 : IVN1=0 17 : IVD1=0 18 : IVGM1=0 19 : IVGE1=0 20 : IVLAM1=0 21 : N11 = N + 1 22 : II=-1 23 : KK=-1 24 : MLREEL=IVDR 25 : JG=N11 26 : SEGINI MLREE1,MLREE2 27 : IVPZ=MLREE2 28 : IVQZ=MLREE1 29 : MXMAT=MCQ 29 : CALL MATVE1(XMAT,PROG,M,N11,MLREE1.PROG,1) 30 : MXMAT=MCP 30 : CALL MATVE1(XMAT,PROG,M,N11,MLREE2.PROG,1) 31 : ALPMAX=1.E25 32 : MLREE1=IVLAMB 33 : IPB=0 34 : DO 1 I=1,M 35 : IF(PROG(I).LT.-1.D-20) THEN 36 : IF(ALPMAX+(MLREE1.PROG(I)/PROG(I)).GT.0.D0) THEN 37 : ALPMAX= -( MLREE1.PROG(I)/PROG(I)) 38 : IPB=I 39 : ENDIF 40 : ENDIF 41 : 1 CONTINUE 42 : IF(IIMPI.EQ.1799) 43 : *WRITE(IOIMP,FMT='('' ALPHAMAX IPB = '',E12.5,2X,I3)')ALPMAX,IPB 44 : NDIS=0 45 : MLENTI=IDVD 46 : DO 7 I=1,NVD 47 : NDIS=NDIS+LECT(I)-1 48 : 7 CONTINUE 49 : JG=(NDIS-NVD)+2*(N11-NVD)+2 50 : SEGINI MLREE6,MLENT1,MLENT2 51 : MLREE6.PROG(1)=0.D0 52 : MLENT1.LECT(1)=-3 53 : MLENT2.LECT(1)=0 54 : * 55 : * CALCUL DES ALPHA CRITIQUES 56 : * ON COMMENCE PAR LES VARIABLES CONTINUES 57 : * 58 : IP=1 59 : MLREEL=IVD 60 : MLREE1=IVN 61 : MLREE2=IVMINU 62 : MLREE3=IVMINL 63 : MLREE4=IVQZ 64 : MLREE5=IVPZ 65 : DO 2 I=NVD+1,N11 66 : CN=PROG(I)*(MLREE2.PROG(I)**2) 67 : CN=CN-(MLREE1.PROG(I)*(MLREE3.PROG(I)**2)) 68 : CD=(MLREE3.PROG(I) ** 2) * MLREE5.PROG(I) 69 : CD=CD-((MLREE2.PROG(I)**2)*MLREE4.PROG(I)) 70 : IP=IP+1 71 : IF(CD.EQ.0.D0) CD=1.D-35 72 : MLREE6.PROG(IP)=CN/CD 73 : MLENT1.LECT(IP)=I-NVD 74 : MLENT2.LECT(IP)=-1 75 : 2 CONTINUE 76 : MLREE2=IVMAXU 77 : MLREE3=IVMAXL 78 : DO 3 I=NVD+1,N11 79 : CN=PROG(I)*(MLREE2.PROG(I)**2) 80 : CN=CN-(MLREE1.PROG(I)*(MLREE3.PROG(I)**2)) 81 : CD=(MLREE3.PROG(I) ** 2) * MLREE5.PROG(I) 82 : CD=CD-((MLREE2.PROG(I)**2)*MLREE4.PROG(I)) 83 : IP=IP+1 84 : IF(CD.EQ.0.D0) CD=1.D-35 85 : MLREE6.PROG(IP)=CN/CD 86 : MLENT1.LECT(IP)=I-NVD 87 : MLENT2.LECT(IP)=-1 88 : 3 CONTINUE 89 : * 90 : * POUR LES VARIABLES DISCRETES 91 : * 92 : IF(NVD.NE.0) THEN 93 : IF(IIMPI.EQ.1799) 94 : * WRITE(IOIMP,FMT='('' IL Y A'',I4,'' VARIABLES DISCRETES'')')NVD 95 : MXMAT=MVDU 96 : MXMA1=MVDL 97 : DO 4 I=1,NVD 98 : DO 4 J = 3,LECT(I) 99 : CN=PROG(I)*XMAT(I,J)*XMAT(I,J-1) 100 : CN=CN-(MLREE1.PROG(I)*MXMA1.XMAT(I,J)*MXMA1.XMAT(I,J-1)) 101 : CD=MLREE4.PROG(I)* MXMA1.XMAT(I,J)*MXMA1.XMAT(I,J-1) 102 : CD=CD-(MLREE5.PROG(I)*XMAT(I,J)*XMAT(I,J-1)) 103 : IP=IP+1 104 : IF(CD.EQ.0.D0) CD=1.D-35 105 : MLREE6.PROG(IP)= CN/CD 106 : IF(ABS(MLREE6.PROG(IP)).LT.1D-10)MLREE6.PROG(IP)=0.D0 107 : MLENT1.LECT(IP)=I 108 : MLENT2.LECT(IP)=J 109 : 4 CONTINUE 110 : ENDIF 111 : MAC2=MLENT1 112 : MAC3=MLENT2 113 : SEGSUP MLREE4,MLREE5 114 : * 115 : * ON ENTRE ALPHA MAX DANS LE TABLEAU DES ALPHA 116 : * 117 : IP=IP+1 118 : MLREE6.PROG(IP)=ALPMAX 119 : * 120 : MLREEL=MLREE6 121 : * 122 : * TRI DES ALPHA 123 : * 124 : JG=PROG(/1) 125 : SEGINI MLREE1 126 : SEGINI MLENT1,MLENT2 127 : MORDRE=MLENT1 128 : DO 5 I=1,JG 129 : MLENT1.LECT(I)=I 130 : 5 CONTINUE 130 : CALL TRIFLO(PROG,MLREE1.PROG,MLENT1.LECT,MLENT2.LECT,JG) 131 : IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' LISTE DES ALPHA APRES TRI 132 : * '',/,(1X,5E12.5))')(PROG(I),I=1,IP) 133 : SEGSUP MLREE1,MLENT2 134 : * 135 : * ELIMINATION DES ALPHA NON ADMISSIBLES 136 : * 137 : IDL=0 138 : IFL=IP+2 139 : DO 6 I=1,IP 140 : IF(PROG(I).LT.1.D-20) IDL=I 141 : IF(PROG(I).GT.ALPMAX)THEN 142 : IF(IFL.EQ.IP+2) IFL=I-1 143 : MLENT1.LECT(I)=-2 144 : ENDIF 145 : 6 CONTINUE 146 : IFL=MIN(IFL,IP) 147 : IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' LISTE DES ALPHA APRES TRI 148 : * ET ELIMINATION '',/,(1X,5E12.5))')(PROG(I),I=IDL,IFL) 149 : * 150 : * RECHERCHE DU PAS ALPHA OPTIMUM 151 : * 152 : MLREE1=IVLAMB 153 : MLREE2=IVDR 154 : SEGINI,MLREE3=MLREE1 155 : IVLAM1=MLREE3 156 : * 157 : * CALCUL DU PREMIER POINT 158 : * 159 : I=IDL 160 : INPINF=IDL 161 : QQZ=PROG(I) 162 : DO 8 J=1,MLREE1.PROG(/1) 163 : MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J) 164 : 8 CONTINUE 165 : IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' VALEUR DU LAMBDA POUR ALPHA 166 : * =0 '',/,(1X,5E12.5))')(MLREE3.PROG(I),I=1,M) 166 : CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU, 167 : * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1) 167 : CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1, 168 : *IVGM1,IVLAM1,IPBASE) 169 : XLP1=0.D0 170 : MLREE4=IVGE1 171 : DO 9 J=1,MLREE1.PROG(/1) 172 : XLP1=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP1 173 : 9 CONTINUE 174 : IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' 1ER ALPHA , PROD SCAL :'' 175 : * ,1X,E12.5,'' , '',E12.5)')QQZ,XLP1 176 : MLREE4=IVXU1 177 : MLREE5=IVU1 178 : MLREE6=IVN1 179 : SEGSUP MLREE4,MLREE5,MLREE6 180 : MLREE4=IVD1 181 : MLREE5=IVGM1 182 : MLREE6=IVGE1 183 : MLREE7=IVXL1 184 : SEGSUP MLREE4,MLREE5,MLREE6,MLREE7 185 : IVXL1=0 186 : IVXU1=0 187 : IVU1=0 188 : IVN1=0 189 : IVD1=0 190 : IVGM1=0 191 : IVGE1=0 192 : * 193 : * CALCUL SUR LE DERNIER POINT 194 : * 195 : INPSUP=IFL 196 : I=IFL 197 : QQZ=PROG(I) 198 : DO 10 J=1,MLREE1.PROG(/1) 199 : MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J) 200 : 10 CONTINUE 200 : CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU, 201 : * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1) 201 : CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1, 202 : *IVGM1,IVLAM1,IPBASE) 203 : XLP2=0.D0 204 : MLREE4=IVGE1 205 : DO 11 J=1,MLREE1.PROG(/1) 206 : XLP2=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP2 207 : 11 CONTINUE 208 : IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' DERNIER ALPHA , PROD SCAL 209 : * :'',1X,E12.5,'' , '',E12.5)')QQZ,XLP2 210 : * 211 : * RECHERCHE PAR DICHOTOMIE SUR LES INDICES pour trouver l'interval 212 : * de alpha interessant 213 : * 214 : IF(XLP1 * XLP2.GE.0.D0) THEN 215 : ALPHA=ALPMAX 216 : IF (ALPMAX.LT.1E25) THEN 217 : KK = -3 218 : II=IPB 219 : ENDIF 220 : GO TO 20 221 : ENDIF 222 : 12 CONTINUE 223 : IK=(INPINF+INPSUP)/2 224 : IF(IK.EQ.INPINF) GO TO 15 225 : MLREE4=IVXU1 226 : MLREE5=IVU1 227 : MLREE6=IVN1 228 : SEGSUP MLREE4,MLREE5,MLREE6 229 : MLREE4=IVD1 230 : MLREE5=IVGM1 231 : MLREE6=IVGE1 232 : MLREE7=IVXL1 233 : SEGSUP MLREE4,MLREE5,MLREE6,MLREE7 234 : IVXL1=0 235 : IVXU1=0 236 : IVU1=0 237 : IVN1=0 238 : IVD1=0 239 : IVGM1=0 240 : IVGE1=0 241 : I=IK 242 : QQZ=PROG(I) 243 : DO 13 J=1,MLREE1.PROG(/1) 244 : MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J) 245 : 13 CONTINUE 245 : CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU, 246 : * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1) 246 : CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1, 247 : *IVGM1,IVLAM1,IPBASE) 248 : XLP3=0.D0 249 : MLREE4=IVGE1 250 : DO 14 J=1,MLREE1.PROG(/1) 251 : XLP3=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP3 252 : 14 CONTINUE 253 : IF(XLP1*XLP3.LE.0.D0) THEN 254 : XLP2=XLP3 255 : INPSUP=I 256 : ELSE 257 : XLP1=XLP3 258 : INPINF=I 259 : ENDIF 260 : GO TO 12 261 : 15 CONTINUE 262 : * 263 : * RECHERCHE PAR DICHOTOMIE pour determiner la valeur de alpha donnant le 264 : * min de l(lambda + alpha * grad) 265 : * 266 : IPP=0 267 : ALPE=PROG(INPINF) 268 : ALGR=PROG(INPSUP) 269 : ALPHA=(ALGR + ALPE)/2.D0 270 : 60 CONTINUE 271 : IPP=IPP + 1 272 : IF(IPP.GT.100) THEN 272 : CALL ERREUR ( 603) 273 : RETURN 274 : ENDIF 275 : ALPHA0=ALPHA 276 : MLREE4=IVXU1 277 : MLREE5=IVU1 278 : MLREE6=IVN1 279 : SEGSUP MLREE4,MLREE5,MLREE6 280 : MLREE4=IVD1 281 : MLREE5=IVGM1 282 : MLREE6=IVGE1 283 : MLREE7=IVXL1 284 : SEGSUP MLREE4,MLREE5,MLREE6,MLREE7 285 : IVXL1=0 286 : IVXU1=0 287 : IVU1=0 288 : IVN1=0 289 : IVD1=0 290 : IVGM1=0 291 : IVGE1=0 292 : * IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' ALPHA '',E12.5)')ALPHA 293 : QQZ=ALPHA 294 : DO 17 J=1,MLREE1.PROG(/1) 295 : MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J) 296 : 17 CONTINUE 296 : CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU, 297 : * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1) 297 : CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1, 298 : *IVGM1,IVLAM1,IPBASE) 299 : XLP3=0.D0 300 : MLREE4=IVGE1 301 : DO 50 J=1,MLREE1.PROG(/1) 302 : XLP3=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP3 303 : 50 CONTINUE 304 : IF(XLP1*XLP3.LE.0.D0) THEN 305 : XLP2=XLP3 306 : ALGR=ALPHA 307 : ELSE 308 : XLP1=XLP3 309 : ALPE=ALPHA 310 : ENDIF 311 : ALPHA =(ALGR + ALPE)/2 312 : EPSI=1.D-10 313 : IF(ALPHA.LT.1.D-12.AND.NDR.GE.1)THEN 314 : TERMIN=.TRUE. 315 : GO TO 20 316 : ENDIF 317 : IF((ABS(ALPHA-ALPHA0)/ALPHA).GE.EPSI) GO TO 60 318 : 20 CONTINUE 319 : MLENT1=MAC2 320 : MLENT2=MAC3 321 : MLENT3=MORDRE 322 : IF(IIMPI.EQ.1799)WRITE(IOIMP,FMT='('' ALPHA '',E12.5)')ALPHA 323 : IF (KK.NE.-3)THEN 324 : * 325 : ****************** TEST SI ARRET SUR PLAN DE DISCONTINUITE ************ 326 : * 327 : ALPHA1=PROG(INPSUP) 328 : IF((ABS(ALPHA1-ALPHA)/ALPHA).LE.EPSI)THEN 329 : KK=MLENT2.LECT(MLENT3.LECT(INPSUP)) 330 : II=MLENT1.LECT(MLENT3.LECT(INPSUP)) 331 : ALPHA=ALPHA1 332 : ENDIF 333 : ALPHA2=PROG(INPINF) 334 : IF(((ABS(ALPHA2-ALPHA)/ALPHA).LE.EPSI.AND.ALPHA2.NE.0.D0) 335 : * .OR.ALPHA.LT.1.D-12)THEN 336 : KK=MLENT2.LECT(MLENT3.LECT(INPINF)) 337 : II=MLENT1.LECT(MLENT3.LECT(INPINF)) 338 : ALPHA=ALPHA2 339 : ENDIF 340 : ENDIF 341 : DO 21 J=1,MLREE1.PROG(/1) 342 : MLREE3.PROG(J)=ALPHA * MLREE2.PROG(J)+MLREE1.PROG(J) 343 : 21 CONTINUE 344 : IF (KK.EQ.-3)THEN 345 : MLREE3.PROG(II)=0. 346 : ENDIF 347 : ALPHA0=ALPHA 348 : MLREE4=IVXU1 349 : MLREE5=IVU1 350 : MLREE6=IVN1 351 : SEGSUP MLREE4,MLREE5,MLREE6 352 : MLREE4=IVD1 353 : MLREE5=IVGM1 354 : MLREE6=IVGE1 355 : MLREE7=IVXL1 356 : SEGSUP MLREE4,MLREE5,MLREE6,MLREE7 357 : IVLAMB=IVLAM1 358 : SEGSUP MLENT1,MLENT2,MLREEL,MLENT3 359 : RETURN 360 : END 361 :
© Cast3M 2003 - All rights reserved.
Disclaimer