1 : C IMPOS1 SOURCE FANDEUR 12/03/15 21:26:25 7314 2 : 3 : C Operateur 'IMPO' 'MAIL' ( | 'MESC' | ) Lig1 Lig2 ( Crit) ; 4 : C 2D | 'FAIB' | 5 : C | 'SYME' | 6 : 7 : SUBROUTINE IMPOS1 8 : 9 : IMPLICIT INTEGER(I-N) 10 : IMPLICIT REAL*8 (A-H,O-Z) 11 : 12 : -INC CCOPTIO 13 : 14 : -INC SMCOORD 15 : -INC SMELEME 16 : 17 : PARAMETER ( X1s3 = 0.33333333333333333333333333333333333333333D0 , 18 : & X1s4 = 0.25000000000000000000000000000000000000000D0 ) 19 : 20 : logical faible 21 : PARAMETER (ncle = 3) 22 : character*4 mcle(ncle) 23 : data mcle / 'MESC','FAIB','SYME' / 24 : * 25 : idimp1 = IDIM + 1 26 : * 27 : * Lecture du mot eventuel SYME(trique) ou MESC ou FAIB(le) 28 : * 29 : mesc = 0 29 : call lirmot(mcle,ncle,mesc,0) 30 : if (ierr.ne.0) return 31 : * 32 : * Lecture des deux lignes en contact 33 : * 33 : call lirobj('MAILLAGE',ipt1,1,iretou) 34 : if (ierr.ne.0) return 34 : call lirobj('MAILLAGE',ipt2,1,iretou) 35 : if (ierr.ne.0) return 36 : * 37 : * Lecture eventuelle de la distance limite (crit) 38 : * 39 : crit=1.e+35 39 : call lirree(crit,0,iretou) 40 : if (ierr.ne.0) return 41 : cri2 = crit*crit 42 : * crit = abs(crit) 43 : *DBG write(ioimp,*) ' critere ',crit,cri2 44 : * 45 : * Petites verifications sur les lignes 46 : * 47 : segact,ipt1,ipt2 47 : if (ipt1.itypel.ne.2) call erreur(16) 47 : if (ipt2.itypel.ne.2) call erreur(16) 48 : if (ierr.ne.0) goto 900 49 : * 50 : * Quelques initialisations 51 : * 52 : nbpts0 = xcoor(/1) / idimp1 53 : * 54 : nel1 = ipt1.num(/2) 55 : nel2 = ipt2.num(/2) 56 : * noe1 = ipt1.num(/1) 57 : noe2 = ipt2.num(/1) 58 : * En fait : noe2 = noe1 = 2 car elements de type SEG2 ! 59 : * 60 : * Increment sur le nombre d'elements de contact retenus pour agrandir 61 : * les longueurs des segments (maillage, mcoord) de maniere adequate 62 : incnel = 10000 63 : * 64 : faible = mesc.eq.2 65 : if (mesc.ne.1) mesc = 0 66 : 67 : *======================================================================= 68 : * Formulation forte ("standard") du contact 69 : *======================================================================= 70 : if (faible) goto 200 71 : * 72 : * Creation du maillage de contact (meleme) 73 : * 74 : nbsous = 0 75 : nbref = 0 76 : nbnn = 4 77 : nbelem = incnel 78 : segini,meleme 79 : itypel = 22 80 : * 81 : * Nombre de points supports (multiplicateurs de Lagrange) associes 82 : * au contact : 1 par element de contact 83 : * Ajustement du segment mcoord 84 : * 85 : nbpts = nbpts0 + (incnel * 1) 86 : segadj,mcoord 87 : * 88 : nblag = nbpts0 89 : nbcoo = (nblag-1)*idimp1 90 : nbelc = 0 91 : * 92 : ipas = 1 93 : * 94 : 5 continue 95 : * 96 : * On boucle sur les elements de la 1ere ligne de contact 97 : * 98 : do 10 iel = 1, nel1 99 : nbeini = nbelc 100 : ip1=ipt1.num(1,iel) 101 : ip2=ipt1.num(2,iel) 102 : xp1=xcoor((ip1-1)*idimp1+1) 103 : yp1=xcoor((ip1-1)*idimp1+2) 104 : xp2=xcoor((ip2-1)*idimp1+1) 105 : yp2=xcoor((ip2-1)*idimp1+2) 106 : xg1 = xp1+xp2 107 : yg1 = yp1+yp2 108 : * 109 : * Boucle sur les noeuds de la 2eme ligne de contact 110 : * 111 : do 20 jel=1,nel2 112 : do 20 jnn=1,noe2 113 : jp=ipt2.num(jnn,jel) 114 : * verification que pas relation du point sur lui meme 115 : if (jp.eq.ip1) goto 20 116 : if (jp.eq.ip2) goto 20 117 : * verification que pas deja la relation 118 : do 21 irela = nbeini+1, nbelc 119 : if (jp.eq.num(4,irela)) goto 20 120 : 21 continue 121 : xp=xcoor((jp-1)*idimp1+1) 122 : yp=xcoor((jp-1)*idimp1+2) 123 : * test distance 124 : d1 = (xp1-xp)**2+(yp1-yp)**2 125 : d2 = (xp2-xp)**2+(yp2-yp)**2 126 : if (cri2.lt.d1.and.cri2.lt.d2) goto 20 127 : * 128 : * Mise a jour de la dimension des segments si necessaire 129 : * 130 : if (nbelc.eq.nbelem) then 131 : nbelem = nbelem + incnel 132 : segadj,meleme 133 : nbpts = nbpts + incnel 134 : segadj,mcoord 135 : endif 136 : * 137 : * on incremente mcoord avec le pt support du multiplicateur 138 : * on range dans le meleme l'element de contact associe 139 : * 140 : nblag = nblag + 1 141 : nbcoo = nbcoo + idimp1 142 : nbelc = nbelc + 1 143 : * 144 : xcoor(nbcoo+1) = (xg1+xp) * X1s3 145 : xcoor(nbcoo+2) = (yg1+yp) * X1s3 146 : xcoor(nbcoo+3) = 0. 147 : * 148 : num(1,nbelc) = nblag 149 : num(2,nbelc) = ip1 150 : num(3,nbelc) = ip2 151 : num(4,nbelc) = jp 152 : * 153 : 20 continue 154 : * 155 : 10 continue 156 : * 157 : * On recommence en inversant les deux lignes en cas de traitement 158 : * symetrique du contact (option MESC non activee) 159 : * 160 : if (ipas.eq.1.and.mesc.ne.1) then 161 : ipas=ipas+1 162 : iaux=ipt1 163 : ipt1=ipt2 164 : ipt2=iaux 165 : nel1 = ipt1.num(/2) 166 : nel2 = ipt2.num(/2) 167 : noe2 = ipt2.num(/1) 168 : goto 5 169 : endif 170 : * 171 : goto 1000 172 : 173 : *======================================================================= 174 : * Formulation "faible" du contact 175 : *======================================================================= 176 : 200 continue 177 : * 178 : * Le nb d'elements maximal a creer pour le maillage de contact FAIBLE : 179 : * 1 entre chaque element de chacune des lignes 180 : * 181 : * Creation du maillage de contact (meleme) 182 : * 183 : nbsous= 0 184 : nbref = 0 185 : nbnn = 5 186 : nbelem = incnel 187 : segini,meleme 188 : itypel = 22 189 : * 190 : * Nomre de points supports (multiplicateurs de Lagrange) associes 191 : * au contact : 1 par element de contact 192 : * Ajustement du segment mcoord 193 : * 194 : nbpts = nbpts0 + (nbelem * 1) 195 : segadj,mcoord 196 : * 197 : nblag = nbpts0 198 : nbcoo = (nblag-1)*idimp1 199 : nbelc = 0 200 : * 201 : * Boucle sur les elements de la 1ere ligne de contact 202 : * 203 : do 210 iel = 1, nel1 204 : ip1 = ipt1.num(1,iel) 205 : ip2 = ipt1.num(2,iel) 206 : ipv = (ip1-1)*idimp1 207 : xp1 = xcoor(ipv+1) 208 : yp1 = xcoor(ipv+2) 209 : ipv = (ip2-1)*idimp1 210 : xp2 = xcoor(ipv+1) 211 : yp2 = xcoor(ipv+2) 212 : xg1 = xp1 + xp2 213 : yg1 = yp1 + yp2 214 : * 215 : * Boucle sur les elements du maillage de la 2e ligne de contact 216 : * 217 : do 220 jel=1,nel2 218 : jp1 = ipt2.num(1,jel) 219 : jp2 = ipt2.num(2,jel) 220 : ipv = (jp1-1)*idimp1 221 : xp3 = xcoor(ipv+1) 222 : yp3 = xcoor(ipv+2) 223 : ipv = (jp2-1)*idimp1 224 : xp4 = xcoor(ipv+1) 225 : yp4 = xcoor(ipv+2) 226 : 227 : * test distance 228 : d1 = (xp1-xp3)**2+(yp1-yp3)**2 229 : d2 = (xp2-xp3)**2+(yp2-yp3)**2 230 : if (d1.le.cri2.or.d2.le.cri2) goto 250 231 : d1 = (xp1-xp4)**2+(yp1-yp4)**2 232 : d2 = (xp2-xp4)**2+(yp2-yp4)**2 233 : if (d1.le.cri2.or.d2.le.cri2) goto 250 234 : goto 220 235 : * 236 : 250 continue 237 : * 238 : * Mise a jour de la dimension des segments si necessaire 239 : * 240 : if (nbelc.eq.nbelem) then 241 : nbelem = nbelem + incnel 242 : segadj,meleme 243 : nbpts = nbpts + incnel 244 : segadj,mcoord 245 : endif 246 : * 247 : * on incremente mcoord avec les pts support des multiplicateurs 248 : * on range dans le meleme l'element de contact correspondant 249 : * 250 : nblag = nblag + 1 251 : nbcoo = nbcoo + idimp1 252 : nbelc = nbelc + 1 253 : * 254 : xcoor(nbcoo+1) = (xg1+xp3+xp4) * X1s4 255 : xcoor(nbcoo+2) = (yg1+yp3+yp4) * X1s4 256 : xcoor(nbcoo+3) = 0. 257 : * 258 : num(1,nbelc) = nblag 259 : num(2,nbelc) = ip1 260 : num(3,nbelc) = ip2 261 : num(4,nbelc) = jp1 262 : num(5,nbelc) = jp2 263 : * 264 : 220 continue 265 : * 266 : 210 continue 267 : * 268 : * goto 1000 269 : * 270 : *======================================================================= 271 : * Fin 272 : *======================================================================= 273 : 1000 continue 274 : *DBG write(ioimp,*) ' nbpts nbpts0 nblag',nbpts,nbpts0,nblag 275 : *DBG write(ioimp,*) ' nbelem nbelc',nbelem,nbelc 276 : * 277 : * On renvoie le resultat mais au prealable on ajuste les longueurs 278 : * 279 : if (nbelem.lt.nbelc) then 279 : call erreur(5) 280 : segsup,meleme 281 : goto 900 282 : else if (nbelem.gt.nbelc) then 283 : nbelem = nbelc 284 : segadj,meleme 285 : nbpts = nblag 286 : segadj,mcoord 287 : endif 288 : segdes,meleme 288 : call ecrobj('MAILLAGE',meleme) 289 : 290 : 900 continue 291 : segdes,ipt1,ipt2 292 : 293 : return 294 : end 295 : 296 :
© Cast3M 2003 - All rights reserved.
Disclaimer