1 : C PROJRG SOURCE BP208322 10/12/01 21:15:14 6798 2 : SUBROUTINE PROJRG(MRIGID,MTAB1,POS,RI1,RI2) 3 : IMPLICIT INTEGER(I-N) 4 : IMPLICIT REAL*8 (A-H,O-Z) 5 : * 6 : *********************************************************** 7 : * PROJECTION D'UNE MATRICE SUR UNE BASE DE MODES * 8 : * _______________________________________________________ * 9 : * * 10 : * CREATION : Nicolas BENECH (11 Avril 1995) * 11 : * MODIFICATIONS : * 12 : * -18/11/2010 Benoit PRABEL : AMELIORATION PERFORMANCE * 13 : * _______________________________________________________ * 14 : * * 15 : * MODULE(S) APPELANT(S) : PJBA * 16 : * * 17 : * MODULE(S) APPELE(S) : ACCTAB, YTMX * 18 : * mucpri, corrsp, xty * 19 : * _______________________________________________________ * 20 : * * 21 : * EN ENTREE : * 22 : * MRIGID : Matrice a projeter * 23 : * MTAB1 : Base de modes, reels ou complexes * 24 : * 'REEL' : indique que l'on utilise le produit * 25 : * scalaire reel (pas de conjugaison) * 26 : * * 27 : * EN SORTIE : * 28 : * RI1 : Matrice projetee (partie reelle) * 29 : * RI2 : Matrice projetee (partie imaginaire) * 30 : * _______________________________________________________ * 31 : * * 32 : * REMARQUE : * 33 : * L'operation realisee est : * 34 : * (MTAB1)t . MRIGID. MTAB1 * 35 : * Dans le cas complexe, la transposition est accompagnee * 36 : * de la conjugaison (si REEL n'est pas mentionne). * 37 : *********************************************************** 38 : * 39 : -INC CCOPTIO 40 : -INC SMCHPOI 41 : -INC SMELEME 42 : -INC SMLCHPO 43 : -INC SMRIGID 44 : -INC SMTABLE 45 : -INC SMLMOTS 46 : pointeur IPLMOX.mlmots, IPLMOY.mlmots 47 : * 48 : * Declarations 49 : * 50 : REAL*8 XVAL, RMAX 51 : CHARACTER*8 LETYPE 52 : CHARACTER*8 TYPMOD,CHARRE 53 : LOGICAL MODCOM 54 : INTEGER I, J, NBMOD, POS, IREEL, IVALRE, IOBRE, isym 55 : REAL*8 XVALRE 56 : LOGICAL LOGRE 57 : * 58 : * par defaut, on considere le vrai produit scalaire avec le conjugué 59 : * transposé y^H.M.x , mais avec POS=1, on prend seulement y^T.M.x 60 : IREEL = -1 61 : IF(POS.EQ.1) IREEL = 1 62 : * 63 : * recup de MTYMAT + symetrie / antisymetrie de la matrice a projeter 64 : segact,MRIGID 65 : LETYPE = MRIGID.MTYMAT 66 : NRIGEL=IRIGEL(/2) 67 : isym=IRIGEL(7,1) 68 : if(NRIGEL.le.1) goto 09 69 : do iel=2,NRIGEL 70 : if (isym.ne.IRIGEL(7,iel)) then 71 : isym=2 72 : goto 09 73 : endif 74 : enddo 75 : 09 CONTINUE 76 : segdes,MRIGID 77 : * 78 : LETYPE = ' ' 78 : CALL ACCTAB(MTAB1,'MOT',0,0.0D0,'MODES',.TRUE.,0, 79 : & 'TABLE',IVALRE,XVALRE,CHARRE,LOGRE,MTAB2) 80 : SEGACT, MTAB2 81 : NBMOD = MTAB2.MLOTAB-2 82 : SEGDES, MTAB2 83 : * 84 : N1 = NBMOD 85 : SEGINI, MLCHP1 86 : SEGINI, MLCHP2 87 : * 88 : * Constitution du maillage support et du segment descriptif 89 : * 90 : NBNN = NBMOD 91 : NBELEM = 1 92 : NBSOUS = 0 93 : NBREF = 0 94 : SEGINI, MELEME 95 : ITYPEL = 1 96 : * 97 : NLIGRD=NBMOD 98 : NLIGRP=NBMOD 99 : nelrig=1 100 : SEGINI,DESCR 101 : * 102 : DO 10, I=1, NBMOD 103 : IPT1 = 0 104 : * 104 : CALL ACCTAB(MTAB2,'ENTIER',I,0.0D0,' ',.TRUE.,0, 105 : & 'TABLE',IVALRE,XVALRE,CHARRE,LOGRE,MTAB3) 105 : CALL ACCTAB(MTAB3,'MOT',0,0.0D0,'SOUSTYPE',.TRUE.,0, 106 : & 'MOT',IVALRE,XVALRE,TYPMOD,LOGRE,IOBRE) 107 : * 108 : * Le calcul est impossible : 109 : * 110 : IF (TYPMOD.EQ.'MODE_ANN') THEN 111 : IF (LETYPE.NE.'ANNULE') THEN 112 : WRITE (*,*) 'Calcul impossible : modes annules !!!' 113 : LETYPE = 'ANNULE' 114 : ENDIF 115 : GOTO 5 116 : ENDIF 117 : * 118 : * Cas reel ou cas complexe ? 119 : * 120 : IF (TYPMOD .EQ. 'MODE_COM') THEN 121 : MODCOM=.TRUE. 121 : CALL ACCTAB(MTAB3,'MOT',0,0.0D0,'DEFORMEE_MODALE_REELLE', 122 : & .TRUE.,0,'CHPOINT',IVALRE,XVALRE,CHARRE,LOGRE,MCHPOI) 123 : MLCHP1.ICHPOI(I) = MCHPOI 123 : CALL ACCTAB(MTAB3,'MOT',0,0.0D0,'DEFORMEE_MODALE_IMAGINAIRE', 124 : & .TRUE.,0,'CHPOINT',IVALRE,XVALRE,CHARRE,LOGRE,MCHPOI) 125 : MLCHP2.ICHPOI(I) = MCHPOI 126 : ELSE 127 : MODCOM = .FALSE. 127 : CALL ACCTAB(MTAB3,'MOT',0,0.0D0,'DEFORMEE_MODALE',.TRUE.,0, 128 : & 'CHPOINT',IVALRE,XVALRE,CHARRE,LOGRE,MCHPOI) 129 : MLCHP1.ICHPOI(I) = MCHPOI 130 : ENDIF 131 : * 131 : CALL ACCTAB(MTAB3,'MOT',0,0.0D0,'POINT_REPERE',.TRUE.,0, 132 : & 'POINT',IVALRE,XVALRE,CHARRE,LOGRE,IPT1) 133 : * 134 : 5 MELEME.NUM(I,1)=IPT1 135 : * 136 : DESCR.LISINC(I) = 'ALFA' 137 : DESCR.LISDUA(I) = 'FALF' 138 : DESCR.NOELEP(I) = I 139 : DESCR.NOELED(I) = I 140 : * 141 : 10 CONTINUE 142 : SEGDES, DESCR 143 : SEGDES, MELEME 144 : * 145 : * Constitution des segments XMATRI 146 : * 147 : NLIGRD=NBMOD 148 : NLIGRP=NBMOD 149 : * 150 : IF (LETYPE .EQ. 'ANNULE') THEN 151 : SEGINI, XMATR1 152 : IF (MODCOM) THEN 153 : SEGINI, XMATR2 154 : SEGDES, XMATR2 155 : ENDIF 156 : SEGDES, XMATR1 157 : GOTO 55 158 : ENDIF 159 : * 160 : * Cas reel 161 : * 162 : SEGINI, XMATR1 163 : c DO 20, I=1, NBMOD 164 : c MCHPO1 = MLCHP1.ICHPOI(I) 165 : c DO 20, J=1, NBMOD 166 : c MCHPO2 = MLCHP1.ICHPOI(J) 166 : c CALL YTMX (MCHPO2, MCHPO1, MRIGID, XVAL) 167 : c XMATR1.RE(I,J,1)=XVAL 168 : c 20 CONTINUE 169 : c bp 19/11/2010 : on ameliore les perf en remplacant ci dessus par : 170 : if (isym.eq.0) then 171 : DO 20 J=1, NBMOD 172 : MCHPO2 = MLCHP1.ICHPOI(J) 172 : CALL MUCPRI(MCHPO2,MRIGID,MCHPO3) 173 : DO I=J, NBMOD 174 : MCHPO1 = MLCHP1.ICHPOI(I) 174 : CALL CORRSP(MRIGID,MCHPO1,MCHPO3,IPLMOX,IPLMOY) 174 : CALL XTY1(MCHPO1,MCHPO3,IPLMOX,IPLMOY,XVAL) 175 : XMATR1.RE(I,J,1)=XVAL 176 : XMATR1.RE(J,I,1)=XVAL 177 : ENDDO 178 : segsup,IPLMOX,IPLMOY 179 : segsup,MCHPO3 180 : 20 CONTINUE 181 : elseif (isym.eq.1) then 182 : DO 21 J=1, NBMOD 183 : XMATR1.RE(J,J,1)=0.D0 184 : if(J.ge.NBMOD) goto 21 185 : MCHPO2 = MLCHP1.ICHPOI(J) 185 : CALL MUCPRI(MCHPO2,MRIGID,MCHPO3) 186 : JP1=J+1 187 : DO I=JP1, NBMOD 188 : MCHPO1 = MLCHP1.ICHPOI(I) 188 : CALL CORRSP(MRIGID,MCHPO1,MCHPO3,IPLMOX,IPLMOY) 188 : CALL XTY1(MCHPO1,MCHPO3,IPLMOX,IPLMOY,XVAL) 189 : XMATR1.RE(I,J,1)=XVAL 190 : XMATR1.RE(J,I,1)=-1.D0*XVAL 191 : ENDDO 192 : segsup,IPLMOX,IPLMOY 193 : segsup,MCHPO3 194 : 21 CONTINUE 195 : else 196 : DO 22, J=1, NBMOD 197 : MCHPO2 = MLCHP1.ICHPOI(J) 197 : CALL MUCPRI(MCHPO2,MRIGID,MCHPO3) 198 : DO I=1, NBMOD 199 : MCHPO1 = MLCHP1.ICHPOI(I) 199 : CALL CORRSP(MRIGID,MCHPO1,MCHPO3,IPLMOX,IPLMOY) 199 : CALL XTY1(MCHPO1,MCHPO3,IPLMOX,IPLMOY,XVAL) 200 : XMATR1.RE(I,J,1)=XVAL 201 : ENDDO 202 : segsup,IPLMOX,IPLMOY 203 : segsup,MCHPO3 204 : 22 CONTINUE 205 : endif 206 : * 207 : * Cas complexe : calcul de termes complementaires 208 : * 209 : IF (MODCOM) THEN 210 : c partie reelle = phiR_i^T.M.phiR_j +/- phiI_i^T.M. phiI_j 211 : c DO 30, I=1, NBMOD 212 : c MCHPO1 = MLCHP2.ICHPOI(I) 213 : c DO 30, J=1, NBMOD 214 : c MCHPO2 = MLCHP2.ICHPOI(J) 214 : c CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL) 215 : c XMATR1.RE(I,J,1)=XMATR1.RE(I,J,1)-IREEL*XVAL 216 : c 30 CONTINUE 217 : DO 30, J=1, NBMOD 218 : MCHPO2 = MLCHP2.ICHPOI(J) 218 : CALL MUCPRI(MCHPO2,MRIGID,MCHPO3) 219 : DO I=1, NBMOD 220 : MCHPO1 = MLCHP2.ICHPOI(I) 220 : CALL CORRSP(MRIGID,MCHPO1,MCHPO3,IPLMOX,IPLMOY) 220 : CALL XTY1(MCHPO1,MCHPO3,IPLMOX,IPLMOY,XVAL) 221 : XMATR1.RE(I,J,1)=XMATR1.RE(I,J,1)-IREEL*XVAL 222 : ENDDO 223 : segsup,IPLMOX,IPLMOY 224 : segsup,MCHPO3 225 : 30 CONTINUE 226 : * 2eme matrice (=partie imaginaire) 227 : c partie imaginaire = phiR_i^T.M.phiI_j -/+ phiI_i^T.M. phiR_j 228 : SEGINI,XMATR2 229 : c DO 40, I=1, NBMOD 230 : c MCHPO1 = MLCHP1.ICHPOI(I) 231 : c DO 40, J=1, NBMOD 232 : c MCHPO2 = MLCHP2.ICHPOI(J) 232 : c CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL) 233 : c XMATR2.RE(I,J,1)=XVAL 234 : c 40 CONTINUE 235 : c DO 50, I=1, NBMOD 236 : c MCHPO1 = MLCHP2.ICHPOI(I) 237 : c DO 50, J=1, NBMOD 238 : c MCHPO2 = MLCHP1.ICHPOI(J) 238 : c CALL YTMX (MCHPO1, MCHPO2, MRIGID, XVAL) 239 : c XMATR2.RE(I,J,1)=XMATR2.RE(I,J,1)+IREEL*XVAL 240 : c 50 CONTINUE 241 : DO 40, J=1, NBMOD 242 : MCHPO2 = MLCHP2.ICHPOI(J) 242 : CALL MUCPRI(MCHPO2,MRIGID,MCHPO3) 243 : DO I=1, NBMOD 244 : MCHPO1 = MLCHP1.ICHPOI(I) 244 : CALL CORRSP(MRIGID,MCHPO1,MCHPO3,IPLMOX,IPLMOY) 244 : CALL XTY1(MCHPO1,MCHPO3,IPLMOX,IPLMOY,XVAL) 245 : XMATR2.RE(I,J,1)=XVAL 246 : ENDDO 247 : segsup,IPLMOX,IPLMOY 248 : segsup,MCHPO3 249 : 40 CONTINUE 250 : DO 50, J=1, NBMOD 251 : MCHPO2 = MLCHP1.ICHPOI(J) 251 : CALL MUCPRI(MCHPO2,MRIGID,MCHPO3) 252 : DO I=1, NBMOD 253 : MCHPO1 = MLCHP2.ICHPOI(I) 253 : CALL CORRSP(MRIGID,MCHPO1,MCHPO3,IPLMOX,IPLMOY) 253 : CALL XTY1(MCHPO1,MCHPO3,IPLMOX,IPLMOY,XVAL) 254 : XMATR2.RE(I,J,1)=XMATR2.RE(I,J,1)+IREEL*XVAL 255 : ENDDO 256 : segsup,IPLMOX,IPLMOY 257 : segsup,MCHPO3 258 : 50 CONTINUE 259 : SEGDES,XMATR2 260 : ENDIF 261 : * 262 : SEGDES,XMATR1 263 : SEGDES,MLCHP1,MLCHP2 264 : * 265 : * Creation des segments IMATRI 266 : * 267 : 55 NELRIG = 1 268 : * SEGINI, IMATR1 269 : * IMATR1.IMATTT(1) = XMATR1 270 : * SEGDES, IMATR1 271 : IF (MODCOM) THEN 272 : * SEGINI, IMATR2 273 : * IMATR2.IMATTT(1) = XMATR2 274 : * SEGDES, IMATR2 275 : ENDIF 276 : * 277 : * Creation des rigidites calculees 278 : * 279 : NRIGE=7 280 : NRIGEL=1 281 : SEGINI, RI1 282 : RI1.MTYMAT = LETYPE 283 : RI1.IFORIG = IFOMOD 284 : RI1.IMGEO1 = 0 285 : RI1.IMGEO2 = 0 286 : RI1.COERIG = 1.D0 287 : RI1.IRIGEL(1,1) = MELEME 288 : RI1.IRIGEL(2,1) = 0 289 : RI1.IRIGEL(3,1) = DESCR 290 : RI1.IRIGEL(4,1) = xMATR1 291 : RI1.IRIGEL(5,1) = NIFOUR 292 : RI1.IRIGEL(6,1) = 0 293 : c RI1.IRIGEL(7,1) = 2 294 : RI1.IRIGEL(7,1) = isym 295 : SEGDES, RI1 296 : IF (MODCOM) THEN 297 : RI1.IRIGEL(7,1) = 2 298 : SEGINI, RI2 = RI1 299 : RI2.IRIGEL(4,1) = xMATR2 300 : SEGDES, RI2 301 : ELSE 302 : RI2 = 0 303 : SEGSUP, MLCHP2 304 : ENDIF 305 : * 306 : RETURN 307 : END 308 : 309 : 310 : 311 : 312 : 313 : 314 : 315 : 316 : 317 : 318 : 319 :
© Cast3M 2003 - All rights reserved.
Disclaimer