1 : C TRANSF SOURCE PV 07/10/05 21:16:58 5890 2 : C CE SOUS PROGRAMME MAILLE A L'INTERIEUR D'UN CONTOUR FERME 3 : C IL EST TIRE DE COCO MAIS UN PEU SIMPLIFIE L'ARCHITECTURE 4 : C S'Y PRETANT MIEUX 5 : C DECEMBRE 1982 : PRISE EN COMPTE DE LA NATURE DES CONTOURS 6 : C (INTERNE OU EXTERNE) 7 : C TEST SI CERTAINS SEGMENTS CREES COUPENT UN CONTOUR 8 : C FEVRIER 1983 ==> QUADRANGLES 9 : C 10 : SUBROUTINE TRANSF(NFI,NFMAX,MAI,MAIMAX,ITOUR,X,MAXPTS,NUM, 11 : # NBNN,NBELEM,NUMELG,NUMNP,XAUX,NUMINI,ICLE,QUAL,INAT,IREGU,IMOYE, 12 : # XMOY,NBMN,IRECHA,ichp) 13 : IMPLICIT INTEGER(I-N) 14 : IMPLICIT REAL*8 (A-H,O-Z) 15 : DIMENSION NFI(1),MAI(1),X(3,1),XAUX(2,1),NUM(NBNN,1),INAT(1) 16 : -INC CCOPTIO 17 : -INC CCREEL 18 : MAXELE=NBELEM 19 : NUMELG=0 20 : NUMNP=MAI(ITOUR+1) 21 : NUMINI=NUMNP+1 22 : IZERO=0 23 : IN=MAI(1)+1 24 : LEMIN=NUMELG+1 25 : KT36=ITOUR 26 : LONG=MAI(ITOUR+1) 27 : C ON FORCE LA DENSITE DE CHAQUE POINT A SA VALEUR REELLE CONNUE 28 : DO 300 IT=1,ITOUR 29 : IT1=MAI(IT-1+1)+1 30 : IT2=MAI(IT+1) 31 : DO 300 I=IT1,IT2 32 : ISUIV=I+1 33 : IF (ISUIV.GT.IT2) ISUIV=IT1 34 : IALTE=0 35 : IPREC=I-1 36 : IF (IPREC.LT.IT1) IPREC=IT2 37 : SUIV=SQRT((X(1,NFI(ISUIV))-X(1,NFI(I)))**2+ 38 : # (X(2,NFI(ISUIV))-X(2,NFI(I)))**2) 39 : SPRE=SQRT((X(1,NFI(IPREC))-X(1,NFI(I)))**2+ 40 : # (X(2,NFI(IPREC))-X(2,NFI(I)))**2) 41 : 300 X(3,NFI(I))=SQRT(SUIV*SPRE) 42 : IPARC=0 43 : 200 IPARC=IPARC+1 44 : C VERIFICATION QUE TOUT EST CORRECT : LA SURFACE EST POSITIVE 45 : SURF=0.D0 46 : DO 3500 IT=1,ITOUR 47 : IT1=MAI(IT-1+1)+1 48 : IT2=MAI(IT+1) 49 : SURF1=0.D0 50 : DO 3501 I=IT1,IT2 51 : ISUIV=I+1 52 : IF (ISUIV.GT.IT2) ISUIV=IT1 53 : SURF1=SURF1+X(1,NFI(I))*X(2,NFI(ISUIV))-X(2,NFI(I))*X(1,NFI(ISUIV 54 : # )) 55 : 3501 CONTINUE 56 : INAT(IT)=1 57 : IF (SURF1.LT.0.D0) INAT(IT)=-1 58 : IF (SURF1.EQ.0.D0) GOTO 1500 59 : SURF=SURF+SURF1 60 : 3500 CONTINUE 61 : IF (SURF.LT.0.D0) GOTO 1500 62 : C FABRICATION DES NOEUDS DECALES 63 : IF=MAI(ITOUR+1) 64 : IF (IMOYE.NE.0) GOTO 935 65 : XMOY=0 66 : DO 310 I=IN,IF 67 : 310 XMOY=XMOY+LOG(X(3,NFI(I))) 68 : XMOY=EXP(XMOY/(IF-IN+1)) 69 : IF (IIMPI.EQ.1) WRITE (IOIMP,8877) XMOY 70 : 8877 FORMAT(' LA VALEUR MOYENNE DE LA MAILLE RETENUE EST: ',G12.5) 71 : 935 CONTINUE 72 : XMOMA=0 73 : XMOMI=XGRAND 74 : DO 3001 I=IN,IF 75 : XMOMI=MIN(XMOMI,X(3,NFI(I))) 76 : XMOMA=MAX(XMOMA,X(3,NFI(I))) 77 : 3001 CONTINUE 78 : COS=1.d0/REAL(IREGU) 79 : DO 311 I=IN,IF 80 : * modif pour prendre en compte un chpt de densite 81 : if (ichp.ne.0) xmoy=denchp(ichp,x(1,nfi(i)),x(2,nfi(i))) 82 : XSU=XMOY+XMOY*COS 83 : XIN=XMOY-MIN(XMOMI,XMOY*COS) 84 : IF (X(3,NFI(I)).GT.XIN) GOTO 920 85 : X(3,NFI(I))=X(3,NFI(I))+MIN(COS*X(3,NFI(I)),XMOMI) 86 : GOTO 311 87 : 920 IF (X(3,NFI(I)).LT.XSU) GOTO 921 88 : X(3,NFI(I))=X(3,NFI(I))*(1.d0-COS) 89 : GOTO 311 90 : 921 X(3,NFI(I))=XMOY 91 : 311 CONTINUE 92 : C ON METS UN ELEMENT DANT TOUS LES COINS (JUSQU'A 100 DEG) 93 : 360 CONTINUE 94 : IT=0 95 : 251 IT=IT+1 96 : IF (IT.GT.ITOUR) GOTO 253 97 : IDI=MAI(IT-1+1)+1 98 : IFI=MAI(IT+1) 99 : I=IDI-1 100 : IAUX=0 101 : VECTOR=-XGRAND 102 : 250 I=I+1 103 : IAUX=IAUX+1 104 : IF (IAUX.GT.IFI-IDI+1) GOTO 1200 105 : IF (I.GT.IFI) I=IDI 106 : IF (I.EQ.IDI.AND.IIMPI.EQ.1) WRITE (IOIMP,9878) IT,IDI,IFI,IAUX 107 : 9878 FORMAT (' IT,IDI,IFI,IAUX ',4I5) 108 : IPRES=I-1 109 : IF (IPRES.LT.IDI) IPRES=IFI 110 : ISUIV=I+1 111 : IF (ISUIV.GT.IFI) ISUIV=IDI 112 : IF (NFI(IPRES).EQ.NFI(ISUIV)) GOTO 250 113 : SCAL=(X(1,NFI(IPRES))-X(1,NFI(I)))*(X(1,NFI(ISUIV))-X(1,NFI(I)))+ 114 : # (X(2,NFI(IPRES))-X(2,NFI(I)))*(X(2,NFI(ISUIV))-X(2,NFI(I))) 115 : VECT=(X(1,NFI(IPRES))-X(1,NFI(I)))*(X(2,NFI(ISUIV))-X(2,NFI(I)))- 116 : # (X(2,NFI(IPRES))-X(2,NFI(I)))*(X(1,NFI(ISUIV))-X(1,NFI(I))) 117 : IF (VECT.EQ.0.D0.AND.SCAL.EQ.0.D0) GOTO 1500 118 : ANG=ATAN2(VECT,SCAL) 119 : * write (6,*) 'atan2 - 4 ',ang,vect,vector,iaux,ifi,idi 120 : IF (ANG.LT.-2.34D0.OR.ANG.GE.1.D-5) GOTO 250 121 : IF (ANG.LT.-1.6D0.AND.NBNN.EQ.3) GOTO 250 122 : IF (VECT.LE.VECTOR*0.999999d0.AND.IAUX.NE.IFI-IDI+1.AND.NBNN.EQ.3) 123 : > GOTO 250 124 : * write (6,*) ' apres atan2 -4 ' 125 : LL1=NFI(IPRES) 126 : LL2=NFI(ISUIV) 127 : C ON VA TESTER SI LE SEGMENT EN COUPE UN AUTRE 127 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 128 : IF (IRECL.NE.0.AND.IFI-IDI.GT.2) GOTO 250 129 : VECTSA=VECT 130 : IF (IAUX.EQ.IFI-IDI+1.OR.NBNN.EQ.4) GOTO 390 131 : VECTOR=VECT 132 : IAUX=0 133 : GOTO 250 134 : 390 CONTINUE 135 : C POUR LES QUA4 ON ESSAYE DE COMPLETER LE TRIANGLE EN UN QUADRANGLE 136 : IF (NBNN.EQ.3) GOTO 6100 137 : IF (IFI-IDI.LE.2) GOTO 6100 138 : IF (ANG.GT.-0.8D0) GOTO 6100 139 : IAV=IPRES-1 140 : IF (IAV.LT.IDI) IAV=IFI 141 : C CALCUL DE L'ANGLE 142 : SCAL=(X(1,NFI(IAV))-X(1,NFI(IPRES)))*(X(1,NFI(I))-X(1,NFI(IPRES))) 143 : # +(X(2,NFI(IAV))-X(2,NFI(IPRES)))*(X(2,NFI(I))-X(2,NFI(IPRES))) 144 : VECT=(X(1,NFI(IAV))-X(1,NFI(IPRES)))*(X(2,NFI(I))-X(2,NFI(IPRES))) 145 : # -(X(2,NFI(IAV))-X(2,NFI(IPRES)))*(X(1,NFI(I))-X(1,NFI(IPRES))) 146 : IF (VECT.EQ.0D0.AND.SCAL.EQ.0.D0) GOTO 1500 147 : ANGA=ATAN2(-VECT,SCAL) 148 : * write (6,*) 'atan2 - 5 ',anga 149 : IF (ANGA.GT.0.8D0.OR.ANGA.LE.1.d-5) GOTO 6006 150 : I=IAV 151 : GOTO 250 152 : 6006 IF (ANGA.GT.2.34D0.OR.ANGA.LT.0.D0) GOTO 6010 153 : LL1=NFI(IAV) 154 : LL2=NFI(ISUIV) 154 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 155 : IF (IRECL.NE.0) GOTO 6010 156 : XLA=(X(1,NFI(IAV))-X(1,NFI(ISUIV)))**2 157 : # +(X(2,NFI(IAV))-X(2,NFI(ISUIV)))**2 158 : XCOMP=X(3,NFI(I))*X(3,NFI(IPRES)) 159 : IF (XLA.LT.0.333D0*XCOMP.OR.XLA.GT.3.D0*XCOMP) GOTO 6010 160 : C SI 5 NOEUDS SUR LE CONTOUR TEST SUPPLEMENTAIRE 161 : IF (IFI-IDI.NE.4) GOTO 6008 162 : ICINQ=IAV-1 163 : IF (ICINQ.LT.IDI) ICINQ=IFI 164 : SUC=(X(1,NFI(IAV))-X(1,NFI(ICINQ)))*(X(2,NFI(ISUIV))-X(2,NFI(ICINQ 165 : #)))-(X(2,NFI(IAV))-X(2,NFI(ICINQ)))*(X(1,NFI(ISUIV))-X(1,NFI(ICINQ 166 : #))) 167 : IF (SUC.GT.XLA*XLA/3.D0) GOTO 6008 168 : GOTO 6010 169 : 6008 NUMELG=NUMELG+1 170 : IF (NUMELG.GE.MAXELE) GOTO 1500 171 : NUM(1,NUMELG)=NFI(IAV) 172 : NUM(2,NUMELG)=NFI(IPRES) 173 : NUM(3,NUMELG)=NFI(I) 174 : NUM(4,NUMELG)=NFI(ISUIV) 175 : IF=IF-1 176 : IFI=IFI-2 177 : X(3,NFI(IAV))=SQRT(X(3,NFI(IAV))*X(3,NFI(IPRES))) 178 : X(3,NFI(ISUIV))=SQRT(X(3,NFI(I))*X(3,NFI(ISUIV))) 179 : DO 6002 J=I,IF 180 : NFI(J)=NFI(J+1) 181 : 6002 CONTINUE 182 : IF=IF-1 183 : IF (IPRES.GT.I) IPRES=IPRES-1 184 : DO 6003 J=IPRES,IF 185 : 6003 NFI(J)=NFI(J+1) 186 : DO 6004 J=IT,ITOUR 187 : 6004 MAI(J+1)=MAI(J+1)-2 188 : I=IAV-1 189 : IDI=MAI(IT-1+1)+1 190 : IFI=MAI(IT+1) 191 : VECTOR=-XGRAND 192 : IAUX=0 193 : IF (IFI-IDI.NE.1) GOTO 250 194 : ITOUR=ITOUR-1 195 : IF (ITOUR.EQ.0) GOTO 1600 196 : DO 6005 I=IT,ITOUR 197 : INAT(I)=INAT(I+1) 198 : 6005 MAI(I+1)=MAI(I+1+1)-2 199 : IF=MAI(ITOUR+1) 200 : ID=MAI(IT-1+1)+1 201 : DO 6007 I=ID,IF 202 : 6007 NFI(I)=NFI(I+2) 203 : IT=IT-1 204 : GOTO 251 205 : C CALCUL DE L'ANGLE 206 : 6010 IAP=ISUIV+1 207 : IF (IAP.GT.IFI) IAP=IDI 208 : SCAL=(X(1,NFI(I))-X(1,NFI(ISUIV)))*(X(1,NFI(IAP))-X(1,NFI(ISUIV))) 209 : # +(X(2,NFI(I))-X(2,NFI(ISUIV)))*(X(2,NFI(IAP))-X(2,NFI(ISUIV))) 210 : VECT=(X(1,NFI(I))-X(1,NFI(ISUIV)))*(X(2,NFI(IAP))-X(2,NFI(ISUIV))) 211 : # -(X(2,NFI(I))-X(2,NFI(ISUIV)))*(X(1,NFI(IAP))-X(1,NFI(ISUIV))) 212 : IF (VECT.EQ.0D0.AND.SCAL.EQ.0.D0) GOTO 1500 213 : ANGA=ATAN2(-VECT,SCAL) 214 : * write (6,*) ' atan2 - 1 ',anga 215 : IF (ANGA.GE.0.D0.AND.ANGA.LT.0.8D0) GOTO 250 216 : IF (ANGA.GT.2.34D0.OR.ANGA.LT.0.D0) GOTO 6020 217 : LL1=NFI(IAP) 218 : LL2=NFI(IPRES) 218 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 219 : IF (IRECL.NE.0) GOTO 6020 220 : XLA=(X(1,NFI(IAP))-X(1,NFI(IPRES)))**2 221 : # +(X(2,NFI(IAP))-X(2,NFI(IPRES)))**2 222 : XCOMP=X(3,NFI(I))*X(3,NFI(ISUIV)) 223 : IF (XLA.LT.0.333D0*XCOMP.OR.XLA.GT.3.D0*XCOMP) GOTO 6020 224 : C SI 5 NOEUDS SUR LE CONTOUR TEST SUPPLEMENTAIRE 225 : IF (IFI-IDI.NE.4) GOTO 6018 226 : ICINQ=IPRES-1 227 : IF (ICINQ.LT.IDI) ICINQ=IFI 228 : SUC=(X(1,NFI(IPRES))-X(1,NFI(ICINQ)))*(X(2,NFI(IAP))-X(2,NFI(ICINQ 229 : #)))-(X(2,NFI(IPRES))-X(2,NFI(ICINQ)))*(X(1,NFI(IAP))-X(1,NFI(ICINQ 230 : #))) 231 : IF (SUC.GT.XLA*XLA/3.D0) GOTO 6018 232 : GOTO 6020 233 : 6018 NUMELG=NUMELG+1 234 : IF (NUMELG.GE.MAXELE) GOTO 1500 235 : NUM(1,NUMELG)=NFI(IPRES) 236 : NUM(2,NUMELG)=NFI(I) 237 : NUM(3,NUMELG)=NFI(ISUIV) 238 : NUM(4,NUMELG)=NFI(IAP) 239 : IFI=IFI-2 240 : IF=IF-1 241 : X(3,NFI(IAP))=SQRT(X(3,NFI(IAP))*X(3,NFI(ISUIV))) 242 : X(3,NFI(IPRES))=SQRT(X(3,NFI(IPRES))*X(3,NFI(I))) 243 : DO 6012 J=I,IF 244 : 6012 NFI(J)=NFI(J+1) 245 : IF=IF-1 246 : IF (ISUIV.GT.I) ISUIV=ISUIV-1 247 : DO 6013 J=ISUIV,IF 248 : 6013 NFI(J)=NFI(J+1) 249 : DO 6014 J=IT,ITOUR 250 : 6014 MAI(J+1)=MAI(J+1)-2 251 : I=IAV 252 : IDI=MAI(IT-1+1)+1 253 : IFI=MAI(IT+1) 254 : VECTOR=-XGRAND 255 : IAUX=0 256 : IF (IFI-IDI.NE.1) GOTO 250 257 : ITOUR=ITOUR-1 258 : IF (ITOUR.EQ.0) GOTO 1600 259 : DO 6015 I=IT,ITOUR 260 : INAT(I)=INAT(I+1) 261 : 6015 MAI(I+1)=MAI(I+1+1)-2 262 : IF=MAI(ITOUR+1) 263 : ID=MAI(IT-1+1)+1 264 : DO 6017 I=ID,IF 265 : 6017 NFI(I)=NFI(I+2) 266 : IT=IT-1 267 : GOTO 251 268 : 6020 CONTINUE 269 : C ON TENTE DE FAIRE UN LOSANGE 270 : NTENT=NUMNP+1 271 : IF (NTENT.GT.MAXPTS) GOTO 1500 272 : X(1,NTENT)=X(1,NFI(IPRES))+X(1,NFI(ISUIV))-X(1,NFI(I)) 273 : X(2,NTENT)=X(2,NFI(IPRES))+X(2,NFI(ISUIV))-X(2,NFI(I)) 274 : LL1=NFI(I) 275 : LL2=NTENT 275 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 276 : * write (6,*) ' transa - irecl 1 ',irecl 277 : IF (IRECL.NE.0) GOTO 6100 278 : LL1=NFI(ISUIV) 279 : 1040 IF (IRECL.NE.0) GOTO 6100 280 : C VERIFICATION DES DISTANCES 281 : XLA=(X(1,NTENT)-X(1,NFI(IPRES)))**2 282 : # +(X(2,NTENT)-X(2,NFI(IPRES)))**2 283 : XCOMP=X(3,NFI(IPRES))*X(3,NFI(I)) 284 : IF (XLA.LT.0.333D0*XCOMP.OR.XLA.GT.3.D0*XCOMP) GOTO 6100 285 : XLB=(X(1,NTENT)-X(1,NFI(ISUIV)))**2 286 : # +(X(2,NTENT)-X(2,NFI(ISUIV)))**2 287 : XCOMP=X(3,NFI(ISUIV))*X(3,NFI(I)) 288 : IF (XLB.LT.0.333D0*XCOMP.OR.XLB.GT.3.D0*XCOMP) GOTO 6100 289 : C VERIFICATION : RECHERCHE DU POINT LE PLUS PRES 290 : INATIT=INAT(IT) 291 : XL=XGRAND 292 : IP=0 293 : ITP=0 294 : DO 6030 ITT=1,ITOUR 295 : IF (INATIT.EQ.1.AND.INAT(ITT).EQ.1.AND.IT.NE.ITT) GOTO 6030 296 : IDIA=MAI(ITT-1+1)+1 297 : IFIA=MAI(ITT+1) 298 : DO 6031 IA=IDIA,IFIA 299 : IF (NFI(IA).EQ.NFI(I)) GOTO 6031 300 : IF (NFI(IA).EQ.NFI(IPRES)) GOTO 6031 301 : IF (NFI(IA).EQ.NFI(ISUIV)) GOTO 6031 302 : DIS=(X(1,NTENT)-X(1,NFI(IA)))**2+(X(2,NTENT)-X(2,NFI(IA)))**2 303 : IF (DIS.GT.XL) GOTO 6031 304 : XL=DIS 305 : IP=IA 306 : ITP=ITT 307 : 6031 CONTINUE 308 : 6030 CONTINUE 309 : * write (6,*) ' ip le plus proche ',ip 310 : C ITP IP EST LE POINT LE PLUS PROCHE 311 : IF (XL.LT.0.3333D0*(XLA*XLB)**0.25D0*X(3,NFI(IP))) GOTO 6100 312 : IF (XL.EQ.0.D0) GOTO 6100 313 : IDIP=MAI(ITP+1-1)+1 314 : IFIP=MAI(ITP+1) 315 : IPA=IP-1 316 : IF (IPA.LT.IDIP) IPA=IFIP 317 : IPS=IP+1 318 : IF (IPS.GT.IFIP) IPS=IDIP 319 : SCALA=(X(1,NFI(IP))-X(1,NTENT))*(X(1,NFI(IPA))-X(1,NTENT)) 320 : # +(X(2,NFI(IP))-X(2,NTENT))*(X(2,NFI(IPA))-X(2,NTENT)) 321 : VECTA=(X(1,NFI(IP ))-X(1,NTENT))*(X(2,NFI(IPA))-X(2,NTENT)) 322 : # -(X(2,NFI(IP ))-X(2,NTENT))*(X(1,NFI(IPA))-X(1,NTENT)) 323 : IF (SCALA.EQ.0.D0.AND.VECTA.EQ.0.D0) GOTO 1500 324 : ANGA=ATAN2(-VECTA,SCALA) 325 : * write (6,*) 'atan2 - 2 ',anga 326 : 327 : IF (ANGA.GT.2.34D0.OR.ANGA.LT.0.D0) GOTO 6100 328 : SCALS=(X(1,NFI(IPS))-X(1,NTENT))*(X(1,NFI(IP ))-X(1,NTENT)) 329 : # +(X(2,NFI(IPS))-X(2,NTENT))*(X(2,NFI(IP ))-X(2,NTENT)) 330 : VECTS=(X(1,NFI(IPS))-X(1,NTENT))*(X(2,NFI(IP ))-X(2,NTENT)) 331 : # -(X(2,NFI(IPS))-X(2,NTENT))*(X(1,NFI(IP ))-X(1,NTENT)) 332 : IF (SCALS.EQ.0.D0.AND.VECTS.EQ.0.D0) GOTO 1500 333 : ANGS=ATAN2(-VECTS,SCALS) 334 : * write (6,*) 'atan2 - 3 ',angs 335 : IF (ANGS.GT.2.34D0.OR.ANGS.LE.1.d-5) GOTO 6100 336 : C OK ON Y VA 337 : NUMELG=NUMELG+1 338 : NUM(1,NUMELG)=NFI(IPRES) 339 : NUM(2,NUMELG)=NFI(I) 340 : NUM(3,NUMELG)=NFI(ISUIV) 341 : NUM(4,NUMELG)=NTENT 342 : NUMNP=NUMNP+1 343 : IF (X(3,NFI(I)).GT.XIN) GOTO 6060 344 : X(3,NTENT)=X(3,NFI(I))+MIN(COS*X(3,NFI(I)),XMOMI) 345 : GOTO 6062 346 : 6060 IF (X(3,NFI(I)).LT.XSU) GOTO 6061 347 : X(3,NTENT)=X(3,NFI(I))*(1.D0-COS) 348 : GOTO 6062 349 : 6061 X(3,NTENT)=XMOY 350 : 6062 CONTINUE 351 : NFI(I)=NUMNP 352 : VECTOR=-XGRAND 353 : IAUX=0 354 : I=IPRES 355 : GOTO 250 356 : 6100 CONTINUE 357 : IF(ANG.LT.-1.9D0) GOTO 250 358 : C SI 4 NOEUDS SUR LE CONTOUR TEST SUPPLEMENTAIRE 359 : IF (IFI-IDI.NE.3) GOTO 6101 360 : IQUAT=IPRES-1 361 : IF (IQUAT.LT.IDI) IQUAT=IFI 362 : SUC=(X(1,NFI(IPRES))-X(1,NFI(IQUAT)))*(X(2,NFI(ISUIV))- 363 : # X(2,NFI(IQUAT)))-(X(2,NFI(IPRES))-X(2,NFI(IQUAT)))* 364 : #(X(1,NFI(ISUIV))-X(1,NFI(IQUAT))) 365 : IF (SUC.GT.-VECTSA/5.D0) GOTO 6101 366 : IAA=I 367 : I=ISUIV 368 : ISUIV=IQUAT 369 : IQUAT=IPRES 370 : IPRES=IAA 371 : GOTO 6101 372 : 6101 IAX=IAUX 373 : IAUX=0 374 : VECTOR=-XGRAND 375 : C3=(X(1,NFI(ISUIV))-X(1,NFI(IPRES)))**2+(X(2,NFI(ISUIV))- 376 : # X(2,NFI(IPRES)))**2 377 : C3=SQRT(C3) 378 : IF(C3.GT.1.5d0*X(3,NFI(I)).AND.MAI(IT+1)-MAI(IT-1+1).NE.3) 379 : > GOTO 261 380 : NUMELG=NUMELG+1 381 : IF (NUMELG.GT.MAXELE) GOTO 1500 382 : NUM(1,NUMELG)=NFI(IPRES) 383 : NUM(2,NUMELG)=NFI(I) 384 : NUM(3,NUMELG)=NFI(ISUIV) 385 : IF (NBNN.EQ.4) NUM(4,NUMELG)=0 386 : IFI=IFI-1 387 : IF=IF-1 388 : X(3,NFI(IPRES))=SQRT(X(3,NFI(IPRES))*X(3,NFI(I))) 389 : X(3,NFI(ISUIV))=SQRT(X(3,NFI(ISUIV))*X(3,NFI(I))) 390 : DO 255 J=I,IF 391 : NFI(J)=NFI(J+1) 392 : 255 CONTINUE 393 : DO 256 J=IT,ITOUR 394 : 256 MAI(J+1)=MAI(J+1)-1 395 : I=IPRES-1 396 : IDI=MAI(IT-1+1)+1 397 : IFI=MAI(IT+1) 398 : IF (MAI(IT+1)-MAI(IT-1+1).NE.2) GOTO 250 399 : ID=MAI(IT-1+1)+1 400 : DO 257 J=ID,IF 401 : NFI(J)=NFI(J+2) 402 : 257 CONTINUE 403 : 244 CONTINUE 404 : ITOUR=ITOUR-1 405 : IF (ITOUR.EQ.0) GOTO 1600 406 : DO 258 I=IT,ITOUR 407 : INAT(I)=INAT(I+1) 408 : 258 MAI(I+1)=MAI(I+1+1)-2 409 : IF=MAI(ITOUR+1) 410 : IT=IT-1 411 : GOTO 251 412 : 261 NUMNP=NUMNP+1 413 : IF (NUMNP.GE.MAXPTS) GOTO 1500 414 : A=X(3,NFI(IPRES))+X(3,NFI(ISUIV)) 415 : B=X(3,NFI(ISUIV))/A 416 : C=X(3,NFI(IPRES))/A 417 : X(1,NUMNP)=C*X(1,NFI(ISUIV))+B*X(1,NFI(IPRES)) 418 : X(2,NUMNP)=C*X(2,NFI(ISUIV))+B*X(2,NFI(IPRES)) 419 : IF (X(3,NFI(I)).GT.XIN) GOTO 6070 420 : X(3,NUMNP)=X(3,NFI(I))+MIN(COS*X(3,NFI(I)),XMOMI) 421 : GOTO 6072 422 : 6070 IF (X(3,NFI(I)).LT.XSU) GOTO 6071 423 : X(3,NUMNP)=X(3,NFI(I))*(1.d0-COS) 424 : GOTO 6072 425 : 6071 X(3,NUMNP)=XMOY 426 : 6072 CONTINUE 427 : NUANC=NUANC+1 428 : NUMELG=NUMELG+1 429 : IF (NUMELG+1.GE.MAXELE) GOTO 1500 430 : NUM(1,NUMELG)=NFI(IPRES) 431 : NUM(2,NUMELG)=NFI(I) 432 : NUM(3,NUMELG)=NUMNP 433 : IF (NBNN.EQ.4) NUM(4,NUMELG)=0 434 : NUMELG=NUMELG+1 435 : NUM(1,NUMELG)=NFI(I) 436 : NUM(2,NUMELG)=NFI(ISUIV) 437 : NUM(3,NUMELG)=NUMNP 438 : IF (NBNN.EQ.4) NUM(4,NUMELG)=0 439 : NFI(I)=NUMNP 440 : I=IPRES-1 441 : GOTO 250 442 : 1200 CONTINUE 443 : IDI=MAI(IT-1+1)+1 444 : IFI=MAI(IT+1) 445 : DO 1 I=IDI,IFI 446 : IPRES=I-1 447 : IF (IPRES.LT.IDI) IPRES=IFI 448 : ISUIV=I+1 449 : IF (ISUIV.GT.IFI) ISUIV=IDI 450 : XL1=SQRT((X(1,NFI(IPRES))-X(1,NFI(I)))**2+ 451 : # (X(2,NFI(IPRES))-X(2,NFI(I)))**2) 452 : XL2=SQRT((X(1,NFI(ISUIV))-X(1,NFI(I)))**2+ 453 : # (X(2,NFI(ISUIV))-X(2,NFI(I)))**2) 454 : XDIR=XL2*(X(2,NFI(IPRES))-X(2,NFI(I)))+XL1*(X(2,NFI(I))-X(2, 455 : # NFI(ISUIV))) 456 : YDIR=XL1*(X(1,NFI(ISUIV))-X(1,NFI(I)))+XL2*(X(1,NFI(I))-X(1, 457 : # NFI(IPRES))) 458 : SDIR=SQRT(XDIR**2+YDIR**2) 459 : IF (SDIR.EQ.0.D0) XDIR=XL2*(X(2,NFI(IPRES))-X(2,NFI(I))) 460 : IF (SDIR.EQ.0.D0) YDIR=-(XL2*(X(1,NFI(IPRES))-X(1,NFI(I)))) 461 : IF (SDIR.EQ.0.D0) SDIR=SQRT(XDIR**2+YDIR**2) 462 : RAP=X(3,NFI(I))/SDIR 463 : IF (NBNN.EQ.3) RAP=RAP*0.85D0 464 : XAUX(1,I)=X(1,NFI(I))+RAP*XDIR 465 : XAUX(2,I)=X(2,NFI(I))+RAP*YDIR 466 : 1 CONTINUE 467 : C TESTS POUR LA SEPARATION OU LA FUSION DES CONTOURS 468 : ITO=0 469 : IP1=0 470 : IP2=0 471 : RAP=10 472 : IDI=MAI(IT-1+1)+1 473 : IFI=MAI(IT+1) 474 : IF (IFI-IDI.LE.3) GOTO 373 475 : IFIP=IFI-3 476 : IDIP=IDI 477 : DO 351 I=IDIP,IFIP 478 : IG=I+3 479 : IPRES=I-1 480 : IF (IPRES.LT.IDI) IPRES=IFI 481 : ISUIV=I+1 482 : IF (ISUIV.GT.IFI) ISUIV=IDI 483 : IFTA=IFI-MAX(0,IDIP+2-I) 484 : IF (IG.GT.IFTA) GOTO 351 485 : DO 352 J=IG,IFTA 486 : JSUIV=J+1 487 : JPRES=J-1 488 : IF (JPRES.LT.IDI) JPRES=IFI 489 : IF (JSUIV.GT.IFI) JSUIV=IDI 490 : IF (NFI(I).EQ.NFI(J)) GOTO 352 491 : IF (NFI(IPRES).EQ.NFI(J)) GOTO 352 492 : IF (NFI(I+1).EQ.NFI(J)) GOTO 352 493 : IF (NFI(I).EQ.NFI(JSUIV)) GOTO 352 494 : IF (NFI(I).EQ.NFI(JPRES)) GOTO 352 495 : DI=SQRT((X(1,NFI(I))-X(1,NFI(J)))**2+(X(2,NFI(I))-X(2,NFI(J)))**2) 496 : XRAP=DI/MAX(X(3,NFI(I)),X(3,NFI(J))) 497 : IF (XRAP.EQ.0.D0) GOTO 1500 498 : C REGLAGE POSSIBLE: A PARTIR DE QUAND COUPER EN 2 LE SEGMENT GENERE 499 : IF (XRAP.GT.1.415D0) GOTO 352 500 : IF (XRAP.GT.RAP*0.99999D0) GOTO 352 501 : C VERIFICATION QUE LES POINTS SE FONT FACE 502 : SCAL=(XAUX(1,I)-X(1,NFI(I)))*(XAUX(1,J)-X(1,NFI(J))) 503 : # +(XAUX(2,I)-X(2,NFI(I)))*(XAUX(2,J)-X(2,NFI(J))) 504 : C REGLAGE POUR DETECTER SI LES POINTS SE FONT FACE LA FUSION MARCHE 505 : C DE LA MEME FACON 506 : IF (SCAL.GE.0.D0) GOTO 352 507 : SCAL=(XAUX(1,I)-X(1,NFI(I)))*(X(1,NFI(J))-X(1,NFI(I))) 508 : # +(XAUX(2,I)-X(2,NFI(I)))*(X(2,NFI(J))-X(2,NFI(I))) 509 : IF (SCAL.LT.0.D0) GOTO 352 510 : SCAL=(XAUX(1,J)-X(1,NFI(J)))*(X(1,NFI(I))-X(1,NFI(J))) 511 : # +(XAUX(2,J)-X(2,NFI(J)))*(X(2,NFI(I))-X(2,NFI(J))) 512 : IF (SCAL.LT.0.D0) GOTO 352 513 : C ON INTERDIT DE COUPER UN COUTOUR EXTERIEUR EN UN EXTERIEUR + 514 : C UN INTERIEUR OU UN INTERIEUR EN DEUX INTERIEUR 515 : NATINI=INAT(IT) 516 : SURF1=0.D0 517 : SURF2=0.D0 518 : LPRES=IFI 519 : DO 1001 L=IDI,IFI 520 : SU=X(1,NFI(LPRES))*X(2,NFI(L))-X(2,NFI(LPRES))*X(1,NFI(L)) 521 : IF (L.LE.I.OR.LPRES.GE.J) SURF1=SURF1+SU 522 : IF (L.GT.I.AND.LPRES.LT.J) SURF2=SURF2+SU 523 : LPRES=L 524 : 1001 CONTINUE 525 : SU=X(1,NFI(I))*X(2,NFI(J))-X(2,NFI(I))*X(1,NFI(J)) 526 : SURF1=SURF1+SU 527 : SURF2=SURF2-SU 528 : IF (NATINI.GE.0.AND.SURF1*SURF2.LT.0.D0) GOTO 352 529 : IF (NATINI.LT.0.AND.SURF1.LT.0.D0.AND.SURF2.LT.0.D0) GOTO 352 530 : C TEST SI LE SEGMENT CANDIDAT EN COUPE UN AUTRE 531 : LL1=NFI(I) 532 : LL2=NFI(J) 532 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 533 : IF (IRECL.NE.0) GOTO 352 534 : LL3=NFI(IPRES) 535 : LL4=NFI(ISUIV) 535 : CALL TRANSB(X,IRECL,LL1,LL2,LL3,LL4) 536 : IF (IRECL.NE.0) GOTO 352 537 : LL1=NFI(J) 538 : LL2=NFI(I) 539 : LL3=NFI(JPRES) 540 : LL4=NFI(JSUIV) 540 : CALL TRANSB(X,IRECL,LL1,LL2,LL3,LL4) 541 : IF (IRECL.NE.0) GOTO 352 542 : 380 CONTINUE 543 : C OK 544 : RAP=XRAP 545 : ITO=IT 546 : IP1=I 547 : IP2=J 548 : SUK1=SURF1 549 : SUK2=SURF2 550 : 352 CONTINUE 551 : 351 CONTINUE 552 : IF (ITO.EQ.0) GOTO 353 553 : IF (IIMPI.EQ.1) WRITE (IOIMP,5566) NFI(IP1),NFI(IP2) 554 : 5566 FORMAT (' SEPARATION DES CONTOURS ',2I6) 555 : IDI=MAI(ITO-1+1)+1 556 : IFI=MAI(ITO+1) 557 : IDEB=MAI(ITOUR+1)+1 558 : ITOUR=ITOUR+1 559 : IF (ITOUR.GE.MAIMAX) GOTO 1500 560 : IF (IDEB.GE.NFMAX) GOTO 1500 561 : NFI(IDEB)=NFI(IP1) 562 : IDEB=IDEB+1 563 : IF (IDEB.GE.NFMAX) GOTO 1500 564 : IP=IP2 565 : 355 IF (IP.GT.IFI) IP=IDI 566 : IF (IP.EQ.IP1) GOTO 354 567 : NFI(IDEB)=NFI(IP) 568 : IDEB=IDEB+1 569 : IF (IDEB.GE.NFMAX) GOTO 1500 570 : IP=IP+1 571 : GOTO 355 572 : 354 IDEB=IDEB-1 573 : MAI(ITOUR+1)=IDEB 574 : INAT(ITOUR)=1 575 : IF (SUK1.LT.0.D0) INAT(ITOUR)=-1 576 : ITOUR=ITOUR+1 577 : IF (ITOUR.GE.MAIMAX) GOTO 1500 578 : IDEB=IDEB+1 579 : IF (IDEB.GE.NFMAX) GOTO 1500 580 : IP=IP1 581 : 357 IF (IP.GT.IP2) GOTO 356 582 : NFI(IDEB)=NFI(IP) 583 : IP=IP+1 584 : IDEB=IDEB+1 585 : IF (IDEB.GE.NFMAX) GOTO 1500 586 : GOTO 357 587 : 356 CONTINUE 588 : IDEB=IDEB-1 589 : MAI(ITOUR+1)=IDEB 590 : INAT(ITOUR)=1 591 : IF (SUK2.LT.0.D0) INAT(ITOUR)=-1 592 : C SUPPRESSION DE ITO 593 : IDEC=MAI(ITO+1)-MAI(ITO-1+1) 594 : ID=MAI(ITO-1+1)+1 595 : ITOUR=ITOUR-1 596 : DO 358 I=ITO,ITOUR 597 : INAT(I)=INAT(I+1) 598 : 358 MAI(I+1)=MAI(I+1+1)-IDEC 599 : IF=MAI(ITOUR+1) 600 : DO 359 I=ID,IF 601 : 359 NFI(I)=NFI(I+IDEC) 602 : IT=IT-1 603 : GOTO 251 604 : 353 CONTINUE 605 : C DEUXIEME TEST PORTANT SUR LES XAUX PARCOUR DECALE TRAITEMENT SIMILAIR 606 : ITO=0 607 : IP1=0 608 : IP2=0 609 : RAP=10 610 : IDI=MAI(IT-1+1)+1 611 : IFI=MAI(IT+1) 612 : IFIP=IFI-2 613 : IDIP=IDI 614 : DO 371 I=IDIP,IFIP 615 : IG=I+2 616 : IFTA=IFI-MAX(0,IDIP+1-I) 617 : IF (IG.GT.IFTA) GOTO 371 618 : DO 372 J=IG,IFTA 619 : IF (I.EQ.IDI.AND.J.EQ.IFI) GOTO 372 620 : DI=SQRT((XAUX(1,I)-XAUX(1,J))**2+(XAUX(2,I)-XAUX(2,J))**2) 621 : XRAP=DI/MAX(X(3,NFI(I)),X(3,NFI(J))) 622 : IF (XRAP.GT.0.707D0) GOTO 372 623 : C ON VERIFIE QUE LES POINTS SONT EN FACE 624 : SCAL=(XAUX(1,I)-X(1,NFI(I)))*(XAUX(1,J)-X(1,NFI(J))) 625 : # +(XAUX(2,I)-X(2,NFI(I)))*(XAUX(2,J)-X(2,NFI(J))) 626 : IF (SCAL.GE.0D0) GOTO 372 627 : XRAP=-SCAL/(X(3,NFI(I))*X(3,NFI(J))) 628 : IF (XRAP.GT.RAP*0.99999D0) GOTO 372 629 : NATINI=INAT(IT) 630 : SURF1=0.D0 631 : SURF2=0.D0 632 : LPRES=IFI 633 : DO 1002 L=IDI,IFI 634 : SU=X(1,NFI(LPRES))*X(2,NFI(L))-X(2,NFI(LPRES))*X(1,NFI(L)) 635 : IF (L.LE.I.OR.LPRES.GE.J) SURF1=SURF1+SU 636 : IF (L.GT.I.AND.LPRES.LT.J) SURF2=SURF2+SU 637 : LPRES=L 638 : 1002 CONTINUE 639 : SU=X(1,NFI(I))*X(2,NFI(J))-X(2,NFI(I))*X(1,NFI(J)) 640 : SURF1=SURF1+SU 641 : SURF2=SURF2-SU 642 : IF (NATINI.GE.0.AND.SURF1*SURF2.LT.0.D0) GOTO 372 643 : IF (NATINI.LT.0.AND.SURF1.LT.0.D0.AND.SURF2.LT.0.D0) GOTO 372 644 : C TEST SI LES DEUX SEGMENTS CANDITATS EN COUPENT UN AUTRE 645 : IF (NUMNP+1.GE.MAXPTS) GOTO 1500 646 : X(1,NUMNP+1)=(X(3,NFI(J))*XAUX(1,I)+X(3,NFI(I))*XAUX(1,J))/ 647 : # (X(3,NFI(I))+X(3,NFI(J))) 648 : X(2,NUMNP+1)=(X(3,NFI(J))*XAUX(2,I)+X(3,NFI(I))*XAUX(2,J))/ 649 : # (X(3,NFI(I))+X(3,NFI(J))) 650 : LL1=NFI(I) 651 : LL2=NUMNP+1 651 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 652 : IF (IRECL.NE.0) GOTO 372 653 : IPRES=I-1 654 : IF (IPRES.LT.IDI) IPRES=IFI 655 : ISUIV=I+1 656 : IF (ISUIV.GT.IFI) ISUIV=IDI 657 : LL3=NFI(IPRES) 658 : LL4=NFI(ISUIV) 658 : CALL TRANSB(X,IRECL,LL1,LL2,LL3,LL4) 659 : IF (IRECL.NE.0) GOTO 372 660 : LL1=NUMNP+1 661 : LL2=NFI(J) 661 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 662 : IF (IRECL.NE.0) GOTO 372 663 : JPRES=J-1 664 : IF (JPRES.LT.IDI) JPRES=IFI 665 : JSUIV=J+1 666 : IF (JSUIV.GT.IFI) JSUIV=IDI 667 : LL1=NFI(J) 668 : LL2=NUMNP+1 669 : LL3=NFI(JPRES) 670 : LL4=NFI(JSUIV) 670 : CALL TRANSB(X,IRECL,LL1,LL2,LL3,LL4) 671 : IF (IRECL.NE.0) GOTO 372 672 : RAP=XRAP 673 : ITO=IT 674 : IP1=I 675 : IP2=J 676 : SUK1=SURF1 677 : SUK2=SURF2 678 : 372 CONTINUE 679 : 371 CONTINUE 680 : IF (ITO.EQ.0) GOTO 373 681 : IF (IIMPI.EQ.1) WRITE (IOIMP,7788) NFI(IP1),NFI(IP2) 682 : 7788 FORMAT (' SEPARATION 2EME TYPE ',2I6) 683 : C ON RAJOUTE UN POINT 684 : NUMNP=NUMNP+1 685 : X(1,NUMNP)=(X(3,NFI(IP2))*XAUX(1,IP1)+X(3,NFI(IP1))*XAUX(1,IP2))/ 686 : # (X(3,NFI(IP1))+X(3,NFI(IP2))) 687 : X(2,NUMNP)=(X(3,NFI(IP2))*XAUX(2,IP1)+X(3,NFI(IP1))*XAUX(2,IP2))/ 688 : # (X(3,NFI(IP1))+X(3,NFI(IP2))) 689 : IDI=MAI(ITO-1+1)+1 690 : IFI=MAI(ITO+1) 691 : IDEB=MAI(ITOUR+1)+1 692 : ITOUR=ITOUR+1 693 : IF (ITOUR.GE.MAIMAX) GOTO 1500 694 : IF (IDEB.GE.NFMAX) GOTO 1500 695 : NFI(IDEB)=NFI(IP1) 696 : IDEB=IDEB+1 697 : IF (IDEB.GE.NFMAX) GOTO 1500 698 : NFI(IDEB)=NUMNP 699 : X(3,NFI(IDEB))=SQRT(X(3,NFI(IP1))*X(3,NFI(IP2))) 700 : IDEB=IDEB+1 701 : IP=IP2 702 : 375 IF (IP.GT.IFI) IP=IDI 703 : IF (IP.EQ.IP1) GOTO 374 704 : IF (IDEB.GE.NFMAX) GOTO 1500 705 : NFI(IDEB)=NFI(IP) 706 : IDEB=IDEB+1 707 : IP=IP+1 708 : GOTO 375 709 : 374 IDEB=IDEB-1 710 : MAI(ITOUR+1)=IDEB 711 : INAT(ITOUR)=1 712 : IF (SUK1.LT.0.D0) INAT(ITOUR)=-1 713 : ITOUR=ITOUR+1 714 : IDEB=IDEB+1 715 : IF (ITOUR.GE.MAIMAX.OR.IDEB.GE.NFMAX) GOTO 1500 716 : IP=IP1 717 : 377 IF (IP.GT.IP2) GOTO 376 718 : NFI(IDEB)=NFI(IP) 719 : IP=IP+1 720 : IDEB=IDEB+1 721 : IF (IDEB.GE.NFMAX) GOTO 1500 722 : GOTO 377 723 : 376 CONTINUE 724 : NFI(IDEB)=NUMNP 725 : X(3,NFI(IDEB))=SQRT(X(3,NFI(IP1))*X(3,NFI(IP2))) 726 : MAI(ITOUR+1)=IDEB 727 : INAT(ITOUR)=1 728 : IF (SUK2.LT.0.D0) INAT(ITOUR)=-1 729 : C SUPPRESSION DE ITO 730 : IDEC=MAI(ITO+1)-MAI(ITO-1+1) 731 : ID=MAI(ITO-1+1)+1 732 : ITOUR=ITOUR-1 733 : DO 378 I=ITO,ITOUR 734 : INAT(I)=INAT(I+1) 735 : 378 MAI(I+1)=MAI(I+1+1)-IDEC 736 : IF=MAI(ITOUR+1) 737 : DO 379 I=ID,IF 738 : 379 NFI(I)=NFI(I+IDEC) 739 : IT=IT-1 740 : GOTO 251 741 : 373 CONTINUE 742 : GOTO 251 743 : 253 CONTINUE 744 : C POUR LE CAS DES MAILLAGES AVEC TROUS RECOLLEMENT EVENTUEL DES CONTOUR 745 : IOPT=0 746 : ITO1=0 747 : ITO2=0 748 : IP1=0 749 : IP2=0 750 : RAP=10 751 : ITOUI=ITOUR-1 752 : IF (ITOUI.EQ.0) GOTO 520 753 : DO 500 IT1=1,ITOUI 754 : IDI1=MAI(IT1-1+1)+1 755 : IFI1=MAI(IT1+1) 756 : ITD=IT1+1 757 : DO 502 IT2=ITD,ITOUR 758 : IF (IT2.EQ.IT1) GOTO 502 759 : C IMPOSSIBLE DE FUSIONNER DEUX CONTOURS EXTERIEURS 760 : IF (INAT(IT1).EQ.1.AND.INAT(IT2).EQ.1) GOTO 502 761 : IDI2=MAI(IT2-1+1)+1 762 : IFI2=MAI(IT2+1) 763 : DO 503 I2=IDI2,IFI2 764 : IPRES2=I2-1 765 : IF (IPRES2.LT.IDI2) IPRES2=IFI2 766 : ISUIV2=I2+1 767 : IF (ISUIV2.GT.IFI2) ISUIV2=IDI2 768 : DO 501 I1=IDI1,IFI1 769 : IPRES1=I1-1 770 : IF (IPRES1.LT.IDI1) IPRES1=IFI1 771 : ISUIV1=I1+1 772 : IF (ISUIV1.GT.IFI1) ISUIV1=IDI1 773 : IF (NFI(I1).EQ.NFI(I2)) GOTO 501 774 : IF (NFI(IPRES1).EQ.NFI(I2)) GOTO 501 775 : IF (NFI(ISUIV1).EQ.NFI(I2)) GOTO 501 776 : IF (NFI(I1).EQ.NFI(IPRES2)) GOTO 501 777 : IF (NFI(I1).EQ.NFI(ISUIV2)) GOTO 501 778 : DI=SQRT((X(1,NFI(I1))-X(1,NFI(I2)))**2+ 779 : # (X(2,NFI(I1))-X(2,NFI(I2)))**2) 780 : XRAP=DI/MAX(X(3,NFI(I1)),X(3,NFI(I2))) 781 : IF (XRAP.EQ.0.D0) GOTO 1500 782 : IF (XRAP.GT.1.415D0) GOTO 530 783 : IF (XRAP.GT.RAP*0.99999D0) GOTO 530 784 : C VERIFICATION QUE LES POINTS SE FONT FACE 785 : SCAL=(XAUX(1,I1)-X(1,NFI(I1)))*(XAUX(1,I2)-X(1,NFI(I2))) 786 : # +(XAUX(2,I1)-X(2,NFI(I1)))*(XAUX(2,I2)-X(2,NFI(I2))) 787 : IF (SCAL.GE.0D0) GOTO 530 788 : SCAL=(XAUX(1,I1)-X(1,NFI(I1)))*(X(1,NFI(I2))-X(1,NFI(I1))) 789 : # +(XAUX(2,I1)-X(2,NFI(I1)))*(X(2,NFI(I2))-X(2,NFI(I1))) 790 : IF (SCAL.LT.0.D0) GOTO 530 791 : SCAL=(XAUX(1,I2)-X(1,NFI(I2)))*(X(1,NFI(I1))-X(1,NFI(I2))) 792 : # +(XAUX(2,I2)-X(2,NFI(I2)))*(X(2,NFI(I1))-X(2,NFI(I2))) 793 : IF (SCAL.LT.0.D0) GOTO 530 794 : LL1=NFI(I1) 795 : LL2=NFI(I2) 795 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 796 : IF (IRECL.NE.0) GOTO 501 797 : 381 CONTINUE 798 : RAP=XRAP 799 : ITO1=IT1 800 : ITO2=IT2 801 : IP1=I1 802 : IP2=I2 803 : IOPT=1 804 : GOTO 501 805 : 530 CONTINUE 806 : C DEUXIEME POSSIBILITE 807 : DI=SQRT((XAUX(1,I1)-XAUX(1,I2))**2+(XAUX(2,I1)-XAUX(2,I2))**2) 808 : XRAP=DI/MAX(X(3,NFI(I1)),X(3,NFI(I2))) 809 : IF (XRAP.GT.0.8D0) GOTO 501 810 : IF (XRAP.GT.RAP*0.99999D0) GOTO 501 811 : C ON VERIFIE QUE LES POINTS SONT EN FACE 812 : SCAL=(XAUX(1,I1)-X(1,NFI(I1)))*(XAUX(1,I2)-X(1,NFI(I2))) 813 : # +(XAUX(2,I1)-X(2,NFI(I1)))*(XAUX(2,I2)-X(2,NFI(I2))) 814 : IF (SCAL.GE.0.D0) GOTO 501 815 : LL1=NFI(I1) 816 : LL2=NFI(I2) 816 : CALL TRANSA(NFI,MAI,ITOUR,X,XMOY,IRECL,LL1,LL2,LL3,LL4) 817 : IF (IRECL.NE.0) GOTO 501 818 : LL3=NFI(IPRES1) 819 : LL4=NFI(ISUIV1) 819 : CALL TRANSB(X,IRECL,LL1,LL2,LL3,LL4) 820 : IF (IRECL.NE.0) GOTO 501 821 : LL1=NFI(I2) 822 : LL2=NFI(I1) 823 : LL3=NFI(IPRES2) 824 : LL4=NFI(ISUIV2) 824 : CALL TRANSB(X,IRECL,LL1,LL2,LL3,LL4) 825 : IF (IRECL.NE.0) GOTO 501 826 : RAP=XRAP 827 : ITO1=IT1 828 : ITO2=IT2 829 : IP1=I1 830 : IP2=I2 831 : IOPT=2 832 : 501 CONTINUE 833 : 503 CONTINUE 834 : 502 CONTINUE 835 : 500 CONTINUE 836 : IF (IOPT.EQ.0) GOTO 520 837 : IF (IOPT.EQ.1) GOTO 531 838 : NUMNP=NUMNP+1 839 : IF (NUMNP.GE.MAXPTS) GOTO 1500 840 : X(1,NUMNP)=(X(3,NFI(IP2))*X(1,NFI(IP1))+X(3,NFI(IP1))* 841 : # X(1,NFI(IP2)))/(X(3,NFI(IP1))+X(3,NFI(IP2))) 842 : X(2,NUMNP)=(X(3,NFI(IP2))*X(2,NFI(IP1))+X(3,NFI(IP1))* 843 : # X(2,NFI(IP2)))/(X(3,NFI(IP1))+X(3,NFI(IP2))) 844 : 531 CONTINUE 845 : IF (IIMPI.EQ.1) WRITE (IOIMP,6655) ITO1,ITO2,NFI(IP1),NFI(IP2) 846 : 6655 FORMAT(' FUSION DES CONTOURS ',4I6) 847 : IDI1=MAI(ITO1-1+1)+1 848 : IFI1=MAI(ITO1+1) 849 : IDI2=MAI(ITO2-1+1)+1 850 : IFI2=MAI(ITO2+1) 851 : IDEB=MAI(ITOUR+1)+1 852 : ITOUR=ITOUR+1 853 : IF (ITOUR.GE.MAIMAX.OR.IDEB.GE.NFMAX) GOTO 1500 854 : IF (IOPT.EQ.1) GOTO 532 855 : NFI(IDEB)=NUMNP 856 : X(3,NFI(IDEB))=SQRT(X(3,NFI(IP1))*X(3,NFI(IP2))) 857 : IDEB=IDEB+1 858 : IF (IDEB.GE.NFMAX) GOTO1500 859 : 532 NFI(IDEB)=NFI(IP1) 860 : IP=IP1 861 : 505 IP=IP+1 862 : IF (IP.GT.IFI1) IP=IDI1 863 : IF (IP.EQ.IP1) GOTO 506 864 : IDEB=IDEB+1 865 : IF (IDEB.GE.NFMAX) GOTO 1500 866 : NFI(IDEB)=NFI(IP) 867 : GOTO 505 868 : 506 IDEB=IDEB+1 869 : IF (IDEB.GE.NFMAX) GOTO 1500 870 : NFI(IDEB)=NFI(IP1) 871 : IDEB=IDEB+1 872 : IF (IDEB.GE.NFMAX) GOTO 1500 873 : IF (IOPT.EQ.1) GOTO 534 874 : NFI(IDEB)=NUMNP 875 : X(3,NFI(IDEB))=SQRT(X(3,NFI(IP1))*X(3,NFI(IP2))) 876 : IDEB=IDEB+1 877 : IF (IDEB.GE.NFMAX) GOTO 1500 878 : 534 CONTINUE 879 : NFI(IDEB)=NFI(IP2) 880 : IP=IP2 881 : 507 IP=IP+1 882 : IF (IP.GT.IFI2) IP=IDI2 883 : IDEB=IDEB+1 884 : IF (IDEB.GE.NFMAX) GOTO 1500 885 : NFI(IDEB)=NFI(IP) 886 : IF (IP.NE.IP2) GOTO 507 887 : MAI(ITOUR+1)=IDEB 888 : INAT(ITOUR)=1 889 : IF (INAT(ITO1).EQ.-1.AND.INAT(ITO2).EQ.-1) INAT(ITOUR)=-1 890 : C SUPPRESSION DE ITO1 ITO2 891 : IDEC=MAI(ITO1+1)-MAI(ITO1-1+1) 892 : ID=MAI(ITO1-1+1)+1 893 : ITOUR=ITOUR-1 894 : DO 508 I=ITO1,ITOUR 895 : INAT(I)=INAT(I+1) 896 : 508 MAI(I+1)=MAI(I+1+1)-IDEC 897 : IF=MAI(ITOUR+1) 898 : DO 509 I=ID,IF 899 : 509 NFI(I)=NFI(I+IDEC) 900 : ITO2=ITO2-1 901 : IDEC=MAI(ITO2+1)-MAI(ITO2-1+1) 902 : ID=MAI(ITO2-1+1)+1 903 : ITOUR=ITOUR-1 904 : DO 510 I=ITO2,ITOUR 905 : INAT(I)=INAT(I+1) 906 : 510 MAI(I+1)=MAI(I+1+1)-IDEC 907 : IF=MAI(ITOUR+1) 908 : DO 511 I=ID,IF 909 : 511 NFI(I)=NFI(I+IDEC) 910 : GOTO 360 911 : 520 CONTINUE 912 : C REPARTITION DES POINTS SUR LE CONTOUR DECALE 913 : C FABRICATION DES TRIANGLES JOIGNANT LES DEUX CONTOURS 914 : IAJOU=ITOUR 915 : IT=0 916 : 320 IT=IT+1 917 : IAJOU=IAJOU+1 918 : IF (IAJOU.GE.MAIMAX) GOTO 1500 919 : IF (IT.GT.ITOUR) GOTO 321 920 : INAT(IAJOU)=INAT(IT) 921 : IDI=MAI(IT-1+1)+1 922 : IFI=MAI(IT+1) 923 : IAJDI=MAI(IAJOU-1+1)+1 924 : IF (NBNN.EQ.3) GOTO 6300 925 : C POUR LES QUADRANGLES ON UTILISE UN AUTRE PROCEDE DE GENERATION 926 : NUMNP=NUMNP+1 927 : IF (NUMNP.GE.MAXPTS) GOTO 1500 928 : IPREM=NUMNP 929 : X(1,NUMNP)=XAUX(1,IDI) 930 : X(2,NUMNP)=XAUX(2,IDI) 931 : X(3,NUMNP)=X(3,NFI(IDI)) 932 : IF (IAJDI.GE.NFMAX) GOTO 1500 933 : NFI(IAJDI)=NUMNP 934 : DO 6301 I=IDI,IFI 935 : ISUIV=I+1 936 : IF (ISUIV.GT.IFI) ISUIV=IDI 937 : XCOMP=X(3,NUMNP)*X(3,NFI(ISUIV)) 938 : XLONG=(XAUX(1,ISUIV)-X(1,NUMNP))**2+(XAUX(2,ISUIV)-X(2,NUMNP))**2 939 : IF (XLONG.LT.0.5D0*XCOMP) GOTO 6302 940 : IF (XLONG.GT.2.D0*XCOMP) GOTO 6303 941 : C ON FAIT UN QUADRILATERE 942 : 6320 CONTINUE 943 : NUMELG=NUMELG+1 944 : IF (NUMNP+1.GE.MAXPTS.OR.NUMELG.GE.MAXELE) GOTO 1500 945 : NUM(1,NUMELG)=NFI(I) 946 : NUM(2,NUMELG)=NFI(ISUIV) 947 : NUM(3,NUMELG)=NUMNP+1 948 : NUM(4,NUMELG)=NUMNP 949 : IF (ISUIV.EQ.IDI) GOTO 6312 950 : NUMNP=NUMNP+1 951 : X(1,NUMNP)=XAUX(1,ISUIV) 952 : X(2,NUMNP)=XAUX(2,ISUIV) 953 : X(3,NUMNP)=X(3,NFI(ISUIV)) 954 : IAJDI=IAJDI+1 955 : IF (IAJDI.GE.NFMAX) GOTO 1500 956 : NFI(IAJDI)=NUMNP 957 : GOTO 6301 958 : 6312 NUM(3,NUMELG)=IPREM 959 : IF (IPREM.NE.NUMNP) GOTO 6301 960 : NUMELG=NUMELG-1 961 : C LE SEGMENT EST TROP PETIT ==> ON FAIT UN TRIANGLE 962 : 6302 X(1,NUMNP)=0.5D0*(X(1,NUMNP)+XAUX(1,ISUIV)) 963 : X(2,NUMNP)=0.5D0*(X(2,NUMNP)+XAUX(2,ISUIV)) 964 : X(3,NUMNP)=SQRT(X(3,NUMNP)*X(3,NFI(ISUIV))) 965 : XAUX(1,I)=X(1,NUMNP) 966 : XAUX(2,I)=X(2,NUMNP) 967 : NUMELG=NUMELG+1 968 : IF (NUMELG.GE.MAXELE) GOTO 1500 969 : NUM(1,NUMELG)=NFI(I) 970 : NUM(2,NUMELG)=NFI(ISUIV) 971 : NUM(3,NUMELG)=NUMNP 972 : NUM(4,NUMELG)=0 973 : IF (ISUIV.NE.IDI) GOTO 6301 974 : IF (NUMNP.EQ.IPREM) GOTO 6301 975 : NUMU=NUMELG+1 976 : 6321 NUMU=NUMU-1 977 : IF (NUM(3,NUMU).NE.NUMNP) GOTO 6322 978 : NUM(3,NUMU)=IPREM 979 : C CORRECTION PROBLEME YALA SYST SUR CRAY 06/18/86 980 : IF (NUM(4,NUMU).NE.IPREM+1) GOTO 6321 981 : C ALORS COUPER NUMU EN 2 TRIANGLES POUR NE PAS AVOIR 2 QUAD SE TOUCHAN 982 : C PAR TROIS POINTS 983 : NUMELG=NUMELG+1 984 : IF (NUMELG.GE.MAXELE) GOTO 1500 985 : NUM(1,NUMELG)=NUM(4,NUMU) 986 : NUM(2,NUMELG)=NUM(2,NUMU) 987 : NUM(3,NUMELG)=NUM(3,NUMU) 988 : NUM(4,NUMELG)=0 989 : NUM(3,NUMU)=NUM(4,NUMU) 990 : NUM(4,NUMU)=0 991 : 6322 CONTINUE 992 : X(1,IPREM)=X(1,NUMNP) 993 : X(2,IPREM)=X(2,NUMNP) 994 : NUMNP=NUMNP-1 995 : IAJDI=IAJDI-1 996 : GOTO 6301 997 : C LE SEGMENT EST TROP GRAND ==> ON FAIT UN TRIANGLE ET UN QUADRANGLE 998 : 6303 IALTE=1-IALTE 999 : IF (NUMNP+2.GE.MAXPTS.OR.NUMELG+2.GE.MAXELE) GOTO 1500 1000 : X(1,NUMNP+1)=0.5D0*(X(1,NUMNP)+XAUX(1,ISUIV)) 1001 : X(2,NUMNP+1)=0.5D0*(X(2,NUMNP)+XAUX(2,ISUIV)) 1002 : X(3,NUMNP+1)=SQRT(X(3,NUMNP)*X(3,NFI(ISUIV))) 1003 : IAJDI=IAJDI+1 1004 : IF (IAJDI+1.GE.NFMAX) GOTO 1500 1005 : NFI(IAJDI)=NUMNP+1 1006 : X(1,NUMNP+2)=XAUX(1,ISUIV) 1007 : X(2,NUMNP+2)=XAUX(2,ISUIV) 1008 : X(3,NUMNP+2)=X(3,NFI(ISUIV)) 1009 : IAJDI=IAJDI+1 1010 : NFI(IAJDI)=NUMNP+2 1011 : NUMELG=NUMELG+1 1012 : NUM(1,NUMELG)=NUMNP+1 1013 : NUM(2,NUMELG)=NUMNP 1014 : NUM(3,NUMELG)=NFI(I) 1015 : NUM(4,NUMELG)=0 1016 : IF (IALTE.EQ.1) NUM(4,NUMELG)=NFI(ISUIV) 1017 : NUMELG=NUMELG+1 1018 : IF (IALTE.NE.1) GOTO 6330 1019 : NUM(1,NUMELG)=NUMNP+1 1020 : NUM(2,NUMELG)=NFI(ISUIV) 1021 : NUM(3,NUMELG)=NUMNP+2 1022 : NUM(4,NUMELG)=0 1023 : GOTO 6331 1024 : 6330 NUM(1,NUMELG)=NFI(I) 1025 : NUM(2,NUMELG)=NFI(ISUIV) 1026 : NUM(3,NUMELG)=NUMNP+2 1027 : NUM(4,NUMELG)=NUMNP+1 1028 : 6331 CONTINUE 1029 : NUMNP=NUMNP+2 1030 : IF (ISUIV.NE.IDI) GOTO 6301 1031 : NUM(3,NUMELG)=IPREM 1032 : NUMNP=NUMNP-1 1033 : IAJDI=IAJDI-1 1034 : GOTO 6301 1035 : 6301 CONTINUE 1036 : GOTO 14 1037 : 6300 CONTINUE 1038 : ICA=IDI 1039 : ICN=ICA 1040 : NUANC=NUMNP 1041 : NUMNP=NUMNP+1 1042 : IF (NUMNP.GE.MAXPTS) GOTO 1500 1043 : X(1,NUMNP)=XAUX(1,ICA) 1044 : X(2,NUMNP)=XAUX(2,ICA) 1045 : NUSAUV=NUMNP 1046 : X(3,NUMNP)=X(3,NFI(ICA)) 1047 : IF (IAJDI.GE.NFMAX) GOTO 1500 1048 : NFI(IAJDI)=NUMNP 1049 : DISTR=0.D0 1050 : 2 ICAS=ICA+1 1051 : DISTA=1.D0 1052 : ISAUT=0 1053 : IF (ICAS.GT.IFI) ICAS=IDI 1054 : 21 ICNS=ICN+1 1055 : ISAUT=ISAUT+1 1056 : IF (ICNS.GT.IFI) ICNS=IDI 1057 : IF (ICN.GT.IFI) GOTO 10 1058 : DISCO=X(3,NFI(ICNS)) 1059 : DISTT=SQRT((XAUX(1,ICNS)-X(1,NUMNP))**2+(XAUX(2,ICNS)-X(2,NUMNP)) 1060 : # **2) 1061 : DISTA=EXP((REAL(ISAUT-1)*LOG(DISTA)+LOG(DISCO))/REAL(ISAUT)) 1062 : IF (DISTT.GE.DISTA) GOTO 20 1063 : NUMELG=NUMELG+1 1064 : IF (NUMELG.GT.MAXELE) GOTO 1500 1065 : NUM(1,NUMELG)=NUSAUV 1066 : NUM(2,NUMELG)=NFI(ICN) 1067 : NUM(3,NUMELG)=NFI(ICNS) 1068 : ICN=ICN+1 1069 : GOTO 21 1070 : 20 CONTINUE 1071 : DISTSU=DISTT-DISTA 1072 : DINCR=SQRT((XAUX(1,ICNS)-XAUX(1,ICN))**2+(XAUX(2,ICNS)-XAUX(2,ICN) 1073 : # )**2) 1074 : RAP=MIN(1.D0,DISTSU/DINCR) 1075 : NUPREC=NUSAUV 1076 : NUMNP=NUMNP+1 1077 : IF (NUMNP.GE.MAXPTS) GOTO 1500 1078 : NUSAUV=NUMNP 1079 : IAJDI=IAJDI+1 1080 : IF (IAJDI.GE.NFMAX) GOTO 1500 1081 : X(3,NUMNP)=DISTA 1082 : NFI(IAJDI)=NUMNP 1083 : X(1,NUMNP)=XAUX(1,ICN)*RAP+(1.D0-RAP)*XAUX(1,ICNS) 1084 : X(2,NUMNP)=XAUX(2,ICN)*RAP+(1.D0-RAP)*XAUX(2,ICNS) 1085 : C VERIFICATION QU'ON NE PASSE PAS TROP PRES DES PTS 1086 : C NFI(ICN) ET NFI(ICNS) 1087 : ZRA=((X(1,NUMNP)-X(1,NFI(ICN)))*(X(2,NUMNP-1)-X(2,NFI(ICN)))-(X(2 1088 : # ,NUMNP)-X(2,NFI(ICN)))*(X(1,NUMNP-1)-X(1,NFI(ICN))))/(DISTA**2) 1089 : IF (ZRA.GT.0.4D0) GOTO 2003 1090 : XQ=X(1,NUMNP) 1091 : YQ=X(2,NUMNP) 1092 : X(1,NUMNP)=XQ+(1-ZRA)*(X(2,NUMNP-1)-YQ) 1093 : X(2,NUMNP)=YQ+(1-ZRA)*(XQ-X(1,NUMNP-1)) 1094 : 2003 CONTINUE 1095 : IF (ISAUT.EQ.1) GOTO 240 1096 : ICNP=ICN-1 1097 : IF (ICNP.LT.IDI) GOTO 240 1098 : C1=(X(1,NUPREC)-X(1,NFI(ICN)))**2+ 1099 : # (X(2,NUPREC)-X(2,NFI(ICN)))**2 1100 : C2=(X(1,NUSAUV)-X(1,NFI(ICNP)))**2+ 1101 : # (X(2,NUSAUV)-X(2,NFI(ICNP)))**2 1102 : IF (C1.LE.C2) GOTO 240 1103 : NUM(1,NUMELG)=NUSAUV 1104 : NUM(2,NUMELG)=NUPREC 1105 : NUM(3,NUMELG)=NFI(ICNP) 1106 : NUMELG=NUMELG+1 1107 : IF (NUMELG.GE.MAXELE) GOTO 1500 1108 : NUM(1,NUMELG)=NUSAUV 1109 : NUM(2,NUMELG)=NFI(ICNP) 1110 : NUM(3,NUMELG)=NFI(ICN) 1111 : ICA=ICN 1112 : GOTO 2 1113 : 240 CONTINUE 1114 : NUMELG=NUMELG+1 1115 : IF (NUMELG.GE.MAXELE) GOTO 1500 1116 : NUM(1,NUMELG)=NUSAUV 1117 : NUM(2,NUMELG)=NUPREC 1118 : NUM(3,NUMELG)=NFI(ICN) 1119 : ICA=ICN 1120 : GOTO 2 1121 : 10 CONTINUE 1122 : * CORRECTION PROBLEME VAGHI 5/2/87 1 SEUL POINT INTERNE 1123 : IF (NUSAUV.EQ.NUANC+1) GOTO 14 1124 : C=SQRT((X(1,NUSAUV)-X(1,NUANC+1))**2+(X(2,NUSAUV)-X(2,NUANC+1)) 1125 : # **2) 1126 : IF (C.GT.X(3,NFI(IFI))) GOTO 11 1127 : X(1,NUANC+1)=0.5D0*(X(1,NUANC+1)+X(1,NUSAUV)) 1128 : X(2,NUANC+1)=0.5D0*(X(2,NUANC+1)+X(2,NUSAUV)) 1129 : NUMU=NUMELG+1 1130 : 12 NUMU=NUMU-1 1131 : IF (NUM(1,NUMU).NE.NUSAUV) GOTO 15 1132 : NUM(1,NUMU)=NUANC+1 1133 : GOTO 12 1134 : 15 CONTINUE 1135 : NUMNP=NUMNP-1 1136 : IAJDI=IAJDI-1 1137 : GOTO 14 1138 : 11 NUMELG=NUMELG+1 1139 : IF (NUMELG.GE.MAXELE) GOTO 1500 1140 : NUM(1,NUMELG)=NUANC+1 1141 : NUM(2,NUMELG)=NUMNP 1142 : NUM(3,NUMELG)=NFI(IDI) 1143 : 14 MAI(IAJOU+1)=IAJDI 1144 : IF (MAI(IAJOU+1)-MAI(IAJOU-1+1).LT.3) IAJOU=IAJOU-1 1145 : IF (NUMELG.GE.MAXELE.OR.NUMNP.GE.MAXPTS) GOTO 1500 1146 : GOTO 320 1147 : 321 CONTINUE 1148 : IAJOU=IAJOU-1 1149 : IDEC=MAI(ITOUR+1)-IN+1 1150 : IF=MAI(IAJOU+1)-IDEC 1151 : DO 204 I=IN,IF 1152 : NFI(I)=NFI(I+IDEC) 1153 : 204 CONTINUE 1154 : ITDEC=IAJOU-ITOUR 1155 : DO 202 IT=1,ITDEC 1156 : INAT(IT)=INAT(IT+ITOUR) 1157 : MAI(IT+1)=MAI(IT+ITOUR+1)-IDEC 1158 : 202 CONTINUE 1159 : ITOUR=ITDEC 1160 : IF (ITOUR.EQ.0) GOTO 1600 1161 : IF (NUMNP.GE.MAXPTS.OR.NUMELG.GE.MAXELE) GOTO 1500 1162 : GOTO 200 1163 : 1500 IRECHA=27 1164 : ICLE=0 1165 : RETURN 1166 : 1600 CONTINUE 1167 : LEMAX=NUMELG 1168 : ICLE=2 1169 : NINF=NUMNP 1170 : RETURN 1171 : END 1172 : 1173 : 1174 : 1175 : 1176 :
© Cast3M 2003 - All rights reserved.
Disclaimer