C REDUAF SOURCE FANDEUR 13/01/30 21:15:00 7686 subroutine reduaf(jchelm,mmodtm,mchel2,istri,iret,KERRE) implicit real*8(a-h,o-z) implicit integer (I-N) * on essaye de faire la reduction du champ par elemnt jchelm * sur le modele mmodel * resultat mchel2 iret=1 ok iret=0 rate *² INFCHE(.,6): LE SMINTE POINTE SUR UN SEGMENT D'INTEGRATION : * =1 AUX NOEUDS * =2 AU CENTRE DE GRAVITE * =3 AUX POINTS DE GAUSS POUR LA RIGIDITE * =4 AUX POINTS DE GAUSS POUR LA MASSE * =5 AUX POINTS DE GAUSS POUR LES CONTRAINTES * -INC CCOPTIO -INC SMCHAML -INC SMMODEL -INC SMCOORD c*-INC SMINTE -INC SMELEME SEGMENT INFO INTEGER INFELL(JG) ENDSEGMENT segment travg integer icpr(xcoor(/1)/(idim+1)) integer ityp (nbele) integer imin(nbele) integer ipart(nbele) integer ipos(nbele) endsegment segment inum(ino+1) segment ielno(ima) segment izone(nz) segment infc(nz) segment locglo(nmacom,nz) segment nokkom character*(8) nokomp(nbco) integer maxptv(nbco),maxpti(nbco) endsegment segment idimma integer imaxdi(n1) endsegment iret=0 KERRE = 0 mchelm=jchelm mmodel=mmodtm segact mchelm,mmodel segini,mchel2=mchelm nz=imache(/1) * Cas tres particulier : IF (nz.EQ.0) THEN * write(ioimp,*) 'CAS PARTICULIER nz = ',nz segdes,mchelm,mmodel * segdes,mchel2 iret = 1 RETURN ENDIF if (istri.eq.1) goto 1954 * test formulation melange: on ajoute les modeles de phases nvim=0 do i=1,kmodel(/1) imodel=kmodel(i) segact imodel if (formod(1).eq.'MELANGE') nvim=nvim+ivamod(/1) enddo if (nvim.ne.0) then nc=kmodel(/1) n1=nc+nvim segini mmode1 do i=1,kmodel(/1) mmode1.kmodel(i)=kmodel(i) enddo do 600 i=1,kmodel(/1) imodel=kmodel(i) if (formod(1).eq.'MELANGE') then do 601 j=1,ivamod(/1) if (tymode(j).ne.'IMODEL') call erreur(21) nc=nc+1 mmode1.kmodel(nc)=ivamod(j) 601 continue endif 600 continue mmodel=mmode1 endif 1954 continue * * test rapide sur les pointeurs des maillages on prend tout ce que l'on peut * n133=0 ia=0 segini izone do 500 i=1,kmodel(/1) imodel=kmodel(i) segact imodel * write (6,*) ' reduad model i imodel nefmod conmod ', * > i,imodel,nefmod,conmod if(nefmod.eq.22) then segdes imodel go to 500 endif do 501 j=1,nz * on parcourt tous les NZ chamelem elementaires * write (6,*) ' reduaf j imamod imache conmod conche ', * > j,imamod,imache(j),conmod,conche(j) if( imamod.eq.imache(j)) then * correspondance du maillage chamelem j avec maillage modele i if( conche(j).eq.conmod.or.conche(j).eq.' ') then ia=ia+1 if (ia.gt.mchel2.imache(/1)) then n1=ia+10 n133=ia l1=titche(/1) n3=infche(/2) segadj mchel2 endif * write (6,*) ' reduaf boucle 501 ia j conche ',ia,j,conche(j) mchel2.conche(ia)=conche(j) mchel2.imache(ia)=imamod mcham3=ichaml(j) segini,mcham2=mcham3 *? do 505 k=1,mcham2.ielval(/1) *? melva3=mcham2.ielval(k) *?** segini,melva2=melva3 *? melva2=melva3 *? mcham2.ielval(k)=melva2 *? segdes melva2,melva3 *? 505 continue mchel2.ichaml(ia)=mcham2 *? segdes mcham2,mcham3 segdes,mcham2 do 504 k=1,infche(/2) mchel2.infche(ia,k)=infche(j,k) 504 continue izone(j)=1 endif endif 501 continue 500 continue * les a-t-on tous pris? si oui, le resultat est le chamelem d'entree *** if(ia.eq.nz) then *** iret = 1 * desactivation des segments *** do 509 i=1,kmodel(/1) *** imodel=kmodel(i) *** segdes imodel **509 continue *** segdes mchelm,mmodel,mchel2 *** segsup izone *** return *** else * sinon, le resultat est un chamelem partiel if(ia.ne.0) then n1=ia n133=ia l1=titche(/1) n3=infche(/2) segadj mchel2 mchel3=mchel2 segini,mchel2=mchelm ib=0 do 506 j=1,nz if( izone(j).eq.1) go to 506 ib=ib+1 * write (6,*) ' reduaf boucle 506 ib j conche ',ib,j,conche(j) mchel2.conche(ib)=conche(j) mchel2.imache(ib)=imache(j) mchel2.ichaml(ib)=ichaml(j) do 507 k=1,infche(/2) mchel2.infche(ib,k)=infche(j,k) 507 continue 506 continue n1=ib l1=titche(/1) n3=infche(/2) segadj mchel2 segdes mchelm mchelm=mchel2 mchel2=0 else mchel3=0 endif *** endif * on doit tester element par element * travail sur les maillages support du mchaml * segsup izone * * le nouveau chmelem ne peut avoir plus que kmodel(/1) * *(nbinfele+conmod+nophas) *nbinfele est le nombre de infele different dans le chamelem a reduire * on commence par compter nbinfele nz=imache(/1) segini infc izo=0 n3=infche(/2) do 200 i=1,nz if(infc(i).ne.0) go to 200 izo=izo+1 infc(i)=izo do 201 j=i+1,nz do 202 k=1,n3 if(k.eq.4) go to 202 if(infche(i,k).ne.infche(j,k)) go to 201 202 continue * la zone j a le meme infche que la zone i * on teste aussi si meme constituant if(conche(i).ne.conche(j))go to 201 infc(j)=izo 201 continue 200 continue * * on va dire que chaque sous zone du modele a izo sous partie * donc nouveau chamelem a izo*kmodel(/1) sous partie : creation. * Au debut chaque sous mchel contient toutes les composantes n1=izo * kmodel(/1) segini,idimma l1=titche(/1) segini mchel2 mchel2.titche=titche mchel2.ifoche=ifoche do 205 i=1,kmodel(/1) imodel=kmodel(i) segact imodel MELE1=NEFMOD if( nefmod.eq.22) go to 205 lesupp=6 if(infele(6).ne.0) then infe6 = infele(6) else CALL ELQUOI(MELE1,0,lesupp,IPINF,IMODEL) info=ipinf segact info infe6 = infell(6) segsup info endif * write(6,*) 'infell',( infell( iou),iou=1,infell(/1)) na=(i-1)*izo do 206 j=1,izo *** imaxdi(na+j)=max(infell(6),infell(2)) * imaxdi(na+j)= infell(6) imaxdi(na+j)= infe6 mchel2.imache(na+j)=imamod mchel2.conche(na+j)=conmod * write (6,*) ' reduad boucle 206 j na conmod ',j,na,conmod 206 continue * segsup info 205 continue *on cree locglo(j,i)=k qui dit que la j ieme inconnue * de la i ieme zone du mchelm est en k ieme position globale * d'abord activation du mchelm et max de n2=nomche(/2)(=ielval(/1)) nz=imache(/1) nmacom=0 do 150 i=1,nz ipt1=imache(i) segact ipt1 mchaml=ichaml(i) segact mchaml nmacom=max(nmacom,ielval(/1)) do 151 j=1,ielval(/1) melval=ielval(j) * write(6,*) 'reduaf activation de melval ith',ith,melval segact melval 151 continue 150 continue nbco=20 nbcv=0 segini nokkom,locglo do 152 i=1,nz mchaml=ichaml(i) do 153 j=1,nomche(/2) do 154 k=1,nbcv if(nokomp(k).eq.nomche(j)) go to 155 154 continue nbcv=nbcv+1 if(nbcv.gt.nbco) then nbco=nbco+20 segadj nokkom endif nokomp(nbcv)=nomche(j) k=nbcv 155 continue melval=ielval(j) if( typche(j).eq.'REAL*8') then * write(6,*) 'i j k nokomp maxpt velche ', * $ i,j, k , nokomp(k),maxptv(k),velche(/1) maxptv(k)=max(maxptv(k),velche(/1)) else maxpti(k)=max(maxpti(k),ielche(/1)) endif locglo(j,i)=k 153 continue 152 continue nbco=nbcv segadj nokkom n2=nbco * * il reste a prendre chaque sous chamelem est a le dispachtcher * dans les differents bouts. Pour faire des tests rapides on cherche pour * chaque element dans les sous zones du model quel est le point de * plus petit numero et on construit un tableau qui dit que les elements * ayant un certain point mini i sont dans les zones j * nbsous=kmodel(/1) nbele=0 * on compte le nombre d'elements do 1 i=1,nbsous imodel=kmodel(i) if(nefmod.eq.22) go to 1 meleme=imamod segact meleme nbele= num(/2)+ nbele 1 continue segact mcoord segini travg * fabrication d'une numerotation locale icpr * fabrication de imin(k) qui donne le plus petit noeud de l'element * fabrication de ipart (k) qui dit dans quel sous paquet est l'element * fabrication de ityp qui donne le type de l element * fabrique de ipos(k): numero d'ordre de l'element dans le sous paquet idej=0 ino=0 do 2 i=1,nbsous imodel = kmodel(i) if( nefmod.eq.22) go to 2 meleme=imamod do 3 j=1,num(/2) ia = num(1,j) do 4 m=1,num(/1) ib=num(m,j) ia = min ( ia,ib) if(icpr(ib).eq.0) then ino=ino+1 icpr(ib)=ino endif 4 continue jidej= j + idej imin(jidej)=ia ipart(jidej)=i ityp(jidej)=itypel ipos(jidej)=j 3 continue idej = idej + num(/2) 2 continue segini inum * on va compter combien (au max)d'elements touchent un noeud mini do 5 i=1,nbele ia=imin(i) ib = icpr(ia) inum(ib)=inum(ib)+1 5 continue * pour chaque noeud, on stocke la fin de rangement dans ielno do i=1,ino inum(i+1)=inum(i+1)+inum(i) enddo ima=inum(ino+1) * write (6,*) ' reduaf ima',ima segini ielno * dans le tableau ielno on dira quels sont les elements ayant * le ino iem noeud local comme plus petit noeud do 7 i=1,nbele ia = imin(i) iloc=icpr(ia) ib=inum(iloc) inum(iloc)=inum(iloc)-1 ielno(ib)=i 7 continue * on a maintenant dans inum le debut de rangement dans ielno * * on travaille sous chaml par sou-chaml * do 8 i=1,imache(/1) ipt1=imache(i) segact ipt1 mchaml=ichaml(i) segact mchaml do 218 j=1,ielval(/1) melval=ielval(j) segact melval 218 continue ityref=ipt1.itypel c* indinf=infc(i) c* iconce= (i-1)*izo do 9 j=1,ipt1.num(/2) imi=ipt1.num(1,j) c* ina=0 do 10 k=2,ipt1.num(/1) imi= min(imi,ipt1.num(k,j)) 10 continue iloc = icpr(imi) if( iloc.eq.0) go to 9 do 11 k=inum(iloc)+1,inum(iloc+1) ie = ielno(k) if(ie.eq.0) write (6,*) ' remplissage incorrect de ielno ' * l element ie est un bon candidat on teste le itypel * puis le constituant if(ityref.ne.ityp(ie)) go to 11 nbs= ipart(ie) imodel=kmodel(nbs) meleme=imamod if(nefmod.eq.22) go to 11 iloo=ipos(ie) if(conche(i).ne.conmod.and.conche(i).ne.' ' ) go to 11 * on regarde si chaque noeud de j appartient aussi a ie do 12 m=1,ipt1.num(/1) if(num(m,iloo).ne.ipt1.num(m,j)) go to 11 12 continue * on a trouve une correspondance on peut dispatcher les valeurs ipla= ipos(ie) ires=(nbs-1)*izo+infc(i) if( mchel2.ichaml(ires) .eq.0) then n2=nbco segini mcham2 mchel2.ichaml(ires)=mcham2 do 21 m=1,n3 mchel2.infche(ires,m)=infche(i,m) 21 continue else mcham2=mchel2.ichaml(ires) endif do 16 m=1,ielval(/1) inoco=locglo(m,i) melval=ielval(m) if(mcham2.ielval(inoco).eq.0) then * nb de point support du chamelem iou=infche(i,6) if( iou.eq.1) then n2ptel=ipt1.num(/1) n1ptel=n2ptel elseif(iou.eq.2) then n2ptel=1 n1ptel=1 else n1ptel=imaxdi(ires) n2ptel=imaxdi(ires) endif ipt3=imamod * n1ptel=velche(/1) * n1ptel=imaxdi(ires) * n2ptel=ielche(/1) * n2ptel=imaxdi(ires) * write(6,*) ' n1ptel ', n1ptel if(typche(m).eq.'REAL*8') then n1el=ipt3.num(/2) n2el=0 n2ptel=0 else n2el=ipt3.num(/2) n1el=0 n1ptel=0 endif segini melva2 mcham2.ielval(inoco)=melva2 mcham2.typche(inoco)=typche(m) mcham2.nomche(inoco)=nokomp(inoco) else * verification du type if( mcham2.typche(inoco).ne.typche(m)) then * incompatibilite du type de composante entre champs KERRE = 472 endif melva2= mcham2.ielval(inoco) endif * write(6,*) ' typche(m)' ,typche(m) if(typche(m).eq.'REAL*8' ) then ** if(ipla.gt.melva2.velche(/2) ) then * write(6,*) ' erreur ipla dim2',ipla,melva2.velche(/2) * endif * if(velche(/1) . gt.melva2.velche(/1)) then * write(6,*) 'bizarre ***********',velche(/1),melva2.velche(/1) * endif ibj=min (velche(/2),j) * write(6,*) 'melva2.velche(/1)',melva2.velche(/1) do 17 n=1,melva2.velche(/1) ibn=min(velche(/1),n) melva2.velche(n,ipla)=velche(ibn,ibj) 17 continue else ibj=min (ielche(/2),j) do 18 n=1,melva2.ielche(/1) ibn=min(ielche(/1),n) melva2.ielche(n,ipla)=ielche(ibn,ibj) 18 continue endif 16 continue 11 continue 9 continue 8 continue * * suppression des segments de trvail * segsup locglo,infc,ielno segsup inum segsup travg segsup inum segsup nokkom * * compactage du resultat * n1=0 do 30 i=1,mchel2.imache(/1) if( mchel2.ichaml(i).ne.0) then n1 = n1+1 * write (6,*) ' reduad boucle 30 i n1 conche ',i,n1, * > mchel2.conche(i) mchel2.conche(n1)=mchel2.conche(i) mchel2.ichaml(n1)=mchel2.ichaml(i) mchel2.imache(n1)=mchel2.imache(i) do 31 j=1,n3 mchel2.infche(n1,j)=mchel2.infche(i,j) 31 continue mcham2=mchel2.ichaml(i) * on compacte les composantes n2=0 do 32 j=1,mcham2.ielval(/1) if( mcham2.ielval(j).ne.0) then n2=n2+1 mcham2.ielval(n2)=mcham2.ielval(j) mcham2.typche(n2)=mcham2.typche(j) mcham2.nomche(n2)=mcham2.nomche(j) melva2=mcham2.ielval(j) segdes melva2 endif 32 continue segadj mcham2 segdes mcham2 endif 30 continue * * desactivation des entrees :modele et chamelem * do 33 i=1,kmodel(/1) imodel=kmodel(i) ipt1=imamod segdes ipt1 segdes imodel 33 continue segdes mmodel do 34 i=1,imache (/1) ipt1=imache(i) segdes ipt1 mchaml=ichaml(i) do 35 j=1,ielval(/1) melval=ielval(j) *pv segdes melval 35 continue segdes mchaml 34 continue segdes mchelm * * ajout de mchel3 * if( n1.eq.0) then segsup mchel2 if(mchel3.eq.0) then * write(6,*) ' on tombe ici' KERRE = 21 return else mchel2=mchel3 go to 1111 endif endif n1dep=n1 n1=n1+n133 segadj mchel2 if(mchel3.ne.0) then segact mchel3 do 40 j=1,n133 mchel2.ichaml(n1dep+j)=mchel3.ichaml(j) mchel2.imache(n1dep+j)=mchel3.imache(j) mchel2.conche(n1dep+j)=mchel3.conche(j) * write (6,*) ' reduaf boucle 40 j n1dep conch', * > j,n1dep,mchel3.conche(j) do 41 k=1,mchel3.infche(/2) mchel2.infche(n1dep+j,k)=mchel3.infche(j,k) 41 continue 40 continue segsup mchel3 endif segdes mchel2 * * sortie * 1111 iret=1 * * on compacte le resultat * segact mchel2 do i=1,mchel2.ichaml(/1) mchaml=mchel2.ichaml(i) segact mchaml*mod * write (6,*) ' reduaf boucle finale i conche ', * > i,mchel2.conche(i) do j=1,ielval(/1) ielin=ielval(j) call comred(ielin) c* ielval(j)=ielin *pv melva1=ielin *pv segdes melva1 enddo segdes mchaml enddo *pv segdes mchel2 segsup idimma return end