1 : C TRJRND SOURCE CHAT 05/01/13 03:51:03 5004 2 : SUBROUTINE TRJRND(IZN3,IZVIT,MELEME,TTEMP,NOE,IEL,NDIM,ITER, 3 : * IELL2,XARI,DTI,ICONT,INOEU) 4 : C 5 : ******************************************************************* 6 : * 7 : * CE SOUS PROGRAMME CHERCHE UN VOISIN ET CALCULE LA TRAJECTOIRE 8 : * LORSQU'ON TOMBE SUR UN NOEUD. 9 : * ENTREES: 10 : * IZN3 SEGMENT QUI STOCKE LES RESULTATS 11 : * IZVIT SEGMENT QUI CONTIENT LES INFOS SUR LES FLUX 12 : * MELEME MAILLAGE 13 : * TTEMP TEMPS COURANT PLUS LE PAS DE LA MAILLE IEL1 14 : * INOEU NOEUD PAR LEQUEL ON SORT DE IEL1 15 : * IEL1 ELEMENT DUQUEL ON SORT 16 : * NDIM 17 : * SORTIES: 18 : * IEL2 ELEMENT DANS LEQUEL PASSE LA TRAJECTOIRE 19 : * XARI POSITION DE SORTIE DE L'ELEMENT 20 : * DTINT TEMPS MIS POUR PARCOURIR IEL2 21 : * ICONT FACE DE SORTIE DE IEL2 22 : * 23 : ******************************************************************* 24 : C 25 : IMPLICIT INTEGER(I-N) 26 : IMPLICIT REAL*8 (A-H,O-Z) 27 : C 28 : -INC CCOPTIO 29 : -INC SMCOORD 30 : -INC SMELEME 31 : -INC SMCHPOI 32 : C 33 : SEGMENT IZN3 34 : INTEGER NAPAR3(NPOS),NUM3(NPOS) 35 : REAL*8 CREF3(NDIM,NPOS),TPAR(NPOS) 36 : ENDSEGMENT 37 : C 38 : SEGMENT IZVIT 39 : REAL*8 TEMTRA(NVIPT) 40 : INTEGER IPUN(NBS),IDUN(NBS),IPVPT(NVIPT),IFORML 41 : ENDSEGMENT 42 : C IDUN(I) nombre d elements avant le sous maillage I 43 : C IPVPT pointeurs de izvpt pour chaque pas de temps 44 : C 45 : SEGMENT IZUN 46 : REAL*8 UN(I1,I2,I3) 47 : ENDSEGMENT 48 : POINTEUR IZUN1.IZUN,IZUN2.IZUN 49 : C 50 : DIMENSION XARI(3),XDEP(3),X(4) 51 : C 52 : DIMENSION XTRI3(2,3),XQUA4(2,4),XTET4(3,4) 53 : DIMENSION XPRI6(3,6),XCUB8(3,8) 54 : C 55 : DATA XTRI3 /0.D0,0.D0, 1.D0,0.D0, 0.D0,1.D0/ 56 : C 57 : DATA XQUA4 /-1.D0,-1.D0, 1.D0,-1.D0, 1.D0,1.D0, -1.D0,1.D0/ 58 : C 59 : DATA XTET4 /0.D0,0.D0,0.D0, 1.D0,0.D0,0.D0, 0.D0,1.D0,0.D0, 60 : * 0.D0,0.D0,1.D0 / 61 : C 62 : DATA XCUB8 /-1.D0,-1.D0,-1.D0, 1.D0,-1.D0,-1.D0, 1.D0,1.D0,-1.D0, 63 : * -1.D0,1.D0,-1.D0, -1.D0,-1.D0,1.D0, 1.D0,-1.D0,1.D0, 64 : * 1.D0,1.D0,1.D0, -1.D0,1.D0,1.D0/ 65 : C 66 : DATA XPRI6 /0.D0,0.D0,-1.D0, 1.D0,0.D0,-1.D0, 0.D0,1.D0,-1.D0, 67 : * 0.D0,0.D0,1.D0, 1.D0,0.D0,1.D0, 0.D0,1.D0,1.D0/ 68 : 69 : C 70 : C*** RATTRAPAGE 71 : C ON VA CHERCHER DANS LE MAILLAGE LES ELEMENTS QUI ONT UN 72 : C NOEUD EN COMMUN AVEC L'ELEMENT DE CALCUL 73 : C 74 : C*** ON CHERCHE DANS LE MAILLAGE UN ELEMENT IEL2 AYANT LE NOEUD 75 : C D'ENTREE EN COMMUN AVEC LE NOEUD DE SORTIE DE IEL1 76 : C 77 : C INITIALISATION 78 : NPAPAR=0 79 : IEL2=0 80 : IND=0 81 : IVPT=1 82 : ITEST=4 83 : NPOS=CREF3(/2) 84 : C 85 : C*** ON ACTIVE LE MAILLAGE ET ON FAIT UNE RECHERCHE PAR SOUS MAILLAGE 86 : C IPT5 EST LE SOUS DOMAINE CONTENANT LA MAILLE PAR LAQUELLE ON 87 : C SORT PAR UN NOEUD 88 : C IPT1 REPRESENTE LE SOUS MAILLAGE DANS LEQUEL ON RECHERCHE UN VOISIN 89 : C 90 : SEGACT MELEME 91 : NBSOUS=LISOUS(/1) 92 : NBS=NBSOUS 93 : IF(NBSOUS.EQ.0) NBS=1 94 : IPT1=MELEME 94 : CALL MELNEL(IEL,MELEME,IPT5,NEL5,1) 95 : SEGACT IPT5 96 : IELE=IEL-NEL5 97 : DO 100 ISOUS=1,NBS 98 : IF(NBSOUS.GT.0)THEN 99 : IPT1=LISOUS(ISOUS) 100 : ENDIF 101 : SEGACT IPT1 102 : NEL=IPT1.NUM(/2) 103 : NBN=IPT1.NUM(/1) 104 : C 105 : C POUR CHAQUE NOEUD DE CHAQUE ELEMENT ON REGARDE SI CE NOEUD 106 : C CORRESPOND AU NOEUD DE L'ELEMENT COURANT PAR LEQUEL PASSE 107 : C LA TRAJECTOIRE.DES QU'ON TROUVE UN ELEMENT AYANT LE NOEUD 108 : C EN COMMUN ON CALCULE LA TRAJECTOIRE DE LA PARTICULE DANS CET 109 : C ELEMENT. 110 : C 111 : C 112 : C BOUCLE SUR LES ELEMENTS 113 : DO 1 IEL2=1,NEL 114 : C ON VA TESTER LES ELEMENTS DIFFERENTS DE L'ELEMENT COURANT 115 : IF((IEL2.NE.IELE).OR.(IPT1.NE.IPT5))THEN 116 : C BOUCLE SUR LES NOEUDS 117 : DO 2 INO=1,NBN 118 : C TEST DE COMPARAISON DES NOEUDS 119 : IF((IPT1.NUM(INO,IEL2)).EQ.(IPT5.NUM(NOE,IELE)))THEN 120 : IELL2=IEL2+IND 121 : ITEST=0 122 : ITYG=NUMGEO(IPT1.ITYPEL) 123 : C 124 : C ON DETERMINE LES COORDONNES DU POINT D'ENTREE DE LA PARTICULE 125 : C A PARTIR DU NOEUD TROUVE.ON CALCULE ENSUITE LA TRAJECTOIRE ET 126 : C ON TESTE SI LES COORDONNEES DE SORTIE SONT VALABLES 127 : C 127 : CALL TRJVEL(IZVIT,IZUN,IELL2,IVPT,TTEMP) 128 : C 129 : C ON TESTE LES TRIANGLES 130 : C 131 : IF(ITYG.EQ.4)THEN 132 : DO 7 I=1,2 133 : XDEP(I)=XTRI3(I,INO) 134 : 7 CONTINUE 134 : CALL TRJMA2(XARI,XDEP,UN(1,1,IEL2),DTI,2, 135 : * ICONT,4,IART,INOEU) 136 : XARI(3)=1.D0-XARI(1)-XARI(2) 137 : DO 5 I=1,3 138 : IF(XARI(I)*(1.D0-XARI(I)).GE.0.D0)ITEST=ITEST+1 139 : 5 CONTINUE 140 : ITEST=ITEST-1 141 : C 142 : ELSEIF(ITYG.EQ.8)THEN 143 : C 144 : C ON TESTE LES QUADRANGLES 145 : C 146 : DO 8 I=1,2 147 : XDEP(I)=XQUA4(I,INO) 148 : 8 CONTINUE 148 : CALL TRJMA2(XARI,XDEP,UN(1,1,IEL2),DTI,2, 149 : * ICONT,8,IART,INOEU) 150 : DO 3 I=1,NDIM 151 : XY=(1.D0+XARI(I))*(1.D0-XARI(I)) 152 : IF(XY.GE.0.D0)ITEST=ITEST+1 153 : 3 CONTINUE 154 : C 155 : C ON EST EN DIMENSION 3 ON TESTE LES CUBES 156 : C 157 : ELSEIF(ITYG.EQ.14)THEN 158 : DO 9 I=1,3 159 : XDEP(I)=XCUB8(I,INO) 160 : 9 CONTINUE 160 : CALL TRJMA2(XARI,XDEP,UN(1,1,IEL2),DTI,3, 161 : * ICONT,14,IART,INOEU) 162 : DO 13 I=1,NDIM 163 : XY=(1.D0+XARI(I))*(1.D0-XARI(I)) 164 : IF(XY.GE.0.D0)ITEST=ITEST+1 165 : 13 CONTINUE 166 : C 167 : C ON TESTE LES PRISMES 168 : ELSEIF(ITYG.EQ.16)THEN 169 : C 170 : DO 10 I=1,3 171 : XDEP(I)=XPRI6(I,INO) 172 : 10 CONTINUE 172 : CALL TRJMA2(XARI,XDEP,UN(1,1,IEL2),DTI,3, 173 : * ICONT,16,IART,INOEU) 174 : X(1)=XARI(1) 175 : X(2)=XARI(2) 176 : X(3)=1.D0-X(1)-X(2) 177 : X(4)=XARI(3) 178 : DO 14 I=1,3 179 : IF(X(I)*(1.D0-X(I)).GE.0.D0)ITEST=ITEST+1 180 : 14 CONTINUE 181 : IF((X(4)+1.D0)*(1.D0-X(4)).GE.0.D0)ITEST=ITEST+1 182 : ITEST=ITEST-1 183 : C 184 : C ON TESTE LES TETRAEDRES 185 : ELSEIF(ITYG.EQ.23)THEN 186 : C 187 : DO 11 I=1,3 188 : XDEP(I)=XTET4(I,INO) 189 : 11 CONTINUE 189 : CALL TRJMA2(XARI,XDEP,UN(1,1,IEL2),DTI,3, 190 : * ICONT,23,IART,INOEU) 191 : DO 20 I=1,3 192 : IF(XARI(I)*(1.D0-XARI(I)).GE.0.D0)ITEST=ITEST+1 193 : 20 CONTINUE 194 : XXX=XARI(1)+XARI(2)+XARI(3) 195 : IF(XXX*(1.D0-XXX).GE.0.D0)ITEST=ITEST+1 196 : ITEST=ITEST-1 197 : C 198 : ENDIF 199 : C 200 : C UNE FOIS QUE L'ON A CALCULER LA TRAJECTOIRE POUR UN ELEMENT 201 : C ON TESTE SI LA TRAJECTOIRE APPARTIENT A CET ELEMENT 202 : C DANS CE CAS ON STOCKE LES VALEURS A CONSERVER 203 : C 204 : IF((ITEST.EQ.NDIM).AND.(INO.NE.INOEU))THEN 205 : TPAR(ITER)=TTEMP 206 : DO 4 ID1=1,NDIM 207 : CREF3(ID1,ITER)=XDEP(ID1) 208 : 4 CONTINUE 209 : NAPAR3(ITER)=IELL2 210 : IF ((NPOS-ITER).LE.1) THEN 211 : NPOS=NPOS+50 212 : SEGADJ IZN3 213 : ITER=ITER+1 214 : ENDIF 215 : ITER=ITER+1 216 : C SI ON A TROUVE LE BON ELEMENT ON SORT 217 : GOTO 17 218 : ENDIF 219 : C 220 : C SINON ON RECOMMENCE AVEC UN AUTRE NOEUD OU UN AUTRE ELEMENT 221 : C 222 : ENDIF 223 : 2 CONTINUE 224 : ENDIF 225 : 1 CONTINUE 226 : IND=IND+NEL 227 : IF(IPT1.NE.IPT5)SEGDES IPT1 228 : 100 CONTINUE 229 : SEGDES IPT5 230 : C 231 : C CAS OU ON ARRIVE EN BORD DE DOMAINE 232 : C SI ON TROUVE PAS DE VOISIN ON RETOURNE DANS TRJAVA 233 : DTI=0.D0 234 : IF(ITEST.EQ.4)GOTO 17 235 : C CAS OU LA TRAJECTOIRE NE PASSE PAS PAR LE VOISIN TROUVE 236 : DO 18 I=1,NDIM 237 : XARI(I)=XDEP(I) 238 : 18 CONTINUE 239 : 17 CONTINUE 240 : RETURN 241 : END 242 : 243 : 244 : 245 : 246 : 247 : 248 : 249 :
© Cast3M 2003 - All rights reserved.
Disclaimer