C ECOU20 SOURCE OF166741 25/11/04 21:15:46 12349 SUBROUTINE ECOU20(MATE,INPLAS,MELE,IPMAIL,NBPTEL,IMAT,ICAR, 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, 4 N2EL,N2PTEL,NBNO,NBPGAU,LW,IVASTF,IVARIF,IVADEP,KERRE) *********************************************************************** * ecoulement inelastique appele par ecoul1 c ppu modif pour les materiaux unidirectionels en plastique * MATERIAUX: - VISCOPLASTIQUES ET FLUAGE INTEGRES PAR CONSTI *********************************************************************** * 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 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 REAL*8 WORK(LW),WORK2(LW2) 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) REAL*8 SIGMA(NSTRS),SIGGD(NSTRS),XMULT(NSTRS),PROD(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) REAL*8 SIGY(NSIGY) 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--------------------------------- WRK8 = 0 WRK22 = 0 minte2 = 0 WRK4 = 0 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 CALL DEFINI(MELE,NCARR,NSTRS,NMATT,CMATE,MATE, . ISTEP,INPLAS,NPINT,IPOTAB,IVADEF, . IPMAIL,IVAMAT, . ITHHER,NUMAT,NUCAR,LOGVIS, . LUNI1,LUNI2,LW,KERRE) IF (KERRE.EQ.999) RETURN 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 CALL RESHPT(1,NBNO,IELE,MELE,NPINT,IPTR1,IRT1) MINTE2=IPTR1 SEGINI WRK22 ENDIF IF (LOGVIS) SEGINI WRK8 SEGINI WRK0,WR00,WRK1,WRK2,WRK3,WRK5 IF(MFR.EQ.7.OR.MFR.EQ.13.OR.LUNI1)THEN SEGINI WRK4 ENDIF 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 * CALL DEFROT(CMATE,MFR,NBNN,IB,MELE,LUNI1,IPTR1, . 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 * CALL DEFVAL(NUMAT1,NBPTEL,NDEF, . 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 * CALL DEFMAT(NMATT,NSTRS,MFR,MELE,INPLAS, . 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 * * >>>>>>>>>> fin du traitement du materiau * * on recupere les caracteristiques geometriques * CALL DEFCAR(NCARR,IB,IGAU,MFR,MELE,IVACAR, . WRK1) * * 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. 17 .OR. 2 (INPLAS .GE. 19 .AND. INPLAS .LE. 25) .OR. 4 INPLAS .EQ. 61 .OR. 4 INPLAS .EQ. 63 .OR. 1 INPLAS .EQ. 53 .OR. INPLAS .EQ. 102 .OR. 8 INPLAS .EQ. 44 .OR. INPLAS .EQ. 76 .OR. 9 INPLAS .EQ. 45 .OR. INPLAS .EQ. 77 .OR. 9 INPLAS .EQ. 84 .OR. INPLAS .EQ. 85 .OR. 9 INPLAS .EQ. 86 .OR. INPLAS .EQ. 70 ) THEN * IF (INPLAS.EQ.44.AND.VAR0(NVARI).EQ.0.0) THEN VAR0(NVARI)=XMAT(20) ENDIF IF (INPLAS.EQ.45.AND.VAR0(NVARI).EQ.0.0) THEN VAR0(NVARI-2)=XMAT(20) VAR0(NVARI-1)=XMAT(21) VAR0(NVARI)=XMAT(27) ENDIF * SEGINI WRK7 SEGINI WRK9 CALL CONSTI(WRK0,WR00,WRK1,WRK5,WRK7,WRK8,WRK9,WTRAV, 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 c write(6,*) istep IF (MFR1.EQ.17.AND.INPLAS.EQ.19) THEN IF (KERREU1.NE.0.AND.NSSINC.EQ.1) THEN CALL ERREUR(KERREU1) ENDIF ENDIF SEGSUP WRK7 SEGSUP WRK9 DTOPTI = MIN(DTOPTI,DTT) NINCMA = MAX(NINCMA,NSSINC) NCOMP = NCOMP + 1 TSOM = TSOM + DTT NSOM = NSOM + NSSINC NINV = NINV + INV TCAR = TCAR + DTT* DTT IF(KERRE.NE.0.AND.KERRE.NE.99) THEN KERR1=1 END IF c ELSE KERRE = 99 ENDIF * * Erreurs * - problèmes de convergence * CALL DEFER1(JNOID,KERR1,KERRE,LOGSUC) * * - autres problèmes * 1990 CALL DEFER2(INPLAS,MFR,MELE,IB,IGAU,KERR1,KERRE) IF (KERRE.NE.0) GOTO 99 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 CALL DEFSIG(MFR,NDEF, . INPLAS,IND,WRK1,WRK5,WTRAV, . IVASTF,IVARIF,IVADEP,COB,XMOB,IB,IGAU, . CMATE,MATE,MELE,KERRER) IF (KERRER.NE.0) GOTO 1000 c c fin de la boucle sur les points de gauss c 1100 continue c c special poutres et tuyaux sauf timoschenko c CALL DEFPOU(MFR,MELE,MELEME,IB,WRK4,IVASTF) 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' * CALL DEFFIN(INPLAS,TSOM,NSOM,NCOMP,NINV,NINCMA, . TCAR,DTOPTI,IPOTAB,KERRE) * 99 CONTINUE SEGSUP WRK0,WR00,WRK1,WRK2,WRK3,WRK5,WTRAV IF (LOGVIS) SEGSUP WRK8 IF (WRK4.NE.0) SEGSUP WRK4 IF (WRK22.NE.0) THEN SEGDES MINTE2 SEGSUP WRK22 ENDIF IF (ITHHER.EQ.1.OR.ITHHER.EQ.2) THEN SEGDES MELVA3,MELVA4,MELVA5 SEGDES MCHAM3,MCHAM4,MCHAM5 ENDIF RETURN END