Télécharger acheck.eso

Retour à la liste

Numérotation des lignes :

  1. C ACHECK SOURCE BP208322 19/04/29 21:15:01 10213
  2. SUBROUTINE ACHECK (IPRIGI,IPMASS,QUAD,SYM,SHIFT,N,FLAG,
  3. & INVER,PIRE,CHOLE,EPSI)
  4.  
  5.  
  6. ***********************************************************************
  7. *
  8. * A C H E C K
  9. *
  10. * FONCTION:
  11. * ---------
  12. *
  13. * VERIFICATION DE LA POSSIBILITE POUR ARPACK DE RESOUDRE +
  14. * CHOIX DE LA MATRICE QUI DEFINIRA LE PRODUIT SCALAIRE DANS ARPACK
  15. * + EVENTUELLE(S) SIMPLIFICATION(S) (CHOLESKY, PROBLEME SYM)
  16. *
  17. *
  18. * PARAMETRES: (E)=ENTREE (S)=SORTIE
  19. * -----------
  20. *
  21. * IPRIGI ENTIER (E) POINTEUR D'UNE RIGIDITE
  22. * IPMASS ENTIER (E) POINTEUR D'UNE MASSE
  23. * QUAD LOGIQUE (E) PROBLEME QUADRATIQUE OU NON
  24. * SYM LOGIQUE (S) PROBLEME SYMETRIQUE OU NON
  25. * SHIFT COMPLEX DP (E) FREQUENCE DE SHIFT
  26. * N ENTIER (E) DIMENSION DU PROBLEME
  27. * FLAG LOGIQUE (S) PROBLEME SOLVABLE OU NON
  28. * INVER LOGIQUE (S) .TRUE. -> PRODUIT SCALAIRE X'KX
  29. * .FALSE. -> PRODUIT SCALAIRE X'MX
  30. * CHOLE LOGIQUE (S) CHOLESKY NON ALTERNATIVE POSSIBLE
  31. * EPSI REEL DP (E) ZERO DE TOLERANCE
  32. *
  33. *
  34. * SOUS-PROGRAMMES APPELES:
  35. * ------------------------
  36. *
  37. * DIAGN1
  38. *
  39. * AUTEUR, DATE DE CREATION:
  40. * -------------------------
  41. *
  42. * PASCAL BOUDA 29 JUIN 2015
  43. *
  44. ***********************************************************************
  45.  
  46. IMPLICIT INTEGER(I-N)
  47. IMPLICIT REAL*8 (A-H,O-Z)
  48.  
  49. -INC CCOPTIO
  50. -INC SMRIGID
  51.  
  52. INTEGER IPRIGI,IPMASS
  53. INTEGER N
  54. LOGICAL QUAD
  55. LOGICAL SYM,INVER,POSITI
  56. LOGICAL FLAG,PIRE
  57. LOGICAL CHOLE
  58. COMPLEX*16 SHIFT
  59.  
  60.  
  61. INTEGER IPCHOI
  62. INTEGER nvp0K,nvp0M
  63. INTEGER NRG,NBR,IANTI
  64. COMPLEX*16 ZERO
  65.  
  66.  
  67. ZERO=CMPLX(0.D0,0.D0)
  68.  
  69. * La decomposition de Cholesky n'est pas codee --> CHOLE toujours false
  70. CHOLE=.FALSE.
  71.  
  72. *** Cas lineaire:
  73. * Si le shift est nul, on peut resoudre tous les problemes (matrice de
  74. * masse ou rig utilisée pour le produit scalaire,
  75. * sinon matrice identité)
  76. * Sinon, par defaut, la matrice utilisee pour le produit scalaire est la
  77. * matrice de masse
  78. * Calcul du nombre de termes diagonaux negatifs:
  79. * -si nul, ok
  80. * -sinon on essaie d'echanger les roles
  81. * (-> K pour le produit scalaire)
  82. * -si nouvel echec,le probleme n'est pas solvable
  83. *
  84. *** Cas quadratique:
  85. * K ou M doit être symetrique semi-definie positive pour le produit
  86. *scalaire.
  87. * Plus precisement, il s'agit de la matrice par blocs
  88. * | M 0 | | K 0 |
  89. * | 0 M | ou | 0 K |
  90. * Il n'y a pas de conditions sur les autres matrices
  91.  
  92. FLAG=.FALSE.
  93. INVER=.FALSE.
  94. PIRE=.FALSE.
  95. POSITI=.TRUE.
  96.  
  97. c --on va tester M--
  98. MRIGID=IPMASS
  99. SEGACT MRIGID
  100. NRG = IRIGEL(/1)
  101. NBR = IRIGEL(/2)
  102.  
  103. IF (NRG .GE. 7) THEN
  104. DO i=1,NBR
  105. IANTI=IRIGEL(7,i)
  106. IF (IANTI .GT. 0) THEN
  107. SEGDES MRIGID
  108. GOTO 101
  109. ENDIF
  110. ENDDO
  111. ENDIF
  112.  
  113. SEGDES MRIGID
  114.  
  115. CALL DIAGN1(IPMASS,nvp0M)
  116. * M def >0 ou semi-def >0
  117. IF (nvp0M .EQ. 0) THEN
  118. FLAG=.TRUE.
  119. GOTO 200
  120. ENDIF
  121.  
  122. 101 CONTINUE
  123.  
  124.  
  125. c cas M non-symetrique ou M non semi-def >0
  126.  
  127. c IF (SHIFT .NE. ZERO) THEN
  128.  
  129. c --on va tester K--
  130. MRIGID=IPRIGI
  131. SEGACT MRIGID
  132. NRG = IRIGEL(/1)
  133. NBR = IRIGEL(/2)
  134.  
  135. IF (NRG .GE. 7) THEN
  136. DO i=1,NBR
  137. IANTI=IRIGEL(7,i)
  138. IF (IANTI .GT. 0) THEN
  139. SEGDES MRIGID
  140. GOTO 102
  141. ENDIF
  142. ENDDO
  143. ENDIF
  144.  
  145. SEGDES MRIGID
  146.  
  147. CALL DIAGN1(IPRIGI,nvp0K)
  148. * K def >0 ou semi-def >0
  149. IF (nvp0K .EQ. 0) THEN
  150. FLAG=.TRUE.
  151. INVER=.TRUE.
  152. GOTO 200
  153. ENDIF
  154.  
  155. c ENDIF
  156.  
  157. 102 CONTINUE
  158. c --cas M et K non-symetrique ou non semi-def >0 --
  159.  
  160.  
  161. *-- cas 'desespere' : aucune matrice n'est bien conditionnee;
  162. * on ne peut resoudre que pour des problemes lineaires a shift nul --
  163. * bp, 2019 : a tester ...
  164. *
  165. * le probleme (matrice A=M^-1*K) n'est (probablement) pas symetrique
  166.  
  167. IF (.NOT. FLAG) THEN
  168. c on peut resoudre avec un shift nul
  169. IF (.NOT. QUAD) THEN
  170. IF (SHIFT .EQ. ZERO) THEN
  171. SYM=.FALSE.
  172. PIRE=.TRUE.
  173. FLAG=.TRUE.
  174. GOTO 300
  175. ENDIF
  176. ENDIF
  177. ENDIF
  178.  
  179.  
  180. 200 CONTINUE
  181. c --on a M ou K symetrique semi-def>0 --
  182.  
  183. * Identification du type de probleme :
  184. * -symetrique -> modes propres reels
  185. * -non symetrique -> modes propres reels ou complexes
  186.  
  187. IF (QUAD) THEN
  188. *Le probleme n'est jamais symetrique
  189. SYM=.FALSE.
  190.  
  191. ELSE
  192. *On regarde la symetrie de la matrice non utilisee pour le ps
  193. SYM=.TRUE.
  194.  
  195. c M utilisee pour le ps : on regarde K
  196. IF (.NOT. INVER) THEN
  197. MRIGID=IPRIGI
  198. SEGACT MRIGID
  199. NRG = IRIGEL(/1)
  200. NBR = IRIGEL(/2)
  201. IF (NRG .GE. 7) THEN
  202. DO i=1,NBR
  203. IANTI=IRIGEL(7,i)
  204. IF (IANTI .GT. 0) THEN
  205. SYM=.FALSE.
  206. ENDIF
  207. ENDDO
  208. ENDIF
  209.  
  210. ELSE
  211. c K utilisee pour le ps : on regarde M
  212. * La matrice M doit etre symetrique si on a inverse les roles
  213. MRIGID=IPMASS
  214. SEGACT MRIGID
  215. NRG = IRIGEL(/1)
  216. NBR = IRIGEL(/2)
  217. IF (NRG .GE. 7) THEN
  218. DO i=1,NBR
  219. IANTI=IRIGEL(7,i)
  220. IF (IANTI .GT. 0) THEN
  221. FLAG=.FALSE.
  222. INVER=.FALSE.
  223. PIRE=.FALSE.
  224. GOTO 102
  225. ENDIF
  226. ENDDO
  227. ENDIF
  228.  
  229. ENDIF
  230.  
  231. SEGDES MRIGID
  232.  
  233. ENDIF
  234.  
  235.  
  236.  
  237. 300 CONTINUE
  238.  
  239. c ERREUR !
  240. IF (.NOT.FLAG) THEN
  241.  
  242. cbp, 2019: "Au moins l'une des matrices doit etre sym. definie positive"
  243. CALL ERREUR(1097)
  244.  
  245. c WRITE(IOIMP,*) 'VIBR ne peut pas resoudre ce probleme :'
  246. IF(QUAD) THEN
  247. c WRITE(IOIMP,*) 'M n est pas symetrique semi-definie positive !'
  248. ELSE
  249. IF(INVER) THEN
  250. IF(SHIFT .NE. ZERO) THEN
  251. c WRITE(IOIMP,*) 'K n est pas symetrique semi-definie positive'
  252. c WRITE(IOIMP,*) 'ou M n est pas symetrique !'
  253. ELSE
  254. c WRITE(IOIMP,*) 'on ne doit pas passer par la !'
  255. CALL ERREUR(5)
  256. ENDIF
  257. ELSE
  258. c WRITE(IOIMP,*) 'on ne doit pas passer par la !'
  259. CALL ERREUR(5)
  260. ENDIF
  261. ENDIF
  262. ENDIF
  263.  
  264.  
  265. END
  266.  
  267.  
  268.  
  269.  
  270.  
  271.  

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