Télécharger qzbasc.eso

Retour à la liste

Numérotation des lignes :

  1. C QZBASC SOURCE CB215821 17/04/20 21:15:27 9406
  2. SUBROUTINE QZBASC (ALFR, ALFI, BETA, MATZ, MELEME, MCOMP,
  3. & IPBC,ERR,NWANTED)
  4. IMPLICIT INTEGER(I-N)
  5. IMPLICIT REAL*8 (A-H,O-Z)
  6. *
  7. ************************************************************************
  8. * CREATION D'UNE BASE DE MODES PROPRES COMPLEXES POUR VIBC
  9. * ______________________________________________________________________
  10. *
  11. * AUTEUR : Nicolas BENECH, 13 avril 1995
  12. * MODIF : BP, 2016-01-15 on reecrit tout
  13. * ______________________________________________________________________
  14. *
  15. * MODULE(S) APPELANT(S) : VIBRAC
  16. *
  17. * MODULE(S) APPELE(S) : CRTABL, ACCTAB, ECCTAB, CREPO1, MUCHPO
  18. * ______________________________________________________________________
  19. *
  20. * EN ENTREE :
  21. * -ALFR : partie reelle des valeurs propres
  22. * -ALFI : partie imaginaire des valeurs propres
  23. * -BETA : denominateur (reel) des valeurs propres
  24. * (denominateur nul --> valeur propre infinie)
  25. * -MATZ : vecteurs propres complexes
  26. * -MELEME : maillage support des chpoints
  27. * -MCOMP : liste des composantes des chpoints
  28. * ______________________________________________________________________
  29. *
  30. * EN SORTIE :
  31. * -IPBC : base de modes complexes
  32. * ______________________________________________________________________
  33. *
  34. * REMARQUE : la variable SEUIL definie au début des données générales
  35. * permet d'identifier les valeurs propres infinies, et
  36. * d'annuler le mode correspondant.
  37. * Sa valeur a ete fixee a 100 x EPSLON(1.0D0) pour le moment
  38. ************************************************************************
  39. *
  40. -INC CCREEL
  41. *-
  42. -INC CCOPTIO
  43. -INC SMELEME
  44. -INC SMLCHPO
  45. -INC SMLMOTS
  46. -INC SMCHPOI
  47. -INC SMRIGID
  48. -INC SMTABLE
  49. -INC SMLREEL
  50. -INC SMLENTI
  51. -INC SMCOORD
  52. * segment de travail = MTRAV modifie pour NNCHPO chpoints
  53. SEGMENT MTRAV
  54. CHARACTER*4 INCO(NNIN)
  55. REAL*8 BB(NNCHPO,NNIN,NNNOE)
  56. INTEGER IBIN(NNIN,NNNOE),IGEO(NNNOE),NHAR(NNIN)
  57. ENDSEGMENT
  58. *
  59. * ******** INCO(NNIN) CONTIENT LE NOMS DES NNIN INCONNUES DIFFERENTES.
  60. *
  61. * ******** BB(I,J) EST LA VALEUR DE LA IEME INCONNUE DU CHAMP POUR
  62. * ******** LE JEME NOEUD DU TABLEAU IGEO.
  63. *
  64. * ******** IBIN(I,J)=1 OU 0. 1 INDIQUE QUE LA I EME INCONNUE DU CHAMP
  65. * ******** EXISTE POUR LE J EME NOEUD DU TABLEAU IGEO.
  66. *
  67. * ******** IGEO(I) EST LE NUMERO A METTRE DANS UN OBJET MELEME POUR
  68. * ******** REFERENCER LE IEME NOEUD
  69. *
  70. * ******** NHAR(I) EST LE NUMERO D'HARMONIQUE SI CALCUL AXI OU
  71. * ******** SIGNIFIE CONTRAINTE PLANE,DEFORMATION PLANE OU DEF PLAN GEN
  72. *
  73. SEGMENT ICPR(XCOOR(/1)/(IDIM+1))
  74. *
  75. REAL*8 XVAL,SEUIL
  76. INTEGER I, J, K, IC, NUMAFF,ERR
  77. LOGICAL MODANN, MODREL, AFFICH,CONV
  78. CHARACTER*4 NOMDDL
  79. *
  80. POINTEUR ALFR.MLREEL, ALFI.MLREEL, BETA.MLREEL
  81. POINTEUR MATZ.XMATRI
  82. POINTEUR MCOMP.MLMOTS
  83. *
  84.  
  85. ************************************************************************
  86. * DONNEES GENERALES *
  87. ************************************************************************
  88.  
  89. * Ecriture des messages pour verification
  90.  
  91. AFFICH = IIMPI.GE.21
  92. c AFFICH = .TRUE.
  93. *
  94. IF (AFFICH)
  95. & WRITE (IOIMP,*) 'QZBASC: Extraction des donnees generales...'
  96. *
  97. * seuil pour le denominateur d'une valeur propre
  98.  
  99. * SEUIL = (EPSLON(1.0D0)*100)
  100. SEUIL = 1.D-99
  101.  
  102. * nombre de modes calcules
  103. NBMOD1=ALFR.PROG(/1)
  104.  
  105.  
  106. **** si on souhaite un nombre de modes inferieur, on devra alors trier
  107.  
  108. IF(NWANTED.GT.0.AND.NWANTED.LT.NBMOD1) THEN
  109. c tri selon 1/module**2 : il faut le calculer explicitement
  110. SEGINI,MLREEL=BETA
  111. I=1
  112. 88 IF(I.GT.NBMOD1) GOTO 22
  113. PROG(I) = (PROG(I)**2) / (ALFI.PROG(I)**2 + ALFR.PROG(I)**2)
  114. c bp : si c'est un complexe, on force l egalite du module
  115. c comme le tri d'ORDON1 est stable, les complexes i et i+1
  116. c resteront a des indices successifs
  117. IF(ABS(ALFI.PROG(I)).NE. 0.0D0) THEN
  118. I=I+1
  119. PROG(I)=PROG(I-1)
  120. ENDIF
  121. I=I+1
  122. GOTO 88
  123. 22 CONTINUE
  124. IORDRE=1
  125. CALL ORDON1(MLREEL,.FALSE.,.TRUE.,IORDRE)
  126. segsup,MLREEL
  127. MLENTI=IORDRE
  128. segact,MLENTI
  129. ELSE
  130. IORDRE=0
  131. NWANTED=NBMOD1
  132. ENDIF
  133.  
  134.  
  135. **** preparation des CHPOINTs deformees *******************************
  136.  
  137. * nombre de ddls (=dimension de MELEME = dimension de MCOMP)
  138. SEGACT,MELEME,MCOMP
  139. NDDL=MELEME.NUM(/2)
  140.  
  141. * on cree ICPR
  142. SEGINI,ICPR
  143. * on dimensionne au maxi MTRAV
  144. NNIN=NDDL
  145. NNNOE=NDDL
  146. NNCHPO=NBMOD1
  147. SEGINI,MTRAV
  148. NNIN=0
  149. NNNOE=0
  150.  
  151. *---- traitement de chaque ddl ----------------------------------
  152. DO 100 J=1,NDDL
  153.  
  154. c - NOEUD
  155. IP=NUM(1,J)
  156. JNOE=ICPR(IP)
  157. c nouveau noeud #IP de numero local #JNOE
  158. IF(JNOE.EQ.0) THEN
  159. NNNOE=NNNOE+1
  160. JNOE=NNNOE
  161. ICPR(IP)=JNOE
  162. IGEO(JNOE)=IP
  163. ENDIF
  164.  
  165. c - COMPOSANTE
  166. NOMDDL=MCOMP.MOTS(J)
  167. IF(NNIN.EQ.0) GOTO 111
  168. DO 110 JIN=1,NNIN
  169. IF(NOMDDL.EQ.INCO(JIN)) GOTO 112
  170. c NOMDDL trouvee dans INCO(JIN)
  171. 110 CONTINUE
  172.  
  173. 111 CONTINUE
  174. c NOMDDL pas trouvee dans INCO(JIN) -> on l'ajoute a la fin
  175. NNIN=NNIN+1
  176. JIN=NNIN
  177. INCO(JIN)=NOMDDL
  178. c bp : pour l'instant on laisse NHAR=NIFOUR, mais il faudrait
  179. c recuperer IRIGEL(5,:) des rigidites d'entree ou autre ...?
  180. NHAR(JIN)=NIFOUR
  181.  
  182. 112 CONTINUE
  183. c NOMDDL trouvee dans INCO(JIN) et noeud #IP dans IGEO(JNOE)
  184.  
  185. c Remplissage de IBIN et BB
  186. IBIN(JIN,JNOE)=J
  187. DO I=1,NBMOD1
  188. BB(I,JIN,JNOE)=MATZ.RE(J,I,1)
  189. ENDDO
  190.  
  191. 100 CONTINUE
  192. *---- fin de la boucle sur les ddl ----------------------------------
  193.  
  194. SEGSUP,MCOMP,ICPR
  195. SEGDES,MELEME
  196. SEGADJ,MTRAV
  197. c write(*,*) 'IGEO=',(IGEO(iou),iou=1,NNNOE)
  198. c write(*,*) 'INCO=',(INCO(iou),iou=1,NNIN)
  199. c on va generer N1 chpoints en 1 passage dans une copie de CRECHP
  200. N1=NBMOD1
  201. SEGINI,MLCHPO
  202. CALL CRECH4(MTRAV,MLCHPO)
  203. SEGSUP,MTRAV
  204.  
  205.  
  206. *******************************************
  207. * Creation de la table BASE_DE_MODES *
  208. *******************************************
  209. *
  210. IF (AFFICH) WRITE(IOIMP,*) 'Creation de la table BASE_DE_MODES...'
  211. *
  212. CALL CRTABL(IPTAB2)
  213. CALL ECCTAB(IPTAB2,'MOT',0,0.0D0,'SOUSTYPE',.TRUE.,0,
  214. & 'MOT',0,0.0D0,'BASE_DE_MODES',.TRUE.,0)
  215. CALL ECCTAB(IPTAB2,'MOT',0,0.0D0,'MAILLAGE',.TRUE.,0,
  216. & 'MAILLAGE',0,0.0D0,' ',.TRUE.,MELEME)
  217.  
  218.  
  219. ************************************************************************
  220. * BOUCLE SUR LES MODES *
  221. ************************************************************************
  222. *
  223. III=1
  224. 80 IF (III .GT. NWANTED) GOTO 20
  225.  
  226. * selection des NWANTED plus petits modes
  227. IF(IORDRE.NE.0) THEN
  228. I = LECT(III)
  229. ELSE
  230. I = III
  231. ENDIF
  232. *
  233. MODREL = (ABS(ALFI.PROG(I)).EQ. 0.0D0)
  234. IF (AFFICH) THEN
  235. IF (MODREL) THEN
  236. WRITE (IOIMP,*) III,': MODE ',I,' REEL',ALFI.PROG(I),
  237. & '-i',ALFR.PROG(I),' / 2pi',BETA.PROG(I)
  238. ELSE
  239. WRITE (IOIMP,*) III,': MODE ',I,' COMPLEXE',ALFI.PROG(I),
  240. & '-i',ALFR.PROG(I),' / 2pi',BETA.PROG(I)
  241. ENDIF
  242. ENDIF
  243. *
  244. * ---- frequence infinie ?
  245. *
  246. MODANN = (ABS(BETA.PROG(I)).LT.SEUIL)
  247.  
  248. * ---- oui
  249. IF (MODANN) THEN
  250. 90 WRITE (IOIMP,*) 'Attention !!! Mode ',I,
  251. & ' annule : frequence infinie.'
  252. CALL CRTABL(MTAB1)
  253. CALL ECCTAB(MTAB1,'MOT',0,0.0D0,'SOUSTYPE', .TRUE.,0,
  254. & 'MOT',0,0.0D0,'MODE_ANNULE',.TRUE.,0)
  255. CALL ECCTAB(IPTAB2,'ENTIER',I,0.0D0,' ',.TRUE.,0,
  256. & 'TABLE', 0, 0.0D0,' ',.TRUE.,MTAB1)
  257. IF (.NOT. MODREL) THEN
  258. MODREL = .TRUE.
  259. III = III + 1
  260. GO TO 90
  261. ENDIF
  262.  
  263. * ---- non
  264. ELSE
  265. *
  266. *------- Deformees reelle (1) et imaginaire (2) du mode i
  267.  
  268. * Recup des CHPOINTs :
  269. MCHPO1=ICHPOI(I)
  270.  
  271. * Cas reel : il faut un chpoint nul -> on crée un chpoint vide
  272. IF (MODREL) THEN
  273. NAT=1
  274. NSOUPO=0
  275. SEGINI,MCHPO2
  276. MCHPO2.IFOPOI=IFOMOD
  277. MCHPO2.MOCHDE='CHPOINT CREE PAR QZBASC'
  278. * Cas complexe : on fait la Recup
  279. ELSE
  280. MCHPO2=ICHPOI(I+1)
  281. ENDIF
  282. c write(*,*) '>>> mode',I,'phiR=',MCHPO1,'phiI=',MCHPO2
  283. *
  284. *
  285. ************************************************
  286. * Creation de la table MODE *
  287. ************************************************
  288. *
  289. *----- valeur propre (reelle ou complexe)
  290. *
  291. * Attention : on enregistre f = w/2pi
  292. * avec iw = landa donc w=-i landa !!!
  293. *
  294. IF (AFFICH) WRITE (*,*) 'Construction de la table MODE ...'
  295. *
  296. CALL CRTABL(MTAB1)
  297. CALL ECCTAB(MTAB1,'MOT',0,0.0D0,'SOUSTYPE',.TRUE.,0,
  298. & 'MOT',0,0.0D0,'MODE_COMPLEXE',.TRUE.,0)
  299. *
  300. CALL ECCTAB(MTAB1,'MOT',0,0.0D0,'NUMERO_MODE',.TRUE.,0,
  301. & 'ENTIER',I,0.0D0,' ',.TRUE.,0)
  302. CALL CREPO1(0.0D0,0.0D0,0.0D0,IPOIN)
  303. CALL ECCTAB(MTAB1,'MOT',0,0.0D0,'POINT_REPERE',.TRUE.,0,
  304. & 'POINT',0,0.0D0,' ',.TRUE.,IPOIN)
  305. XVAL=(-1*ALFR.PROG(I))/(2.*XPI*BETA.PROG(I))
  306. CALL ECCTAB(MTAB1,'MOT',0,0.0D0,'FREQUENCE_IMAGINAIRE',.TRUE.,0
  307. $ ,'FLOTTANT',0,XVAL,' ',.TRUE.,0)
  308. XVAL=(ALFI.PROG(I))/(2.*XPI*BETA.PROG(I))
  309. CALL ECCTAB(MTAB1,'MOT',0,0.0D0,'FREQUENCE_REELLE',.TRUE.,0,
  310. & 'FLOTTANT',0,XVAL,' ',.TRUE.,0)
  311. *
  312. CALL ECCTAB(MTAB1,'MOT',0,0.0D0,'DEFORMEE_MODALE_REELLE',
  313. & .TRUE.,0,'CHPOINT',0,0.0D0,' ',.TRUE.,MCHPO1)
  314. CALL ECCTAB(MTAB1,'MOT',0,0.0D0,'DEFORMEE_MODALE_IMAGINAIRE',
  315. & .TRUE.,0,'CHPOINT',0,0.0D0,' ',.TRUE.,MCHPO2)
  316. *
  317. CALL ECCTAB(IPTAB2,'ENTIER',III,0.0D0,' ',.TRUE.,0,'TABLE',
  318. & 0,0.0D0,' ',.TRUE.,MTAB1)
  319. *
  320. IF (MODREL) GOTO 70
  321. *
  322. *----- Valeur propre complexe conjuguee
  323. *
  324. *
  325. SEGINI, MTAB2=MTAB1
  326. SEGDES, MTAB2
  327. III = III + 1
  328. IF (AFFICH) WRITE (*,*) 'Mode conjugue, no : ',I
  329. CALL ECCTAB(MTAB2,'MOT',0,0.0D0,'NUMERO_MODE',.TRUE.,0,
  330. & 'ENTIER',I,0.0D0,' ',.TRUE.,0)
  331. CALL ECCTAB(MTAB2,'MOT',0,0.0D0,'FREQUENCE_REELLE',.TRUE.,0,
  332. & 'FLOTTANT',0,-1.*XVAL,' ',.TRUE.,0)
  333. XVAL=-1.D0
  334. CALL MUCHPO(MCHPO2,XVAL,IRET,1)
  335. CALL ECCTAB(MTAB2,'MOT',0,0.0D0,'DEFORMEE_MODALE_REELLE',
  336. & .TRUE.,0,'CHPOINT',0,0.0D0,' ',.TRUE.,MCHPO1)
  337. CALL ECCTAB(MTAB2,'MOT',0,0.0D0,'DEFORMEE_MODALE_IMAGINAIRE',
  338. & .TRUE.,0,'CHPOINT',0,0.0D0,' ',.TRUE.,IRET)
  339. *
  340. CALL ECCTAB(IPTAB2,'ENTIER',III,0.0D0,' ',.TRUE.,0,'TABLE',
  341. & 0,0.0D0,' ',.TRUE.,MTAB2)
  342. *
  343. IF (AFFICH) WRITE (*,*) 'Construction du conjugue ok'
  344. 70 CONTINUE
  345. ENDIF
  346.  
  347.  
  348. III = III + 1
  349. GOTO 80
  350. ************************************************************************
  351. * FIN DE BOUCLE SUR LES MODES *
  352. ************************************************************************
  353. 20 CONTINUE
  354.  
  355.  
  356. ***** MENAGE ***********************************************************
  357.  
  358. SEGSUP, ALFR, ALFI, BETA
  359. SEGSUP, MATZ
  360. IF(IORDRE.GT.0) SEGSUP,MLENTI
  361.  
  362. *
  363. ********************************************
  364. * Creation de la table BASE_MODALE *
  365. ********************************************
  366. *
  367. IF (AFFICH) WRITE (*,*) 'Creation de la table BASE_MODALE...'
  368. *
  369. CALL CRTABL(IPBC)
  370. CALL ECCTAB(IPBC,'MOT',0,0.0D0,'SOUSTYPE',.TRUE.,0,
  371. & 'MOT',0,0.0D0,'BASE_MODALE',.TRUE.,0)
  372. CALL ECCTAB(IPBC,'MOT',0,0.0D0,'MODES',.TRUE.,0,
  373. & 'TABLE',0,0.0D0,' ',.TRUE.,IPTAB2)
  374. cbp : on ajoute un indicateur de succes (en plus du message dans vibrac)
  375. CONV=ERR.eq.0
  376. CALL ECCTAB(IPBC,'MOT',0,0.0D0,'CONVERGENCE',.TRUE.,0,
  377. & 'LOGIQUE ',0,0.0D0,' ',CONV,0)
  378. *
  379. RETURN
  380. END
  381.  
  382.  

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