ecou29
C ECOU29 SOURCE CB215821 24/04/12 21:15:42 11897 1 NUMAT,NUCAR,IVASTR,IVARI,IVADEF,IVADET, 1 IVADS,IVAMAT,IVACAR, 2 IPH1,IPH2,IPH3,ITHHER,LHOOK,NSTRS,NVARI,NMATT,NCARR, 3 CMATE,PRECIS,JECHER,IPOTAB,ISTEP,NPINT,JNOID,LOGSUC, *********************************************************************** * ecoulement inelastique appele par ecoul1 c ppu modif pour les materiaux unidirectionels en plastique * MATERIAUX ENDOMMAGEABLES DE LEMAITRE *********************************************************************** * entrees : * * mate = numero de materiau elastique * inplas = numero de materiau inelastique * mele = numero element fini * ipmail = pointeur du maillage * nbptel = nombre de points par element * imat = pointeur sur un segment mptval de materiau (utilise par calsig) * icar = pointeur sur un segment mptval de caracteristiques * geometriques (utilise par calsig) * numat = nb de composantes du melval de imat * nucar = nb de composantes du melval de icar * ivastr =pointeur sur un segment mptval de contraintes * ivari =pointeur sur un segment mptval de variables internes * ivadef =pointeur sur un segment mptval de deformations * ivads =pointeur sur un segment mptval de contraintes (increments) * ivamat =pointeur sur un segment mptval de materiau * ivacar =pointeur sur un segment mptval de cacarteristiques geometrique * iph1 = pointeur sur un mchaml de temperatures au debut du pas * iph2 = pointeur sur un mchaml de temperatures a la fin du pas * iph3 = pointeur sur un mchaml de temperatures de reference * ithher = 0 si pas de chargement thermique * = 1 si chargement thermique mais materiau constant * = 2 si chargement thermique et mat. dependant de la temperature * ipch1,ipch2,ipch3,ithher ne servent que pour les materiaux * endommageables de lemaitre quand ils dependent de la temperature * lhook =taille de la matrice de hooke * nstrs =nombre de composantes de contraintes * nvari =nombre de composantes de variables internes * nmatt =nombre de composnates de proprietes de materiau * ncarr =nombre de composnates de caracteristiques geometriques * cmate =nom du materiau * precis =precision dans les iterations internes * jecher =0 ou 1 pour action dans ecoule * jnoid =0 ou 1 pour action dans ecoule * ipotab =pointeur sur segment table * istep =indicateur d'action pour calcul nonlocal * =0 dans le cas d'un calcul local (normal) * =1 ou 2 dans le cas d'un calcul nonlocal * =1 pour calcul des fonctions seuil uniquement * =2 pour calcul des variables dissipatives a partir * des fonctions seuil moyennees prealablement par nloc * * sorties : * ivastf =pointeur sur un segment mptval de contraintes * ivarif =pointeur sur un segment mptval de variables internes * ivadep =pointeur sur un segment mptval de deformations inelastiques * kerre =indicateur d'erreur * * p dowlatyari fev. 1992 * * c. la borderie fev 92 restructuration et reecriture de certains * passages pour une meilleure lisibilite * * avril 92 ajout istep pour le non local * dec 92 modif pour poutres timoschenko * ************************************************************************ IMPLICIT INTEGER(I-N) IMPLICIT REAL*8(A-H,O-Z) * -INC PPARAM -INC CCOPTIO -INC SMCHAML -INC SMELEME -INC SMCOORD -INC SMMODEL -INC SMINTE -INC CCHAMP -INC CECOU c======================================================================= c la variable kerre regit les impressions d erreurs dans plast c toutes erreurs de ecoule gerees dans ce sous programme c kerre=0 tout ok c de 1 a 6 s aligner sur valeurs donnees par ecoinc c = 7 un element tuyau a une epaisseur nulle c = 21 on ne trouve pas d intersection avec la surface de charge c = 22 sig0 a l exterieur de la surface de charge c c anomalies avec la courbe de traction c = 30 limite elastique nulle c = 31 trop de points c = 32 pas assez de points c = 33 pente incorrecte c = 34 module d'young nul c = 35 manque l'origine c = 36 pente a l'origine non egale a e c = 37 manque la courbe de traction c = 38 nu devrait etre nul c c = 48 donnees erronnees pour drucker-prager c = 49 matrice singuliere dans iter internes drucker-prager c = 51 pb dans drucker prager option non disponible c = 52 pb dans drucker prager donnees incompatibles c = 53 pb dans drucker prager solution impossible c = 54 les valeurs admissibles pour istep sont 0 1 ou 2 c = 55 modele non implante en non local c = 56 probleme dans l'integration du modele mazars c = 57 .... c = 58 .... c = 59 .... c = 60 pb donnees du cam-clay c c = 99 cas non encore disponible c======================================================================= * SEGMENT MPTVAL INTEGER IPOS(NS) ,NSOF(NS) INTEGER IVAL(NCOSOU) CHARACTER*16 TYVAL(NCOSOU) ENDSEGMENT * SEGMENT WRK0 REAL*8 XMAT(NCXMAT) ENDSEGMENT * SEGMENT WR00 CHARACTER*16 TYMAT(NCXMAT) REAL*8 XMAT1(NCXMAT),XMAT2(NCXMAT) ENDSEGMENT * SEGMENT WRK1 REAL*8 DDHOOK(LHOOK,LHOOK),SIG0(NSTRS),DEPST(NSTRS) REAL*8 SIGF(NSTRS),VAR0(NVARI),VARF(NVARI) REAL*8 DEFP(NSTRS),XCAR(ICARA) ENDSEGMENT * SEGMENT WRK2 REAL*8 TRAC(LTRAC) ENDSEGMENT * SEGMENT WRK22 REAL*8 XXE(3,NBNN) ENDSEGMENT * SEGMENT WRK3 ENDSEGMENT * SEGMENT WRK4 REAL*8 XE(3,NBBB) ENDSEGMENT * SEGMENT WRK5 REAL*8 EPIN0(NSTRS),EPINF(NSTRS),EPST0(NSTRS) ENDSEGMENT * SEGMENT WRK6 REAL*8 BB(NSTRS,NNVARI),R(NSTRS),XMU(NSTRS) REAL*8 S(NNVARI),QSI(NNVARI),DDR(NSTRS),BBS(NSTRS) ENDSEGMENT * SEGMENT WRK7 REAL*8 F(NCOURB,2),W(NCOURB),TRUC(NCOURB) ENDSEGMENT * SEGMENT WRK8 REAL*8 DD(NSTRS,NSTRS),DDV(NSTRS,NSTRS),DDINV(NSTRS,NSTRS) ENDSEGMENT * SEGMENT WRK9 REAL*8 YOG(NYOG),YNU(NYNU),YALFA(NYALFA),YSMAX(NYSMAX) REAL*8 YN(NYN),YM(NYM),YKK(NYKK),YALFA1(NYALF1) REAL*8 YBETA1(NYBET1),YR(NYR),YA(NYA),YKX(NYKX),YRHO(NYRHO) INTEGER NKX(NNKX) ENDSEGMENT * SEGMENT WR10 INTEGER IABLO1(NTABO1) REAL*8 TABLO2(NTABO2) ENDSEGMENT * SEGMENT WR11 INTEGER IABLO3(NTABO3) REAL*8 TABLO4(NTABO4) ENDSEGMENT * SEGMENT WTRAV REAL*8 DDAUX(LHOOK,LHOOK),VALMAT(NUMAT) REAL*8 VALCAR(NUCAR),DSIGT(NSTRS) REAL*8 TXR(IDIM,IDIM),DDHOMU(LHOOK,LHOOK) REAL*8 XLOC(3,3),XGLOB(3,3) REAL*8 D1HOOK(LHOOK,LHOOK),ROTHOO(LHOOK,LHOOK) ENDSEGMENT * SEGMENT WPOUT REAL*8 X(2),Y(2),Z(2) ENDSEGMENT * LOGICAL LOGVIS,LOGIN,LOGRE,LOGSUC LOGICAL LUNI1,LUNI2 DIMENSION BID(6),BID2(6),CRIGI(12),CMASS(12) * CHARACTER*72 CHARRE CHARACTER*8 CMATE * * * mise à disposition des temperatures tini tfin tref * aux points de gauss * TETA1=-1.E35 TETA2=-1.E35 TETREF=-1.E35 TREFA=-1.E35 IF (ITHHER.EQ.1.OR.ITHHER.EQ.2) THEN MCHAM3=IPH1 MCHAM4=IPH2 MCHAM5=IPH3 SEGACT MCHAM3 SEGACT MCHAM4 SEGACT MCHAM5 MELVA3=MCHAM3.IELVAL(1) MELVA4=MCHAM4.IELVAL(1) MELVA5=MCHAM5.IELVAL(1) SEGACT MELVA3 SEGACT MELVA4 SEGACT MELVA5 ENDIF c c Initialisations de variables c--------------------------------- c - mise à zéro des variables du commun NECOU si besoin c - modèles viscoplastiques: c . on récupère le pas de temps c . on récupère le nombre maximal de sous-pas c . on met IND=1 c - initialisation des dimensions des tableaux des segments c Sorties: en plus du commun NECOU, on range les autres données c initialisées dans les COMMON IECOU et XECOU c Sauf pour KERRE,LW,LOGVIS,LUNI1 et LUNI2 qui sont sortis comme c argument de DEFINI c . ISTEP,INPLAS,NPINT,IPOTAB,IVADEF, . IPMAIL,IVAMAT, . ITHHER,NUMAT,NUCAR,LOGVIS, . LUNI1,LUNI2,LW,KERRE) IF (KERRE.EQ.999) RETURN c c Initialisations des segments de travail c IF ((CMATE.EQ.'ORTHOTRO'.OR.CMATE.EQ.'ANISOTRO'.OR. 1 CMATE.EQ.'UNIDIREC').AND.(MFR.EQ.1.OR.MFR.EQ.31 1 .OR.MFR.EQ.33)) THEN MINTE2=IPTR1 SEGACT MINTE2 SEGINI WRK22 ENDIF c IF (LOGVIS) SEGINI WRK8 * * initialisation des segments de travail * SEGINI WRK0,WR00,WRK1,WRK2,WRK3,WRK5 IF(MFR.EQ.7.OR.MFR.EQ.13.OR.LUNI1)THEN SEGINI WRK4 ENDIF IF(INPLAS.EQ.26)THEN SEGINI WRK6 ENDIF c SEGINI WTRAV * * * boucle sur les elements * DO 1000 IB=1,NBELEM * * Matériaux orthotropes, anisotropes et unidirectionnels * en formulation massive: * - on cherche les coordonnees des noeuds de l element ib * - calcul des axes locaux * Cas particulier de l'ACIER_UNI * . MELEME,WRK4,WRK22,WTRAV) * * * boucle sur les points de gauss * DO 1100 IGAU=1,NBPTEL * * -recuperation de valmat et de valcar * -on recupere les contraintes initiales * -on recupere les variables internes * -on recupere les deformations inelastiques initiales si besoin * -on recupere les increments de deformations totales * -on cherche la section de l'element ib * -prise en compte de l'epaisseur et de l'excentrement * dans le cas des coques minces avec ou sans cisaillement * transverse * . IMAT,IVACAR,ICAR,IVASTR,IVARI,IVADEF,IVADET, . IVADS,MFR,CMATE,INPLAS,IB,IGAU,IND, . WTRAV,WRK1,WRK5,SECT,EPAIST) * * on recupere les constantes du materiau * calcul des contraintes effectives en milieu poreux * . IVAMAT,IB,IGAU,CMATE,MATE,LUNI1,LUNI2, . WRK1,WRK5,WRK0,WR00,WTRAV,CMASS,CRIGI,COB,XMOB, . BID,BID2,KERR0) IF (KERR0.EQ.99) THEN KERRE=99 GOTO 1000 ELSE IF (KERR0.EQ.10) THEN GOTO 1000 ENDIF * IF ((INPLAS.EQ.29.).OR.(INPLAS.EQ.26)) THEN * * pour les materiaux endommageables de lemaitre traitement special * car ils peuvent dependre de la temperature * NTABO1 = 0 NTABO2 = 0 SEGINI WR10 DO 2200 JC=1,NMATT IF (TYMAT(JC)(1:8).EQ.'REAL*8 ') THEN NTABO1=NTABO1+1 NTABO2=NTABO2+1 SEGADJ WR10 IABLO1(NTABO1)=1 TABLO2(NTABO2)=XMAT(JC) ELSE IF (TYMAT(JC)(9:16).EQ.'EVOLUTIO') THEN IF (KERRE.NE.0) GOTO 1990 NTABO1=NTABO1+1 NTABO=NTABO2 NTABO2=NTABO2+(2*NCOURB) SEGADJ WR10 IABLO1(NTABO1)=2*NCOURB DO 2050 JCC=1,NCOURB TABLO2(NTABO+(2*JCC-1))=TRAC(2*JCC-1) TABLO2(NTABO+(2*JCC))=TRAC(2*JCC) 2050 continue ELSE IF (TYMAT(JC)(9:16).EQ.'NUAGE ') THEN NTABO3 = 0 NTABO4 = 0 SEGINI WR11 IF (KERRE.NE.0) THEN SEGSUP WR10 SEGSUP WR11 KERR1=2 GOTO 1990 ENDIF * segadj wr11 NTABO=NTABO1 NTABOO=NTABO2 NTABO1=NTABO1+NTABO3+1 NTABO2=NTABO2+NTABO4 SEGADJ WR10 IABLO1(NTABO+1)=NTABO3 DO 2075 JCC=1,NTABO3 2075 iablo1(ntabo+1+jcc)=iablo3(jcc) DO 2125 JCC=1,NTABO4 2125 tablo2(ntaboo+jcc)=tablo4(jcc) SEGSUP WR11 ENDIF 2200 continue ENDIF * * >>>>>>>>>> fin du traitement du materiau * * on recupere les caracteristiques geometriques * . WRK1) * CALL DEFCAR(NCARR,ICARA,IB,IGAU,MFR,MELE,IVACAR, * . XCAR) * * * quelques impressions si iimpi = 99 * * IF(IIMPI.EQ.99) THEN * WRITE(IOIMP,66770) IB,IGAU *66770 format(////////2x,'element ',i6,2x,'point ',i3//) * WRITE(IOIMP,66771) MATE,INPLAS *66771 format('0 mate=',i4,2x,'inplas=',i4/) * WRITE(IOIMP,66772) (SIG0(I),I=1,NSTRS) *66772 format(2x,' sig0 '/(6(1x,1pe12.5))) * WRITE(IOIMP,66773) (VAR0(I),I=1,NVARI) *66773 format(2x,' var0 '/(6(1x,1pe12.5))) * WRITE(IOIMP,66774) (DEPST(I),I=1,NSTRS) *66774 format(2x,' depst '/(6(1x,1pe12.5))) * WRITE(IOIMP,66775) (XMAT(I),I=1,NMATT) *66775 format(2x,' xmat '/(6(1x,1pe12.5))) * IF(IVACAR.NE.0)THEN * WRITE(IOIMP,66776) (XCAR(I),I=1,ICARA) *66776 format(2x,' xcar '/(6(1x,1pe12.5))) * ENDIF * ENDIF * * mise à disposition des temperatures tini tfin tref * aux points de gauss * IF (ITHHER.EQ.1.OR.ITHHER.EQ.2) THEN IBMN=MIN(IB,MELVA3.VELCHE(/2)) IGMN=MIN(IGAU,MELVA3.VELCHE(/1)) TETA1=MELVA3.VELCHE(IGMN,IBMN) IBMN=MIN(IB,MELVA4.VELCHE(/2)) IGMN=MIN(IGAU,MELVA4.VELCHE(/1)) TETA2=MELVA4.VELCHE(IGMN,IBMN) IBMN=MIN(IB,MELVA5.VELCHE(/2)) IGMN=MIN(IGAU,MELVA5.VELCHE(/1)) TETREF=MELVA5.VELCHE(IGMN,IBMN) ENDIF * * *--------------------------------------------------------------------- * * ecoulement selon les modeles * *--------------------------------------------------------------------- * * * modeles de viscoplasticite integres par consti * IF ( INPLAS .EQ. 29) THEN * NYOG=IABLO1(1) NYNU=IABLO1(2) NYALFA=IABLO1(3) NYSMAX=IABLO1(4) NYN=IABLO1(5) NYM=IABLO1(6) NYKK=IABLO1(7) NYALF1=IABLO1(8) NYBET1=IABLO1(9) NYR=IABLO1(10) NYA=IABLO1(11) IF ((MFR.EQ.1.OR.MFR.EQ.31.OR.MFR.EQ.33). + AND.IFOUR.EQ.-2) THEN INTMAT=15 ELSE INTMAT=14 ENDIF IF (NTABO1.EQ.INTMAT) THEN NNKX=1 NYKX=IABLO1(12) ELSE NNKX=IABLO1(12) NYKX=0 DO 1881 I=1,NNKX 1881 NYKX=NYKX+(2*IABLO1(12+I)) NYKX=NYKX+NNKX ENDIF NYRHO=IABLO1(NTABO1) NSIGY=1 SEGINI WRK9 SEGSUP WR10 IF (ITHHER.EQ.0.OR.ITHHER.EQ.1) THEN NCOURB=2*NKX(1) ELSE NCOURB=NKX(1) DO 1882 I=1,NNKX 1882 IF (NKX(I).GE.NCOURB) NCOURB=NKX(I) NCOURB=2*NCOURB ENDIF SEGINI WRK7 IF (INPLAS.EQ.29.AND.VAR0(3).GE.0.96) THEN DO 1883 I=1,NVARI VARF(I) = VAR0(I) 1883 CONTINUE VARF(3) = 1.0 DO 1884 I=1,NSTRS EPINF(I) = EPIN0(I) 1884 CONTINUE SEGSUP WRK7 SEGSUP WRK9 ELSE 1 INPLAS,MFR1,DT,NSTRSS,NVARI,NMATT,PRECIS,MSOUPA,JECHER,DTT, 2 NSSINC,INV,KERRE,ICARA,IFOURB,NYOG,NYNU,NYALFA,NYSMAX,NYN, 3 NYM,NYKK,NYALF1,NYBET1,NYR,NYA,NYKX,NNKX,NYRHO,NSIGY,TETA1, 5 TETA2,TREFA,TLIFE,ITHHER,NCOURB,CMATE,N2EL,N2PTEL,IB,IGAU, 6 EPAIST,NBPGAU,MELE,NPINT,NBGMAT,NELMAT,SECT,LHOOK,CRIGI,KERREU1) c c SEGSUP WRK7 SEGSUP WRK9 IF (INPLAS.EQ.29.AND.TLIFE.GE.0.D0) THEN INTERR(1)=IB INTERR(2)=IGAU REAERR(1)=TLIFE ENDIF DTOPTI = MIN(DTOPTI,DTT) NINCMA = MAX(NINCMA,NSSINC) NCOMP = NCOMP + 1 TSOM = TSOM + DTT NSOM = NSOM + NSSINC TCAR = TCAR + DTT* DTT IF(KERRE.NE.0) THEN KERR1=1 END IF END IF c c modele plastique d'endommagement de lemaitre c ++++++++++++++++++++++++++++++++++++++++++++ c traitement du materiau qui depend eventuellement de la temperature c ------------------------------------------------------------------ ELSE IF (INPLAS.EQ.26) THEN NYOG=IABLO1(1) NYNU=IABLO1(2) NYRHO=IABLO1(3) NYALFA=IABLO1(4) IF ((MFR.EQ.1.OR.MFR.EQ.31.OR.MFR.EQ.33). + AND.IFOUR.EQ.-2) THEN INTMAT=10 ELSE INTMAT=9 ENDIF IF (NTABO1.EQ.INTMAT) THEN NNKX=1 NYKX=IABLO1(5) IEPS=0 ELSE NNKX=IABLO1(5) NYKX=0 DO 1789 I=1,NNKX 1789 NYKX=NYKX+(2*IABLO1(5+I)) NYKX=NYKX+NNKX IEPS=1 ENDIF IORIGI=6+(IEPS*NNKX) NYN=IABLO1(IORIGI) NYM=IABLO1(IORIGI+1) NYKK=IABLO1(IORIGI+2) NYSMAX=0 NYALF1=0 NYBET1=0 NYR=0 NYA=0 NSIGY=0 SEGINI WRK9 SEGSUP WR10 c c *** si le pt. de gauss est ruine, les contr. sont annulees et c *** on n' ecoule pas c IF (VAR0(3).GE.1.D0.OR.VAR0(3).GE.DC) THEN DO 1115 IEN=1,NVARI VARF(IEN)=VAR0(IEN) 1115 continue VARF(3)=1.D0 SEGSUP WRK9 ELSE c ---------------------------------------------------------------------- c nnvari est le nbr. de var. int. pilotant les eq. du modele soit r et d c p est en supplement c ---------------------------------------------------------------------- NNVARI=2 IF (ITHHER.EQ.0.OR.ITHHER.EQ.1) THEN NCOURB=2*NKX(1) ELSE NCOURB=NKX(1) DO 1119 I=1,NNKX 1119 if (nkx(i).ge.ncourb) ncourb=nkx(i) NCOURB=4*NCOURB ENDIF IF (KERRE.EQ.0) THEN SEGINI WRK7 1 NMATT,ICARA,INPLAS,NVARI,PRECIS,MFR1,IFOURB,KERRE,NNVARI, 2 NYOG,NYNU,NYRHO,NYALFA,NNKX,NYKX,NCOURB,NYN,NYM,NYKK,TETA1, 3 TETA2,TREFA,ITHHER,CMATE,N2EL,N2PTEL,IB,IGAU,EPAIST, 4 MELE,NPINT,NBGMAT,NELMAT,SECT,LHOOK,CRIGI) SEGSUP WRK7 SEGSUP WRK9 IF(KERRE.GT.200) THEN KERR1=1 END IF END IF END IF ELSE KERRE = 99 ENDIF * * Erreurs * - problèmes de convergence * * * - autres problèmes * . KERR1,KERRE) 1998 IF (KERRE.GT.0) THEN IF (LOGVIS) SEGSUP WRK8 SEGSUP WRK0,WR00,WRK1,WRK2,WRK3,WRK5,WTRAV IF(MFR.EQ.7.OR.MFR.EQ.13.OR.LUNI1) THEN SEGSUP WRK4 ENDIF IF ((CMATE.EQ.'ORTHOTRO'.OR.CMATE.EQ.'ANISOTRO'.OR. 1 CMATE.EQ.'UNIDIREC').AND.(MFR.EQ.1.OR.MFR.EQ.31 1 .OR.MFR.EQ.33)) THEN SEGDES MINTE2 SEGSUP WRK22 ENDIF IF (ITHHER.EQ.1.OR.ITHHER.EQ.2) THEN SEGDES MELVA3 SEGDES MELVA4 SEGDES MELVA5 SEGDES MCHAM3 SEGDES MCHAM4 SEGDES MCHAM5 ENDIF RETURN ENDIF c c c remplissage du segment contenant les contraintes a la fin * ( rearrangement pour milieu poreux ), c les variables internes finales c et les increments de deformations plastiques c . INPLAS,IND,WRK1,WRK5,WTRAV, . IVASTF,IVARIF,IVADEP,COB,XMOB,IB,IGAU, . CMATE,MATE,MELE,KERRER) IF (KERRER.GT.0) GOTO 1000 c c c fin de la boucle sur les points de gauss c 1100 continue c c special poutres et tuyaux sauf timoschenko c c c fin de la boucle sur les elements c 1000 continue c * FIN: modèles visqueux, on stocke le pas de temps * optimal en indice 'dtopti' * . TCAR,DTOPTI,IPOTAB,KERRE) IF (LOGVIS) SEGSUP WRK8 * * SEGSUP WRK0,WR00,WRK1,WRK2,WRK3,WRK5,WTRAV IF(MFR.EQ.7.OR.MFR.EQ.13.OR.LUNI1) THEN SEGSUP WRK4 END IF IF(INPLAS.EQ.26) THEN SEGSUP WRK6 SEGSUP WRK8 END IF IF ((CMATE.EQ.'ORTHOTRO'.OR.CMATE.EQ.'ANISOTRO'.OR. 1 CMATE.EQ.'UNIDIREC').AND.(MFR.EQ.1.OR.MFR.EQ.31 1 .OR.MFR.EQ.33)) THEN SEGDES MINTE2 SEGSUP WRK22 ENDIF * IF (ITHHER.EQ.1.OR.ITHHER.EQ.2) THEN SEGDES MELVA3 SEGDES MELVA4 SEGDES MELVA5 SEGDES MCHAM3 SEGDES MCHAM4 SEGDES MCHAM5 ENDIF * RETURN END
© Cast3M 2003 - Tous droits réservés.
Mentions légales