Télécharger arpitq.eso

Retour à la liste

Numérotation des lignes :

  1. C ARPITQ SOURCE BP208322 20/02/06 21:15:09 10512
  2. SUBROUTINE ARPITQ (IPRTRA,IPMAUP,SIGMA,INVER,EPSI,OUT)
  3.  
  4.  
  5. ***********************************************************************
  6. *
  7. * A R P I T Q
  8. *
  9. * FONCTION:
  10. * ---------
  11. *
  12. * STEP DE LA FACTORISATION D'ARNOLDI POUR UN PROBLEME QUADRATIQUE.
  13. *
  14. * REMARQUE:
  15. * ---------
  16. *
  17. * LE PROBLEME EST DE DIMENSION DOUBLE PAR RAPPORT AU NOMBRE
  18. * D'INCONNUES: CEPENDANT LES FONCTIONS DISPONIBLES POUR LA
  19. * RESOLUTION DE PROBLEMES D'ORDRE 1 SERONT UTILISES CAR LES
  20. * OPERATIONS MATRICE-VECTEUR SE FONT PAR BLOCS
  21. *
  22. * RIGI(1) : RIGIDITE
  23. * RIGI(3) : AMORTISSEMENT
  24. * RIGI(2) : MASSE
  25. *
  26. *
  27. * PARAMETRES: (E)=ENTREE (S)=SORTIE
  28. * -----------
  29. *
  30. *
  31. * IPRTRA ENTIER (E) OPERATEURS DE TRAVAIL
  32. *
  33. * IPMAUP ENTIER (E/S) POINTEUR VARIABLES ARPACK
  34. *
  35. * SIGMA COMPLEX DP (E) VALEUR DU SHIFT (PEUT ETRE NUL)
  36. *
  37. * INVER LOGIQUE (E) .TRUE. -> PRODUIT SCALAIRE X'KX
  38. * .FALSE. -> PRODUIT SCALAIRE X'MX
  39. *
  40. * EPSI REEL DP (E) TOLERANCE
  41. *
  42. * OUT LOGIQUE (S) FLAG DE CONVERGENCE
  43. *
  44. *
  45. * SOUS-PROGRAMMES APPELES:
  46. * ------------------------
  47. *
  48. * DNAUPD,ADCHPO,ARPCH2,MUCPRI,RESOU1,LDMT,DTCHPO
  49. *
  50. * AUTEUR,DATE DE CREATION:
  51. * -------------------------
  52. *
  53. * PASCAL BOUDA 17 JUILLET 2015
  54. *
  55. * LANGAGE:
  56. * --------
  57. *
  58. * FORTRAN 77 & 90
  59. *
  60. ***********************************************************************
  61.  
  62. -INC CCOPTIO
  63. -INC SMLCHPO
  64. -INC TARWORK
  65. c -INC TARTRAK
  66.  
  67. SEGMENT idemen(0)
  68.  
  69.  
  70. INTEGER IPRTRA
  71. INTEGER IPMAUP
  72. LOGICAL INVER
  73. COMPLEX*16 SIGMA
  74. REAL*8 EPSI
  75. LOGICAL OUT
  76.  
  77.  
  78. INTEGER IPRIGI,IPMASS,IPAMOR,IPSCAL
  79. INTEGER TEST
  80. INTEGER IPLCHP
  81. INTEGER IPCTRA(9)
  82. INTEGER NOID,NOEN
  83. INTEGER ndim,ncv,lworkl
  84.  
  85.  
  86. SEGINI IDEMEN
  87. IDEMEN(**)=0
  88. NOID=0
  89. NOEN=1
  90.  
  91. OUT=.FALSE.
  92.  
  93.  
  94. MRITRA=IPRTRA
  95. SEGACT MRITRA
  96. IPRIGI=RIGI(1)
  97. IPMASS=RIGI(2)
  98. IPAMOR=RIGI(3)
  99. IPCSIM=RIGI(5)
  100. IPKCSM=RIGI(6)
  101.  
  102.  
  103. MAUP=IPMAUP
  104. SEGACT MAUP*MOD
  105.  
  106. *Récupération de la dimension des tableaux
  107. ndim=resid(/1)
  108. ncv=v(/2)
  109. lworkl=workl(/1)
  110.  
  111. *Le probleme est non symetrique. Appel de la routine ARPACK specifique
  112. CALL DNAUPD (ido,bmat,ndim,which,nev,EPSI,resid,
  113. & ncv,v,ldv,iparam,ipntr,workd,workl,lworkl,ITRAK,info)
  114.  
  115.  
  116. *Reverse communication: On récupère les paramètres de sortie et on
  117. *effectue des actions en fonction de leurs valeurs
  118. TEST=ido
  119.  
  120. IPMAUP=MAUP
  121. c SEGDES MAUP
  122.  
  123.  
  124. *On verifie si on doit stopper le processus
  125. CALL ARPERR (IPMAUP,.FALSE.,OUT)
  126. IF (OUT) RETURN
  127.  
  128.  
  129. *Initialisation des chpoints de travail
  130. DO i=1,9
  131. IPCTRA(i)=0
  132. ENDDO
  133.  
  134.  
  135.  
  136. IF (TEST .EQ. -1 .OR. TEST .EQ. 1) THEN
  137.  
  138. * &--------------------------------------------------&
  139. * | Calcul du produit matrice vecteur |
  140. * | Y <---- RE( inv(A-SIGMA*B)*B*X ) |
  141. * | |
  142. * | X : workd(ipntr(1)) |
  143. * | Y : workd(ipntr(2)) |
  144. * &--------------------------------------------------&
  145.  
  146. *Note: Initialement, ARPACK est prevu pour distinguer test=-1 ou on
  147. *fournit inv(A-SIGMA*B)*B*X et test=1 ou on fournit inv(A-SIGMA*B)*X
  148. *(economie du produit par la matrice masse deja effectue).
  149. *Cependant,la "logique" cast3m nous retire la garantie de fournir le
  150. *produit pour tout type de matrice.
  151. *En outre, meme si la reussite de l'operation etait garantie, elle
  152. *serait plus lourde a effectuer, cette derniere impliquant l'inversion
  153. *d'une matrice de masse
  154.  
  155. IF (AIMAG(SIGMA) .EQ. 0.D0) THEN
  156.  
  157. *Shift reel : RE( inv(A-SIGMA*B)*B*X ) = inv(A-REAL(SIGMA)*B)*B*X
  158.  
  159. CALL ARPCH2 (IPRIGI,IPRIGI,IPMAUP,IPLCHP,1,3)
  160. MLCHPO=IPLCHP
  161. SEGACT MLCHPO
  162.  
  163.  
  164. IPCTRA(3)=ICHPOI(1)
  165.  
  166. CALL MUCPRI (IPCTRA(3),IPMASS,IPCTRA(6))
  167.  
  168. IPCTRA(1)=ICHPOI(2)
  169.  
  170. CALL MUCPRI (IPCTRA(1),IPCSIM,IPCTRA(7))
  171.  
  172. CALL ADCHPO (IPCTRA(6),IPCTRA(7),IPCTRA(2),-1.D0,-1.D0)
  173.  
  174. SEGDES MLCHPO
  175.  
  176. IDEMEN(1)=IPCTRA(2)
  177.  
  178. IF (SYME(6) .EQ. 0) THEN
  179. CALL RESOU1 (IPKCSM,IDEMEN,NOID,NOEN,1.D-18,0)
  180. ELSEIF (SYME(6) .EQ. 1) THEN
  181. CALL LDMT (IPKCSM,IDEMEN,NOID,NOEN,1.D-18)
  182. ENDIF
  183. IF (IERR.NE.0) RETURN
  184. C RESOU desactive IDEMEN on le reactive
  185. SEGACT,IDEMEN
  186.  
  187.  
  188. IPCTRA(4)=IDEMEN(1)
  189.  
  190. CALL ADCHPO(IPCTRA(4),IPCTRA(1),IPCTRA(5),REAL(SIGMA),1.D0)
  191.  
  192. MLCHPO=IPLCHP
  193. SEGACT MLCHPO*MOD
  194. ICHPOI(1)=IPCTRA(5)
  195. ICHPOI(2)=IPCTRA(4)
  196. IPLCHP=MLCHPO
  197. SEGDES MLCHPO
  198.  
  199. CALL ARPCH2 (IPRIGI,IPRIGI,IPMAUP,IPLCHP,2,1)
  200.  
  201. ELSE
  202.  
  203. *Shift complexe : RE( inv(A-SIGMA*B)*B*X ) a fournir (non implemente)
  204.  
  205. ENDIF
  206.  
  207.  
  208. ELSEIF (TEST .EQ. 2) THEN
  209.  
  210. * &------------------------------------------&
  211. * | Calcul du produit matrice vecteur |
  212. * | |
  213. * | Si INVER |
  214. * | Y <---- DIAG(IPRIGI,IPRIGI)*X |
  215. * | |
  216. * | Sinon |
  217. * | Y <---- DIAG(IPMASS,IPMASS)*X |
  218. * | |
  219. * | X : workd(ipntr(1)) |
  220. * | Y : workd(ipntr(2)) |
  221. * &------------------------------------------&
  222.  
  223. IF (INVER) THEN
  224. IPSCAL=IPRIGI
  225. ELSE
  226. IPSCAL=IPMASS
  227. ENDIF
  228.  
  229. CALL ARPCH2 (IPRIGI,IPRIGI,IPMAUP,IPLCHP,1,3)
  230.  
  231. MLCHPO=IPLCHP
  232. SEGACT MLCHPO
  233.  
  234. IPCTRA(1)=ICHPOI(1)
  235.  
  236. CALL MUCPRI (IPCTRA(1),IPSCAL,IPCTRA(3))
  237.  
  238. IPCTRA(2)=ICHPOI(2)
  239.  
  240. CALL MUCPRI (IPCTRA(2),IPSCAL,IPCTRA(4))
  241.  
  242. SEGDES MLCHPO
  243.  
  244. MLCHPO=IPLCHP
  245. SEGACT MLCHPO*MOD
  246. ICHPOI(1)=IPCTRA(3)
  247. ICHPOI(2)=IPCTRA(4)
  248.  
  249. IPLCHP=MLCHPO
  250. SEGDES MLCHPO
  251.  
  252. CALL ARPCH2 (IPRIGI,IPRIGI,IPMAUP,IPLCHP,2,2)
  253.  
  254.  
  255. ENDIF
  256.  
  257.  
  258. *Destruction des chpoints de travail
  259. DO i=1,9
  260. IF (IPCTRA(i) .NE. 0) THEN
  261. CALL DTCHPO (IPCTRA(i))
  262. ENDIF
  263. ENDDO
  264.  
  265. SEGDES MRITRA
  266.  
  267.  
  268. END
  269.  
  270.  
  271.  
  272.  
  273.  

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