1 : C PYRA3 SOURCE PV 07/10/05 21:16:25 5890 2 : C---------------------------------------------------------------------| 3 : C | 4 : SUBROUTINE PYRA3(II,JJ,IF1,IF2,IGAGNE) 5 : C | 6 : C CETTE SUBROUTINE TENTE DE CREER DEUX PYRAMIDES A PARTIR | 7 : C DES QUADRANGLES IF1 et IF2. | 8 : C - IGAGNE=1 EN CAS DE SUCCES | 9 : C - IGAGNE=0 EN CAS D'ECHEC | 10 : C | 11 : C---------------------------------------------------------------------| 12 : C 13 : IMPLICIT INTEGER(I-N) 14 : IMPLICIT REAL*8(A-H,O-Z) 15 : -INC TDEMAIT 16 : -INC CCOPTIO 17 : LOGICAL REPONS,FACET,SOLPYR,SOLTET,DIAGO,IN2,VERDIV,IN 18 : C 19 : * WRITE(6,1000) 20 : 1000 FORMAT(' ----->>> PYRA3 <<<-----') 21 : C 22 : nptini=nptmax 23 : nfcini=nfcmax 24 : ICTF=0 25 : ICTV=0 26 : IP1=IPRED(IF1,II) 27 : IP2=IPRED(IF1,IP1) 28 : JP1=ISUCC(IF2,II) 29 : JP2=ISUCC(IF2,JP1) 30 : * CREATION POINT MILIEU 31 : NPTMAX=NPTMAX+1 32 : XYZ(4,NPTMAX)=(XYZ(4,IP1)+XYZ(4,IP2)+XYZ(4,JP1)+XYZ(4,JP2))/4. 33 : * deplacement du point pour l'eloigner de ii jj 34 : xn1=(xyz(2,jj)-xyz(2,ii))*((xyz(3,ip1)+xyz(3,ip2))/2-xyz(3,ii))- 35 : > (xyz(3,jj)-xyz(3,ii))*((xyz(2,ip1)+xyz(2,ip2))/2-xyz(2,ii)) 36 : yn1=(xyz(3,jj)-xyz(3,ii))*((xyz(1,ip1)+xyz(1,ip2))/2-xyz(1,ii))- 37 : > (xyz(1,jj)-xyz(1,ii))*((xyz(3,ip1)+xyz(3,ip2))/2-xyz(3,ii)) 38 : zn1=(xyz(1,jj)-xyz(1,ii))*((xyz(2,ip1)+xyz(2,ip2))/2-xyz(2,ii))- 39 : > (xyz(2,jj)-xyz(2,ii))*((xyz(1,ip1)+xyz(1,ip2))/2-xyz(1,ii)) 40 : sn1=sqrt(xn1**2+yn1**2+zn1**2) 41 : xn1=xn1/sn1 42 : yn1=yn1/sn1 43 : zn1=zn1/sn1 44 : xn2=((xyz(2,jp1)+xyz(2,jp2))/2-xyz(2,ii))*(xyz(3,jj)-xyz(3,ii))- 45 : > ((xyz(3,jp1)+xyz(3,jp2))/2-xyz(3,ii))*(xyz(2,jj)-xyz(2,ii)) 46 : yn2=((xyz(3,jp1)+xyz(3,jp2))/2-xyz(3,ii))*(xyz(1,jj)-xyz(1,ii))- 47 : > ((xyz(1,jp1)+xyz(1,jp2))/2-xyz(1,ii))*(xyz(3,jj)-xyz(3,ii)) 48 : zn2=((xyz(1,jp1)+xyz(1,jp2))/2-xyz(1,ii))*(xyz(2,jj)-xyz(2,ii))- 49 : > ((xyz(2,jp1)+xyz(2,jp2))/2-xyz(2,ii))*(xyz(1,jj)-xyz(1,ii)) 50 : sn2=sqrt(xn2**2+yn2**2+zn2**2) 51 : xn2=xn2/sn2 52 : yn2=yn2/sn2 53 : zn2=zn2/sn2 54 : xn=(xn1+xn2)/2 55 : yn=(yn1+yn2)/2 56 : zn=(zn1+zn2)/2 57 : sn=sqrt(xn**2+yn**2+zn**2) 58 : xn=xn/sn 59 : yn=yn/sn 60 : zn=zn/sn 61 : * xmil=(xyz(1,ii)+xyz(1,jj))/2 62 : * ymil=(xyz(2,ii)+xyz(2,jj))/2 63 : * zmil=(xyz(3,ii)+xyz(3,jj))/2 64 : * 65 : xv=xyz(1,jj)-xyz(1,ii) 66 : yv=xyz(2,jj)-xyz(2,ii) 67 : zv=xyz(3,jj)-xyz(3,ii) 68 : sv=sqrt(xv**2+yv**2+zv**2) 69 : xv=xv/sv 70 : yv=yv/sv 71 : zv=zv/sv 72 : xli1=xv*(xyz(1,ip1)-xyz(1,ii))+yv*(xyz(2,ip1)-xyz(2,ii))+ 73 : > zv*(xyz(3,ip1)-xyz(3,ii)) 74 : xli2=xv*(xyz(1,ip2)-xyz(1,ii))+yv*(xyz(2,ip2)-xyz(2,ii))+ 75 : > zv*(xyz(3,ip2)-xyz(3,ii)) 76 : xlj1=xv*(xyz(1,jp1)-xyz(1,ii))+yv*(xyz(2,jp1)-xyz(2,ii))+ 77 : > zv*(xyz(3,jp1)-xyz(3,ii)) 78 : xlj2=xv*(xyz(1,jp2)-xyz(1,ii))+yv*(xyz(2,jp2)-xyz(2,ii))+ 79 : > zv*(xyz(3,jp2)-xyz(3,ii)) 80 : xl=(xli1+xli2+xlj1+xlj2+2*sv+2*0)/8 81 : xl=0.5*sv 82 : xmil=xyz(1,ii)+xl*xv 83 : ymil=xyz(2,ii)+xl*yv 84 : zmil=xyz(3,ii)+xl*zv 85 : expf = xyz(4,nptmax) 86 : xyz(1,nptmax)=xmil+xn*expf*expfac 87 : xyz(2,nptmax)=ymil+yn*expf*expfac 88 : xyz(3,nptmax)=zmil+zn*expf*expfac 89 : * write (6,*) ' pyra3 creation de 2 elements et pt ',nptmax 90 : * write (6,*) (xyz(i,ii),i=1,4) 91 : * write (6,*) (xyz(i,jj),i=1,4) 92 : * write (6,*) (xyz(i,ip1),i=1,4) 93 : * write (6,*) (xyz(i,ip2),i=1,4) 94 : * write (6,*) (xyz(i,jp1),i=1,4) 95 : * write (6,*) (xyz(i,jp2),i=1,4) 96 : * write (6,*) (xyz(i,nptmax),i=1,4) 97 : * CREATION DES ELEMENTS 98 : IPTT=NPTMAX 98 : CALL DIST(iptt,nptaux,GL,IOK,ii,jj,ip1,ip2,jp1,jp2,0,0,0,0) 99 : IF (IOK.EQ.0) THEN 100 : NPTMAX=nptini 101 : WRITE (6,*) ' pyra3 DIST ',nptaux 102 : return 103 : ENDIF 104 : IF (gl.lt.xyz(4,iptt)/4) then 105 : write (6,*) 'pyra3 GL-1' 106 : nptmax=nptini 107 : return 108 : endif 109 : IF (.NOT.IN2(ii,IPTT,nfcini).or..NOT.IN2(jj,IPTT,nfcini)) THEN 110 : write (6,*) ' in2 echec ' 111 : NPTMAX=nptini 112 : return 113 : ENDIF 114 : * creations des faces 115 : NFCMAX=NFCMAX+1 116 : IF3=NFCMAX 117 : NFC(1,IF3)=II 118 : NFC(2,IF3)=iptt 119 : NFC(3,IF3)=IP1 120 : NFC(4,IF3)=0 121 : * 122 : NFCMAX=NFCMAX+1 123 : IF4=NFCMAX 124 : NFC(1,IF4)=IP1 125 : NFC(2,IF4)=IPTT 126 : NFC(3,IF4)=IP2 127 : NFC(4,IF4)=0 128 : * 129 : NFCMAX=NFCMAX+1 130 : IF5=NFCMAX 131 : NFC(1,IF5)=IP2 132 : NFC(2,IF5)=IPTT 133 : NFC(3,IF5)=JJ 134 : NFC(4,IF5)=0 135 : * la face commune 136 : NFCMAX=NFCMAX+1 137 : IF6=NFCMAX 138 : NFC(1,IF6)=JJ 139 : NFC(2,IF6)=IPTT 140 : NFC(3,IF6)=II 141 : NFC(4,IF6)=0 142 : * 143 : NFCMAX=NFCMAX+1 144 : IF7=NFCMAX 145 : NFC(1,IF7)=JJ 146 : NFC(2,IF7)=iptt 147 : NFC(3,IF7)=jp2 148 : NFC(4,IF7)=0 149 : * 150 : NFCMAX=NFCMAX+1 151 : IF8=NFCMAX 152 : NFC(1,IF8)=jp2 153 : NFC(2,IF8)=iptt 154 : NFC(3,IF8)=jp1 155 : NFC(4,IF8)=0 156 : * 157 : NFCMAX=NFCMAX+1 158 : IF9=NFCMAX 159 : NFC(1,IF9)=jp1 160 : NFC(2,IF9)=iptt 161 : NFC(3,IF9)=II 162 : NFC(4,IF9)=0 163 : * si necessaire verification diago 164 : IF (DIAGO(iptt,ii,diacri)) GOTO 275 165 : IF (DIAGO(iptt,jj,diacri)) GOTO 275 166 : IF (DIAGO(iptt,ip1,diacri)) GOTO 275 167 : IF (DIAGO(iptt,ip2,diacri)) GOTO 275 168 : IF (DIAGO(iptt,jp1,diacri)) GOTO 275 169 : IF (DIAGO(iptt,jp2,diacri)) GOTO 275 170 : goto 276 171 : 275 continue 172 : nptmax=nptini 173 : nfcmax=nfcini 174 : write (6,*) ' pyra3 echec diago' 175 : return 176 : 276 continue 177 : C 177 : CALL REPSUB(IF1) 177 : CALL REPSUB(IF3) 177 : CALL REPSUB(IF4) 177 : CALL REPSUB(IF5) 177 : CALL REPSUB(IF6) 178 : * verification du premier element element 179 : IF (.NOT.SOLPYR(IF1,IF3,IF4,IF5,IF6)) then 180 : write (6,*) 'pyra3-1 solpyr invalide' 181 : GOTO 160 182 : endif 183 : IF (.NOT.FACET(IF3)) then 184 : write(6,*) ' pyra3 facet if3 invalide' 185 : GOTO 160 186 : ENDIF 187 : IF (.NOT.FACET(IF4)) then 188 : write(6,*) ' pyra3 facet if4 invalide' 189 : GOTO 160 190 : ENDIF 191 : IF (.NOT.FACET(IF5)) then 192 : write(6,*) ' pyra3 facet if5 invalide' 193 : GOTO 160 194 : ENDIF 194 : CALL REPSUB(IF2) 194 : CALL REPSUB(IF6) 194 : CALL REPSUB(IF7) 194 : CALL REPSUB(IF8) 194 : CALL REPSUB(IF9) 195 : * verification du deuxieme element element 196 : IF (.NOT.SOLPYR(IF2,IF6,IF7,IF8,IF9)) then 197 : write (6,*) 'pyra3-2 solpyr invalide' 198 : GOTO 165 199 : endif 200 : IF (.NOT.FACET(IF7)) then 201 : write(6,*) ' pyra3 facet if7 invalide' 202 : GOTO 165 203 : ENDIF 204 : IF (.NOT.FACET(IF8)) then 205 : write(6,*) ' pyra1 facet if8 invalide' 206 : GOTO 165 207 : ENDIF 208 : IF (.NOT.FACET(IF9)) then 209 : write(6,*) ' pyra1 facet if9 invalide' 210 : GOTO 165 211 : ENDIF 212 : * OK pour creation elements 213 : NVOL=NVOL+1 214 : IF (NFV(1,IF1).EQ.0) NFV(1,IF1)=NVOL 215 : IF (NFV(1,IF1).NE.NVOL) NFV(2,IF1)=NVOL 216 : IF (NFV(1,IF3).EQ.0) NFV(1,IF3)=NVOL 217 : IF (NFV(1,IF3).NE.NVOL) NFV(2,IF3)=NVOL 218 : IF (NFV(1,IF4).EQ.0) NFV(1,IF4)=NVOL 219 : IF (NFV(1,IF4).NE.NVOL) NFV(2,IF4)=NVOL 220 : IF (NFV(1,IF5).EQ.0) NFV(1,IF5)=NVOL 221 : IF (NFV(1,IF5).NE.NVOL) NFV(2,IF5)=NVOL 222 : IF (NFV(1,IF6).EQ.0) NFV(1,IF6)=NVOL 223 : IF (NFV(1,IF6).NE.NVOL) NFV(2,IF6)=NVOL 224 : IVOL(9,NVOL)=35 225 : DO 280 I=1,4 226 : IVOL(I,NVOL)=NFC(I,IF1) 227 : 280 CONTINUE 228 : IVOL(5,NVOL)=iptt 229 : if (iimpi.eq.1) write (6,1102) nfacet,(ivol(i,nvol),i=1,5) 230 : 1102 FORMAT(' PYRA3-1 facettes ',i5,' pyr5 ',8i5) 231 : NVOL=NVOL+1 232 : IVOL(9,NVOL)=35 233 : IF (NFV(1,IF2).EQ.0) NFV(1,IF2)=NVOL 234 : IF (NFV(1,IF2).NE.NVOL) NFV(2,IF2)=NVOL 235 : IF (NFV(1,IF6).EQ.0) NFV(1,IF6)=NVOL 236 : IF (NFV(1,IF6).NE.NVOL) NFV(2,IF6)=NVOL 237 : IF (NFV(1,IF7).EQ.0) NFV(1,IF7)=NVOL 238 : IF (NFV(1,IF7).NE.NVOL) NFV(2,IF7)=NVOL 239 : IF (NFV(1,IF8).EQ.0) NFV(1,IF8)=NVOL 240 : IF (NFV(1,IF8).NE.NVOL) NFV(2,IF8)=NVOL 241 : IF (NFV(1,IF9).EQ.0) NFV(1,IF9)=NVOL 242 : IF (NFV(1,IF9).NE.NVOL) NFV(2,IF9)=NVOL 243 : DO 282 I=1,4 244 : IVOL(I,NVOL)=NFC(I,IF2) 245 : 282 CONTINUE 246 : IVOL(5,NVOL)=iptt 247 : if (iimpi.eq.1) write (6,1101) nfacet,(ivol(i,nvol),i=1,5) 248 : 1101 FORMAT(' PYRA3-2 facettes ',i5,' pyr5 ',8i5) 249 : IGAGNE=1 250 : RETURN 251 : 165 continue 251 : CALL REPSUB(IF9) 251 : CALL REPSUB(IF8) 251 : CALL REPSUB(IF7) 251 : CALL REPSUB(IF6) 251 : CALL REPSUB(IF2) 252 : 160 continue 252 : CALL REPSUB(IF6) 252 : CALL REPSUB(IF5) 252 : CALL REPSUB(IF4) 252 : CALL REPSUB(IF3) 252 : CALL REPSUB(IF1) 253 : nptmax=nptini 254 : nfcmax=nfcini 255 : goto 9000 256 : 9000 CONTINUE 257 : END 258 : 259 :
© Cast3M 2003 - All rights reserved.
Disclaimer