1 : C REMPLX SOURCE CHAT 10/05/28 21:15:20 6677 2 : 3 : SUBROUTINE remplx(irig1,ipchp1,ipchp2 ) 4 : c==================================================================== 5 : c remplissage des LX apres reso avec condensation 6 : c entrees: 7 : c mrigid rigidité de dependance ( irigel(8,irig) # 0) 8 : C (elle contient les multiplicateurs) 9 : c ipchp1 chpoint de force cree par 10 : C (rigtot * soltot) - ftot 11 : c sorties: 12 : C ipchp2 chpoint des multiplicateurs a reintroduire 13 : C dans la solution globale 14 : C 15 : c==================================================================== 16 : * on suppose pour le moment qu'il n'y a qu'une seule harmonique de fourier 17 : 18 : IMPLICIT INTEGER(I-N) 19 : IMPLICIT REAL*8(A-H,O-Z) 20 : CHARACTER*4 CMOT 21 : -INC SMCHPOI 22 : -INC SMRIGID 23 : -INC SMELEME 24 : -INC SMCOORD 25 : -INC CCOPTIO 26 : -INC CCGEOME 27 : -INC CCHAMP 28 : character*4 NOMM,nom2 29 : C 30 : segment iccoun(ncomp) 31 : C 32 : segment iicpr 33 : integer icpr(nbpts),ivers(nbpts),iinb 34 : real*8 cvers(nbpts) 35 : endsegment 36 : segment siver2 37 : integer invnum(nbpts) 38 : integer iver2(nptc) 39 : real*8 cver2(nptc) 40 : endsegment 41 : SEGMENT SNOMIN 42 : CHARACTER*4 NOMIN(0) 43 : ENDSEGMENT 44 : SEGMENT SNOMDU 45 : CHARACTER*4 NODUA(0) 46 : ENDSEGMENT 47 : mrigid = irig1 48 : SEGACT mrigid 49 : nbpts = xcoor(/1)/(idim+1) 50 : segini snomin,snomdu 51 : C identification des inconnues liees 52 : DO 1501 I=1,IRIGEL(/2) 53 : MELEME=IRIGEL(1,I) 54 : C SEGACT MELEME 55 : DESCR=IRIGEL(3,I) 56 : ifo=irigel(5,i) 57 : SEGACT DESCR 58 : DO 1502 J=1,LISINC(/2) 59 : IF(LISINC(J).EQ.'LX '.AND.J.LE.1) GO TO 1502 60 : IF(NOMIN(/2).EQ.0) THEN 61 : NOMIN(**)=LISINC(J) 62 : NODUA(**)=LISDUA(J) 63 : ELSE 64 : DO 1506 K=1,NOMIN(/2) 65 : IF(LISINC(J).EQ.NOMIN(K)) GO TO 1505 66 : 1506 CONTINUE 67 : NOMIN(**)=LISINC(J) 68 : NODUA(**)=LISDUA(J) 69 : 1505 CONTINUE 70 : ENDIF 71 : 1502 CONTINUE 72 : SEGDES DESCR 73 : 1501 CONTINUE 74 : *C construction du tableau des duales 75 : * segini snomdu 76 : * do 325 il=1,nomin(/2) 77 : * NOMM =nomin(IL) 78 : * do 326 in = 1,lnomdd 79 : * if (NOMM.EQ.NOMDD(in)) then 80 : ** NODUA(il) =NOMDU(in) 81 : * goto 327 82 : * endif 83 : *326 continue 84 : *327 continue 85 : *325 continue 86 : 87 : C write(6,*) 'primales bloquees' ,(nomin(k),k=1,nomin(/2)) 88 : C write(6,*) 'leurs duales ' ,(nodua(k),k=1,nomin(/2)) 89 : 90 : ncomp = nomin(/2) 91 : segini iccoun 92 : 93 : do 2000 ic = 1,ncomp 94 : segini iicpr 95 : 96 : iccoun(ic) = iicpr 97 : C write(6,*) ' composante et son iicpr ',ic,iicpr 98 : 2000 continue 99 : 100 : ntplx = 0 101 : C quels sont les noeuds concernes par des liaisons 102 : do 1700 i=1,IRIGEL(/2) 103 : MELEME=IRIGEL(1,I) 104 : SEGACT MELEME 105 : ipt4=meleme 106 : descr = irigel(3,i) 107 : segact descr 108 : xmatri= irigel(4,i) 109 : segact xmatri 110 : if (lisinc(/2).le.1) goto 1701 111 : nomm = lisinc(2) 112 : C on l identifie par un numero ipo 113 : C pour recuperer les composantes ad hoc 114 : 115 : C quelle est sa position dans nomin 116 : 117 : do 3634 ipo=1,nomin(/2) 118 : if(NOMM.eq.nomin(ipo)) goto 3635 119 : 3634 continue 120 : 3635 continue 121 : C recup et comptage des noeuds supports de multiplicateurs 122 : ideb=2 123 : iicpr = iccoun(ipo) 124 : do 1703 iel=1,num(/2) 125 : * xmatri=imattt(iel) 126 : * segact xmatri 127 : ip =num(noelep(ideb),iel) 128 : * write (6,*) ' remplissage ip icpr(ip) ',ip,num(1,iel) 129 : icpr(ip)=num(1,iel) 130 : iinb = iinb+1 131 : ivers(iinb)= ip 132 : cvers(iinb)= re(1,2,iel)*coerig(i) 133 : ntplx = ntplx+1 134 : * write(6,*) 'ip icpr(ip) ipo' ,ip,icpr(ip),ipo 135 : * segdes xmatri 136 : 1703 continue 137 : 1701 continue 138 : segdes meleme,descr,xmatri 139 : 1700 continue 140 : 141 : * write(6,*) ' nombre de pts supports de LX ' ,ntplx 142 : C--------------------------------------------------------- 143 : * write(6,*) 'ivers ' ,(ivers(k),k=1,iinb) 144 : C on ouvre le chpoint des reactions calculees par KU-F 145 : mchpoi = ipchp1 146 : segact mchpoi 147 : nsoup1 = ipchp(/1) 148 : 149 : 150 : C initialisation du chpo de sortie support de ntplx points 151 : nat=1 152 : nsoupo = 1 153 : segini mchpo1 154 : ipchp2= mchpo1 155 : mchpo1.jattri(1) = 1 156 : mchpo1.ifopoi= ifopoi 157 : nc= 1 158 : segini msoup1 159 : mchpo1.ipchp(1)=msoup1 160 : C 161 : nbelem = ntplx 162 : nbnn = 1 163 : nbref = 0 164 : nbsous=0 165 : segini meleme 166 : msoup1.igeoc =meleme 167 : ipt1=meleme 168 : itypel= 1 169 : 170 : n = ntplx 171 : msoup1.nocomp(1) = 'LX ' 172 : msoup1.noharm(1)= ifo 173 : segini mpova1 174 : msoup1.ipoval = mpova1 175 : 176 : inpp = 0 177 : do 1000 isous=1,nsoup1 178 : msoupo= ipchp(isous) 179 : segact msoupo 180 : meleme = igeoc 181 : segact meleme 182 : nc = nocomp(/2) 183 : mpoval=ipoval 184 : segact mpoval 185 : nptc = vpocha(/1) 186 : do 1010 icomp=1,nc 187 : nom2= nocomp(icomp) 188 : C recherche de l inconnue duale 189 : do 334 ipo=1,nomin(/2) 190 : * write(6,*) ' nom de la composante et nodua ',nom2,NODUA(ipo) 191 : if(NOM2.eq.NODUA(ipo)) goto 335 192 : 334 continue 193 : C ce n est pas une composante depedante 194 : goto 1009 195 : 335 continue 196 : C write(6,*) 'composante ->vpocha(1) ipo ' ,nom2 ,nptc,ipo 197 : iicpr = iccoun(ipo) 198 : C boucle sur les noeuds 199 : C write(6,*) ' composante et son iicpr nbre ii ',ipo,iicpr,iinb 200 : segini siver2 201 : C write(6,*) 'num(1,il)',(num(1,il),il=1,nptc) 202 : 203 : C **************************************************** 204 : C [-DAVID-22/07/2004-] Inversion du tableau num suppression de la boucle 520 205 : * do 500 k=1,iinb 206 : * ip = ivers(k) 207 : * do 520 lk=1,nptc 208 : * if(ip.eq.num(1,lk)) then 209 : * iver2(lk)=ip 210 : * cver2(lk)=cvers(k) 211 : * goto 500 212 : * endif 213 : *520 continue 214 : *500 continue 215 : do lk = 1, nptc 216 : if ((num(1,lk).ge.1).and.(num(1,lk).le.nbpts)) then 217 : invnum(num(1,lk))=lk 218 : endif 219 : enddo 220 : do 500 k=1,iinb 221 : ip = ivers(k) 222 : lk = invnum(ip) 223 : if ((lk.ge.1).and.(lk.le.nptc)) then 224 : iver2(lk)=ip 225 : cver2(lk)=cvers(k) 226 : endif 227 : 500 continue 228 : 229 : C write(6,*) 'iver2 ' ,(iver2(k),k=1,nptc) 230 : 231 : do 1020 j=1,nptc 232 : ip = iver2(j) 233 : if(ip.eq.0.or.icpr(ip).eq.0.or.iver2(j).eq.0) goto 1020 234 : * on ne dicise plus par deux puisque simple multiplicateur 235 : xfo=-vpocha(j,icomp)/1.D0/cver2(j) 236 : * write (6,*) ' ip inpp icpr ',ip,inpp,icpr(ip) 237 : ipt1.num(1,inpp+1) = abs(icpr(ip)) 238 : * pour avoir les LX au total 239 : icpr(ip)=-abs(icpr(ip)) 240 : mpova1.vpocha(inpp+1,1)=xfo 241 : C write(6,*) 'ip icpr xfo ',ip,icpr(ip), 242 : C & xfo,lk,(inpp+1) 243 : inpp = inpp+1 244 : 1020 continue 245 : segsup siver2 246 : 1009 continue 247 : 1010 continue 248 : segdes mpoval 249 : segdes msoupo,meleme 250 : 1000 continue 251 : * on complete les LX 252 : do 1022 ic=1,nomin(/2) 253 : iicpr = iccoun(ic) 254 : do 1021 j=1,icpr(/1) 255 : if (icpr(j).le.0) goto 1021 256 : * write (6,*) ' complete j icpr(j) inpp num(/2) ', 257 : * > j,icpr(j),inpp,ipt1.num(/2) 258 : ipt1.num(1,inpp+1) = icpr(j) 259 : inpp=inpp+1 260 : 1021 continue 261 : 1022 continue 262 : * write(6,*) ' on a ecrit ',inpp, ' LX' 263 : segdes mchpoi 264 : nbnn=ipt1.num(/1) 265 : nbelem=inpp 266 : nbsous=0 267 : nbref=0 268 : * write (6,*) ' num mpova1 ',ipt1.num(/2),mpova1.vpocha(/1) 269 : if (nbelem.ne.ipt1.num(/2)) then 270 : * write (6,*) ' remplx ajustement LX ' 271 : segadj ipt1 272 : n=nbelem 273 : nc=1 274 : segadj mpova1 275 : endif 276 : segdes ipt1,mpova1,msoup1,mchpo1 277 : 278 : do 3000 ic=1,nomin(/2) 279 : iicpr = iccoun(ic) 280 : segsup iicpr 281 : 3000 continue 282 : segsup iccoun,snomin,snomdu 283 : segdes mrigid 284 : RETURN 285 : END 286 : 287 : 288 : 289 : 290 : 291 : 292 : 293 : 294 : 295 : 296 : 297 : 298 : 299 : 300 : 301 : 302 : 303 : 304 : 305 : 306 : 307 : 308 : 309 : 310 : 311 : 312 : 313 :
© Cast3M 2003 - All rights reserved.
Disclaimer