1 : C ARC SOURCE CHAT 07/10/29 21:15:23 5940 2 : C 1234567890124567890124567898012456789012345678901234567890123456 3 : 4 : C SOURCE L. DI VALENTIN 20/06/97 5 : 6 : SUBROUTINE ARC(POINT,CENTRE,NORMAL,ANGLE,INBR, 7 : & DEN1OX,DEN2OX,POINT2) 8 : 9 : 10 : C Cette routine, appele par prtran.eso, est chargee de construire 11 : C un arc de cercle en accord avec les elements utilises, de centre 12 : C CENTRE, de normale NORMAL, d'extremites POINT et POINT2 et 13 : C d'angle ANGLE. 14 : C Il y a une redondance entre ces 3 derniers parametres. En fait, 15 : C l'ANGLE est prioritaire devant POINT2 a moins que POINT2 = POINT 16 : C (cf plus loin). La procedure suppose que POINT2 est la 2eme 17 : C extremite de l'arc de cercle. Il est donc important que POINT2 18 : C reference ce point geometrique. 19 : 20 : C INBR est le nombre de decoupage de l'arc de cercle. 21 : C Si INBR = 0, les densites DEN1OX et DEN2OX sont prises en 22 : C compte pour la decoupe. 23 : C 24 : C Si POINT2 = POINT alors quelque soit l'angle, la subroutine 25 : C construira un cercle. Dans ce cas, ANGLE sera change egal a 26 : C 2*pi (=(2.D0*XPI)) 27 : C 28 : C Le maillage resultant est sauve dans la pile a la fin. 29 : 30 : 31 : IMPLICIT INTEGER(I-N) 32 : IMPLICIT REAL*8 (A-H,O-Z) 33 : -INC CCREEL 34 : -INC CCOPTIO 35 : -INC CCGEOME 36 : -INC SMELEME 37 : -INC SMCOORD 38 : 39 : 40 : SEGMENT TABPAR(NBELEM) 41 : 42 : INTEGER POINT,CENTRE,NORMAL,INBR,POINT2 43 : REAL*8 ALPHA,ANGLE,DEN1,DEN2 44 : 45 : 46 : DEN1 = DEN1OX 47 : DEN2 = DEN2OX 48 : 49 : C Vecteur Normal : 50 : 51 : IF (IDIM.EQ.3) THEN 52 : Xn = XCOOR((NORMAL-1)*(IDIM+1) + 1) 53 : Yn = XCOOR((NORMAL-1)*(IDIM+1) + 2) 54 : Zn = XCOOR((NORMAL-1)*(IDIM+1) + 3) 55 : 56 : ALPHA = SQRT(Xn*Xn+Yn*Yn+Zn*Zn) 57 : 58 : Xn = Xn / ALPHA 59 : Yn = Yn / ALPHA 60 : Zn = Zn / ALPHA 61 : 62 : ELSE 63 : Xn = 0.D0 64 : Yn = 0.D0 65 : Zn = 1.D0 66 : 67 : ENDIF 68 : 69 : 70 : C On recupere les coordonnees du CENTRE 71 : 72 : XCEN = XCOOR((CENTRE-1)*(IDIM+1) + 1) 73 : YCEN = XCOOR((CENTRE-1)*(IDIM+1) + 2) 74 : IF (IDIM.EQ.3) THEN 75 : ZCEN = XCOOR((CENTRE-1)*(IDIM+1) + 3) 76 : ENDIF 77 : 78 : C On recupere les coordonnees du vecteur (CENTRE -> POINT) 79 : 80 : X1 = XCOOR((POINT-1)*(IDIM+1) + 1) - XCEN 81 : Y1 = XCOOR((POINT-1)*(IDIM+1) + 2) - YCEN 82 : 83 : IF (IDIM.EQ.3) THEN 84 : Z1 = XCOOR((POINT-1)*(IDIM+1) + 3) - ZCEN 85 : ELSE 86 : Z1 = 0.D0 87 : ENDIF 88 : 89 : RAYON = SQRT(X1*X1+Y1*Y1+Z1*Z1) 90 : 91 : 92 : C On calcule le vecteur definissant avec le vecteur precedent 93 : C le 2eme vecteur de la base du cercle 94 : 95 : X2 = Z1*Yn - Y1*Zn 96 : Y2 = X1*Zn - Z1*Xn 97 : Z2 = Y1*Xn - X1*Yn 98 : 99 : 100 : C Cas du cercle : 101 : 102 : IF (POINT.EQ.POINT2) THEN 103 : ANGLE = (2.D0*XPI) 104 : ENDIF 105 : 106 : 107 : NBPTS = XCOOR(/1)/(IDIM + 1) 108 : 109 : 110 : DLONG = ABS(ANGLE * RAYON) 111 : C DEN1 = 1.D0/INBR 112 : DEN1 = DEN1/DLONG 113 : DEN2 = DEN2/DLONG 114 : 114 : CALL DECOUP(INBR,DEN1,DEN2,APROG,NBELEM,DENI,DECA,DLONG) 115 : DLONG = ANGLE * RAYON 116 : 117 : C NBELEM = INBR 118 : NX = NBELEM - 1 119 : 120 : NBSOUS=0 121 : NBREF=0 121 : IF (ILCOUR.EQ.0) CALL ERREUR(16) 121 : IF (KDEGRE(ILCOUR).EQ.0) CALL ERREUR(16) 122 : NBNN=NBNNE(KDEGRE(ILCOUR)) 122 : IF (NBNN.EQ.0) CALL ERREUR(16) 123 : 124 : 125 : SEGINI MELEME 126 : SEGINI TABPAR 127 : 128 : ITYPEL = KDEGRE(ILCOUR) 129 : 130 : DIN = DEN1 131 : 132 : IPOO = XCOOR(/1)/(IDIM + 1 ) 133 : NUM(1,1) = POINT 134 : 135 : 136 : IF (NX.EQ.0) GOTO 110 137 : 138 : DO 102 I=1,NX 139 : DIN = DIN * APROG 140 : IF (ITYPEL.EQ.2) GOTO 103 141 : IPOO=IPOO+1 142 : NUM(2,I)=IPOO 143 : NUM(3,I)=IPOO+1 144 : 103 IPOO=IPOO+1 145 : TABPAR(I) = DIN 146 : NUM(1,I+1)=IPOO 147 : IF (ITYPEL.EQ.3) GOTO 102 148 : NUM(2,I)=IPOO 149 : 150 : 102 CONTINUE 151 : 152 : 153 : 110 DIN = DIN * APROG 154 : IF (ITYPEL.NE.3) GOTO 115 155 : IPOO=IPOO+1 156 : NUM(2,NBELEM)=IPOO 157 : TABPAR(NBELEM) = DIN 158 : 115 CONTINUE 159 : 160 : 161 : NUM(ITYPEL,NBELEM)=POINT2 162 : 163 : IADR=XCOOR(/1)/(IDIM+1) 164 : 165 : 166 : C Creation des points 167 : 168 : DPAR=0 169 : 170 : IF (NX.EQ.0) GOTO 300 171 : NBPTS=IADR+NX*(ITYPEL-1) 172 : 173 : 174 : SEGADJ MCOORD 175 : 176 : DO 301 I=1,NX 177 : 178 : DIN = TABPAR(I) 179 : DPAR=DPAR+DIN 180 : IF (ITYPEL.EQ.2) GOTO 1103 181 : DPAR1=DPAR-DIN*0.5D0 182 : 183 : PARA=ANGLE*DPAR1 184 : COSA=COS(PARA) 185 : SINA=SIN(PARA) 186 : 187 : XCOOR(IADR*(IDIM+1)+1)=XCEN+COSA*X1+SINA*X2 188 : XCOOR(IADR*(IDIM+1)+2)=YCEN+COSA*Y1+SINA*Y2 189 : IF (IDIM.EQ.3) THEN 190 : XCOOR(IADR*(IDIM+1)+3)=ZCEN+COSA*Z1+SINA*Z2 191 : ENDIF 192 : 193 : XCOOR((IADR+1)*(IDIM+1))=DENI+DECA*DPAR1 194 : IADR=IADR+1 195 : 196 : 1103 PARA=ANGLE*DPAR 197 : COSA=COS(PARA) 198 : SINA=SIN(PARA) 199 : 200 : XCOOR(IADR*(IDIM+1)+1)=XCEN+COSA*X1+SINA*X2 201 : XCOOR(IADR*(IDIM+1)+2)=YCEN+COSA*Y1+SINA*Y2 202 : IF (IDIM.EQ.3) THEN 203 : XCOOR(IADR*(IDIM+1)+3)=ZCEN+COSA*Z1+SINA*Z2 204 : ENDIF 205 : 206 : XCOOR((IADR+1)*(IDIM+1))=DENI+DECA*DPAR 207 : IADR=IADR+1 208 : 209 : 301 CONTINUE 210 : 211 : NBPTS = IADR 212 : 213 : 214 : SEGADJ MCOORD 215 : 216 : 217 : 300 DIN = TABPAR(NBELEM) 218 : DPAR=DPAR+DIN 219 : 220 : IF (ITYPEL.NE.3) GOTO 303 221 : 222 : NBPTS = IADR+1 223 : SEGADJ MCOORD 224 : DPAR1 = DPAR-DIN*0.5D0 225 : 226 : PARA = ANGLE * DPAR1 227 : COSA = COS(PARA) 228 : SINA = SIN(PARA) 229 : 230 : XCOOR(IADR*(IDIM+1)+1) = XCEN + X1*COSA + X2*SINA 231 : XCOOR(IADR*(IDIM+1)+2) = YCEN + Y1*COSA + Y2*SINA 232 : IF (IDIM.EQ.3) THEN 233 : XCOOR(IADR*(IDIM+1)+3) = ZCEN + Z1*COSA + Z2*SINA 234 : ENDIF 235 : 236 : XCOOR((IADR+1)*(IDIM+1)) = DENI+DECA*DPAR1 237 : IADR=IADR + 1 238 : 239 : 240 : 303 DO 309 I=1,NUM(/2) 241 : ICOLOR(I)=IDCOUL 242 : 309 CONTINUE 243 : SEGSUP TABPAR 244 : SEGDES MELEME 245 : 245 : CALL ECROBJ('MAILLAGE',MELEME) 246 : 247 : RETURN 248 : 249 : END 250 : 251 : 252 : 253 : 254 : 255 : 256 : 257 : 258 : 259 :
© Cast3M 2003 - All rights reserved.
Disclaimer