Télécharger ntape6.eso

Retour à la liste

Numérotation des lignes :

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

© Cast3M 2003 - Tous droits réservés.
Mentions légales