1 : C CCGNOR SOURCE GOUNAND 11/04/29 21:15:18 6947 2 : SUBROUTINE CCGNOR(LCOF,NOMLOI, 3 : $ FC, 4 : $ IMPR,IRET) 5 : IMPLICIT REAL*8 (A-H,O-Z) 6 : IMPLICIT INTEGER (I-N) 7 : C*********************************************************************** 8 : C NOM : CCGNOR 9 : C DESCRIPTION : Calcul des composantes d'un vecteur normal 10 : C 11 : C 12 : C LANGAGE : ESOPE 13 : C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF) 14 : C mél : gounand@semt2.smts.cea.fr 15 : C*********************************************************************** 16 : C APPELES : 17 : C APPELE PAR : 18 : C*********************************************************************** 19 : C ENTREES : 20 : C ENTREES/SORTIES : 21 : C SORTIES : - 22 : C TRAVAIL : 23 : C*********************************************************************** 24 : C VERSION : v1, 09/03/07, version initiale 25 : C HISTORIQUE : v1, 09/03/07, création 26 : C HISTORIQUE : 27 : C HISTORIQUE : 28 : C*********************************************************************** 29 : C Prière de PRENDRE LE TEMPS de compléter les commentaires 30 : C en cas de modification de ce sous-programme afin de faciliter 31 : C la maintenance ! 32 : C*********************************************************************** 33 : -INC CCOPTIO 34 : -INC CCREEL 35 : CBEGININCLUDE SMCHAEL 36 : SEGMENT MCHAEL 37 : POINTEUR IMACHE(N1).MELEME 38 : POINTEUR ICHEVA(N1).MCHEVA 39 : ENDSEGMENT 40 : SEGMENT MCHEVA 41 : REAL*8 VELCHE(NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM) 42 : ENDSEGMENT 43 : SEGMENT LCHEVA 44 : POINTEUR LISCHE(NBCHE).MCHEVA 45 : ENDSEGMENT 46 : CENDINCLUDE SMCHAEL 47 : INTEGER NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM,N1 48 : POINTEUR FC.MCHEVA 49 : POINTEUR LCOF.LCHEVA 50 : POINTEUR JMAJAC.MCHEVA 51 : POINTEUR JMIJAC.MCHEVA 52 : POINTEUR JDTJAC.MCHEVA 53 : CHARACTER*8 NOMLOI 54 : INTEGER ICOF 55 : * 56 : -INC TMXMAT 57 : * Objets temporaires 58 : POINTEUR JAC.MXMAT,JT.MXMAT,JP.MXMAT 59 : POINTEUR G.MXMAT,IG.MXMAT 60 : * 61 : SEGMENT MIMAT2 62 : INTEGER IMAT2(2,2) 63 : ENDSEGMENT 64 : POINTEUR EIJ.MIMAT2 65 : SEGMENT MIMAT3 66 : INTEGER IMAT3(3,3,3) 67 : ENDSEGMENT 68 : POINTEUR EIJK.MIMAT3 69 : * 70 : LOGICAL LBID 71 : INTEGER LAXSP 72 : * 73 : INTEGER IMPR,IRET 74 : * 75 : * Executable statements 76 : * 77 : * IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans ccgnor' 78 : * WRITE(IOIMP,*) 'Entrée dans ccgnor' 79 : * WRITE(IOIMP,*) 'NOMLOI=',NOMLOI 80 : C IF (.NOT.(IDIM.EQ.1)) THEN 81 : C WRITE(IOIMP,*) 'IDIM=',IDIM,' ?' 82 : C GOTO 9999 83 : C ENDIF 84 : NLFC=FC.VELCHE(/6) 85 : NPFC=FC.VELCHE(/5) 86 : ICOF=0 87 : * 88 : ICOF=ICOF+1 89 : JMAJAC=LCOF.LISCHE(ICOF) 90 : NLJA=JMAJAC.VELCHE(/6) 91 : NPJA=JMAJAC.VELCHE(/5) 92 : IREF=JMAJAC.VELCHE(/4) 93 : IREL=JMAJAC.VELCHE(/3) 94 : * SEGPRT,JMAJAC 95 : * 96 : * WRITE(IOIMP,*) 'IREL=',IREL 97 : * WRITE(IOIMP,*) 'IREF=',IREF 98 : IF (IREL.NE.IDIM) THEN 99 : WRITE(IOIMP,*) 'Erreur dims JMAJAC' 100 : GOTO 9999 101 : ENDIF 102 : IF (IREL.NE.IREF+1) THEN 103 : WRITE(IOIMP,*) 'Le maillage donne nest pas une surface' 104 : GOTO 9999 105 : ENDIF 106 : IF ((IREL.NE.2).AND.(IREL.NE.3)) THEN 107 : WRITE(IOIMP,*) 'Ne marche quen dimension despace 2 ou 3' 108 : GOTO 9999 109 : ENDIF 110 : * 111 : * Objets temporaires et à préconditionner 112 : * 113 : LDIM1=IREL 114 : LDIM2=IREF 115 : SEGINI,JAC 116 : LDIM1=IREF 117 : LDIM2=IREL 118 : SEGINI,JT 119 : * SEGINI,JP 120 : LDIM1=IREF 121 : LDIM2=IREF 122 : SEGINI,G 123 : SEGINI,IG 124 : * 125 : * Initialisation des tenseurs de permutation 126 : * 127 : SEGINI,EIJ 128 : IMULT=1 129 : ICPT=0 130 : DO I=1,2 131 : DO J=1,2 132 : IF (I.NE.J) THEN 133 : ICPT=ICPT+1 134 : IF (ICPT.EQ.2) THEN 135 : ICPT=0 136 : IMULT=IMULT*(-1) 137 : ENDIF 138 : EIJ.IMAT2(I,J)=IMULT 139 : ENDIF 140 : ENDDO 141 : ENDDO 142 : SEGINI,EIJK 143 : IMULT=1 144 : ICPT=0 145 : DO I=1,3 146 : DO J=1,3 147 : IF (I.NE.J) THEN 148 : DO K=1,3 149 : IF ((K.NE.I).AND.(K.NE.J)) THEN 150 : ICPT=ICPT+1 151 : IF (ICPT.EQ.2) THEN 152 : ICPT=0 153 : IMULT=IMULT*(-1) 154 : ENDIF 155 : EIJK.IMAT3(I,J,K)=IMULT 156 : ENDIF 157 : ENDDO 158 : ENDIF 159 : ENDDO 160 : ENDDO 161 : * SEGPRT,EIJ 162 : * SEGPRT,EIJK 163 : * 164 : DO ILFC=1,NLFC 165 : IF (NLJA.EQ.1) THEN 166 : ILJA=1 167 : ELSE 168 : ILJA=ILFC 169 : ENDIF 170 : DO IPFC=1,NPFC 171 : IF (NPJA.EQ.1) THEN 172 : IPJA=1 173 : ELSE 174 : IPJA=IPFC 175 : ENDIF 176 : * 177 : * Copie du jacobien 178 : * 178 : CALL MAMA(JMAJAC.VELCHE(1,1,1,1,IPJA,ILJA),IREL,IREF, 179 : $ 'COPIE ', 180 : $ JAC.XMAT,IREL,IREF, 181 : $ IMPR,IRET) 182 : IF (IRET.NE.0) GOTO 9999 183 : * Calcul de Jt 183 : CALL MAMA(JAC.XMAT,IREL,IREF, 184 : $ 'TRANSPOS',JT.XMAT,IREF,IREL, 185 : $ IMPR,IRET) 186 : IF (IRET.NE.0) GOTO 9999 187 : * Calcul de G=JtJ 187 : CALL MAMAMA(JT.XMAT,IREF,IREL,JAC.XMAT,IREL,IREF, 188 : $ 'FOIS ',G.XMAT,IREF,IREF, 189 : $ IMPR,IRET) 190 : IF (IRET.NE.0) GOTO 9999 191 : * Calcul de l'inverse, du déterminant et trace de l'inverse de g 192 : LBID=.FALSE. 192 : CALL GEOLI2(IREF,1,1,G.XMAT,IG.XMAT,DETG,LBID,IMPR,IRET) 193 : IF (IRET.NE.0) GOTO 9999 194 : XM=SQRT(DETG) 195 : ** 196 : ** Calcul de la pseudo-inverse J+ = g-1 Jt 197 : ** 197 : * CALL MAMAMA(IG.XMAT,IREF,IREF,JT.XMAT,IREF,IREL, 198 : * $ 'FOIS ',JP.XMAT,IREF,IREL, 199 : * $ IMPR,IRET) 200 : * IF (IRET.NE.0) GOTO 9999 201 : 202 : IF (NOMLOI(1:4).EQ.'VNOR') THEN 202 : CALL CH2INT(NOMLOI(5:5),I,IMPR,IRET) 203 : IF (IRET.NE.0) GOTO 9999 204 : ** 205 : ** Calcul de la pseudo-inverse J+ = g-1 Jt 206 : ** 206 : * CALL MAMAMA(IG.XMAT,IREF,IREF,JT.XMAT,IREF,IREL, 207 : * $ 'FOIS ',JP.XMAT,IREF,IREL, 208 : * $ IMPR,IRET) 209 : * IF (IRET.NE.0) GOTO 9999 210 : 211 : * 212 : IF (IREL.EQ.2) THEN 213 : CONTRI=0.D0 214 : DO J=1,2 215 : CONTRI=CONTRI-(EIJ.IMAT2(I,J)* 216 : $ JMAJAC.VELCHE(1,1,J,1,IPJA,ILJA)) 217 : ENDDO 218 : CONTRI=CONTRI/XM 219 : ELSEIF (IREL.EQ.3) THEN 220 : C XNUM=XM 221 : C XDENO=0.D0 222 : C* C'est louche parce que II ne varie pas de 1 à IREL 223 : C DO II=1,IREF 224 : C DO IJ=1,IREF 225 : C* DO IN=1,IREF 226 : C* DO IO=1,IREF 227 : C XDENO=XDENO+ 228 : C $ (EIJ.IMAT2(II,IJ)* 229 : C $ JAC.XMAT(II,1)*JAC.XMAT(IJ,2))**2 230 : CC $ (EIJ.IMAT2(II,IJ)*EIJ.IMAT2(IN,IO)* 231 : CC $ JAC.XMAT(II,1)*JAC.XMAT(IJ,2)* 232 : CC $ JAC.XMAT(IN,1)*JAC.XMAT(IO,2)) 233 : CC ENDDO 234 : CC ENDDO 235 : C ENDDO 236 : C ENDDO 237 : C XAL=XNUM/XDENO 238 : C CONTRI=0.D0 239 : C DO J=1,IREL 240 : C DO K=1,IREL 241 : C CONTRI=CONTRI+ 242 : C $ (EIJK.IMAT3(I,J,K)* 243 : C $ JAC.XMAT(J,1)*JAC.XMAT(K,2)) 244 : C ENDDO 245 : C ENDDO 246 : C CONTRI=CONTRI*XAL 247 : CONTRI=0.D0 248 : DO IA=1,IREF 249 : DO IB=1,IREF 250 : DO J=1,IREL 251 : DO K=1,IREL 252 : CONTRI=CONTRI+(EIJ.IMAT2(IA,IB)* 253 : $ EIJK.IMAT3(I,J,K)* 254 : $ JMAJAC.VELCHE(1,1,J,IA,IPJA,ILJA)* 255 : $ JMAJAC.VELCHE(1,1,K,IB,IPJA,ILJA)) 256 : C $ JP.XMAT(IA,J)*JP.XMAT(IB,K)) 257 : ENDDO 258 : ENDDO 259 : ENDDO 260 : ENDDO 261 : * SEGPRT,JAC 262 : * SEGPRT,JP 263 : * WRITE(IOIMP,*) 'XM=',XM 264 : CONTRI=CONTRI/(2.D0*XM) 265 : * WRITE(IOIMP,*) 'CONTRI=',CONTRI 266 : ELSE 267 : WRITE(IOIMP,*) 'Erreur grave IREL=',IREL 268 : GOTO 9999 269 : ENDIF 270 : ELSEIF (NOMLOI(1:4).EQ.'VNOJ') THEN 270 : CALL CH2INT(NOMLOI(5:5),I,IMPR,IRET) 271 : IF (IRET.NE.0) GOTO 9999 271 : CALL CH2INT(NOMLOI(6:6),J,IMPR,IRET) 272 : IF (IRET.NE.0) GOTO 9999 272 : CALL CH2INT(NOMLOI(7:7),L,IMPR,IRET) 273 : IF (IRET.NE.0) GOTO 9999 274 : IF (IREL.EQ.2) THEN 275 : CONTRI=0.D0 276 : CONTRI=CONTRI-(EIJ.IMAT2(I,J)* 277 : $ JMAJAC.VELCHE(1,1,L,1,IPJA,ILJA)) 278 : CONTRI=CONTRI/XM 279 : ELSEIF (IREL.EQ.3) THEN 280 : CONTRI=0.D0 281 : DO IA=1,IREF 282 : DO IB=1,IREF 283 : * DO J=1,IREL 284 : DO K=1,IREL 285 : CONTRI=CONTRI+(EIJ.IMAT2(IA,IB)* 286 : $ EIJK.IMAT3(I,J,K)* 287 : $ JMAJAC.VELCHE(1,1,L,IA,IPJA,ILJA)* 288 : $ JMAJAC.VELCHE(1,1,K,IB,IPJA,ILJA)) 289 : CONTRI=CONTRI+(EIJ.IMAT2(IA,IB)* 290 : $ EIJK.IMAT3(I,K,J)* 291 : $ JMAJAC.VELCHE(1,1,K,IA,IPJA,ILJA)* 292 : $ JMAJAC.VELCHE(1,1,L,IB,IPJA,ILJA)) 293 : ENDDO 294 : * ENDDO 295 : ENDDO 296 : ENDDO 297 : * SEGPRT,JAC 298 : * SEGPRT,JP 299 : * WRITE(IOIMP,*) 'XM=',XM 300 : CONTRI=CONTRI/(2.D0*XM) 301 : * WRITE(IOIMP,*) 'CONTRI=',CONTRI 302 : ELSE 303 : WRITE(IOIMP,*) 'Erreur grave IREL=',IREL 304 : GOTO 9999 305 : ENDIF 306 : ELSE 307 : WRITE(IOIMP,*) 'Erreur grave NOMLOI=',NOMLOI 308 : GOTO 9999 309 : ENDIF 310 : FC.VELCHE(1,1,1,1,IPFC,ILFC)= 311 : $ FC.VELCHE(1,1,1,1,IPFC,ILFC)+ 312 : $ CONTRI 313 : ENDDO 314 : ENDDO 315 : SEGSUP,EIJK 316 : SEGSUP,EIJ 317 : * SEGPRT,FC 318 : SEGSUP,JAC 319 : * SEGSUP,JP 320 : SEGSUP,JT 321 : SEGSUP,G 322 : SEGSUP,IG 323 : * 324 : * Normal termination 325 : * 326 : IRET=0 327 : RETURN 328 : * 329 : * Format handling 330 : * 331 : * 332 : * Error handling 333 : * 334 : 9999 CONTINUE 335 : IRET=1 336 : WRITE(IOIMP,*) 'An error was detected in subroutine ccgnor' 337 : RETURN 338 : * 339 : * End of subroutine CCGNOR 340 : * 341 : END 342 : 343 : 344 : 345 : 346 : 347 : 348 :
© Cast3M 2003 - All rights reserved.
Disclaimer