C FCOUL2 SOURCE CB215821 23/01/25 21:15:13 11573 SUBROUTINE FCOUL2(DEPSI,INFIBR,MELE,IPMAIL,MINTE,NBPTEL,IVASTR, 1 IVARI,IVAMAT,IVACAR,NSTRS,NVARI,NMATT,NCARR,TIME0,TIMEF, 2 SIGMA,IVASTF,IVARIF,EPSUP,EPINF,DAMAG,NSTRS2) *********************************************************************** * ECOULEMENT INELASTIQUE POUR LES MODELES A FIBRES * travail sur chaque les element de chaque ss-zone du modele * de section ********************************************************************** * Pierre Pegon (ISPRA) Juillet/Aout 1993 *********************************************************************** * ENTREES : * * DEPSI(6) INCREMENT DE DEFORMATION POUR LA FIBRE CENTRALE * INFIBR = NUMERO DE MATERIAU INELASTIQUE * MELE = NUMERO ELEMENT FINI * IPMAIL = POINTEUR DU MAILLAGE * NBPTEL = NOMBRE DE POINTS PAR ELEMENT * IVASTR = POINTEUR SUR UN SEGMENT MPTVAL DE CONTRAINTES * IVAMAT = POINTEUR SUR UN SEGMENT MPTVAL DE MATERIAU * IVACAR = POINTEUR SUR UN SEGMENT MPTVAL DE CARACT. GEOMETRIQUES * 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 * TIME0 = INSTANT INITIAL * TIMEF = INSTANT FINAL * * SORTIES : * SIGMA(6) EFFORT SUR LA FIBRE MOYENNE * IVASTF = POINTEUR SUR UN SEGMENT MPTVAL DE CONTRAINTES * IVARIF = POINTEUR SUR UN SEGMENT MPTVAL DE VARIABLES INTERNES * ************************************************************************ 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 * SEGMENT MPTVAL INTEGER IPOS(NS) ,NSOF(NS) INTEGER IVAL(NCOSOU) CHARACTER*16 TYVAL(NCOSOU) ENDSEGMENT * SEGMENT WWRK0 REAL*8 XMAT(NCXMAT),XCAR(NCXCAR) ENDSEGMENT * SEGMENT WWRK1 REAL*8 SIG0(NSTRS),SIGF(NSTRS) REAL*8 VAR0(NVARI),VARF(NVARI) ENDSEGMENT * SEGMENT WWRK2 REAL*8 XE(3,NBBB),SHP(6,NBBB) ENDSEGMENT * SEGMENT WRK2 REAL*8 TRAC(LTRAC) ENDSEGMENT * DIMENSION DEPSI(NSTRS2),SIGMA(NSTRS2) DIMENSION DEPS(3),DEPSB(3),SIG0B(3),SIGFB(3) * C+PP IST_DES=0 IST_TOT=0 * write(6,*) ' infibr ',infibr C+PP MFR =NUMMFR(MELE) MELEME=IPMAIL NBNN=NUM(/1) NBELEM=NUM(/2) NDEF=NSTRS * * SEGMENT D'INTEGRATION * C* SEGACT,MINTE <- ACTIF EN E/S * * INITIALISATION DES SEGMENTS DE TRAVAIL * NCXMAT=NMATT + 1 NCXCAR=NCARR NBBB=NBNN SEGINI WWRK0,WWRK1,WWRK2 LTRAC=260 SEGINI WRK2 * * BOUCLE SUR LES ELEMENTS * SEGACT,MCOORD DO 1000 IB=1,NBELEM * * ON CHERCHE LES COORDONNEES DES NOEUDS DE L ELEMENT IB * CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE) * * BOUCLE SUR LES POINTS DE GAUSS * DO 1100 IGAU=1,NBPTEL * * ON CHERCHE LA POSITION DU POINT DE LA SECTION (X->Y) (Y->Z) * YY=0.D0 ZZ=0.D0 DO IE1=1,NBNN CGAUSS=SHPTOT(1,IE1,IGAU) YY=YY+XE(1,IE1)*CGAUSS ZZ=ZZ+XE(2,IE1)*CGAUSS END DO * * ON REMPLIT LES SHP ET ON CALCUL LE JACOBIEN * DO IE2=1,NBNN DO IE1=1,6 SHP(IE1,IE2)=SHPTOT(IE1,IE2,IGAU) END DO END DO C PPf CALL JACOBI(XE,SHP,2,NBNN,DJAC) * * ON EN DEDUIT L'INCREMENT DE DEFORMATION * IF (NSTRS2.EQ.3) THEN DEPS(1)=DEPSI(1)-YY*DEPSI(3) DEPS(2)=DEPSI(2) DEPS(3)=0.D0 ELSE DEPS(1)=DEPSI(1)+ZZ*DEPSI(5)-YY*DEPSI(6) DEPS(2)=DEPSI(2)-ZZ*DEPSI(4) DEPS(3)=DEPSI(3)+YY*DEPSI(4) ENDIF * * ON RECUPERE LES CONTRAINTES INITIALES * MPTVAL=IVASTR DO IC=1,NSTRS MELVAL=IVAL(IC) IBMN=MIN(IB,VELCHE(/2)) IGMN=MIN(IGAU,VELCHE(/1)) SIG0(IC)=VELCHE(IGMN,IBMN) END DO * * ON RECUPERE LES VARIABLES INTERNES * MPTVAL=IVARI DO IC=1,NVARI MELVAL=IVAL(IC) IBMN=MIN(IB,VELCHE(/2)) IGMN=MIN(IGAU,VELCHE(/1)) VAR0(IC)=VELCHE(IGMN,IBMN) END DO * * ON RECUPERE LES CONSTANTES DU MATERIAU * MPTVAL=IVAMAT DO IC=1,NMATT MELVAL=IVAL(IC) IF(IC.LT.3)THEN IIC=IC ELSEIF(IC.LT.(NMATT-1))THEN IIC=IC+2 ELSEIF(IC.LE.(NMATT))THEN IIC=4+IC-NMATT ELSE ENDIF C IF(MELVAL.NE.0)THEN IF(TYVAL(IC)(1:8).NE.'POINTEUR')THEN IBMN=MIN(IB,VELCHE(/2)) IGMN=MIN(IGAU,VELCHE(/1)) XMAT(IIC)=VELCHE(IGMN,IBMN) ELSE IBMN=MIN(IB,IELCHE(/2)) IGMN=MIN(IGAU,IELCHE(/1)) XMAT(IIC)=IELCHE(IGMN,IBMN) ENDIF ELSE XMAT(IIC)=0.D0 IF(TYVAL(IC)(1:8).EQ.'POINTEUR') THEN XMAT(IIC)=0.D0 END IF ENDIF END DO * * ON RECUPERE LES CARACTERISTIQUES GEOMETRIQUES * MPTVAL=IVACAR DO IC=1,NCARR MELVAL=IVAL(IC) * si c'est une caracteristique facultative non remplie melval vaut 0 if (melval.ne.0) then IBMN=MIN(IB,VELCHE(/2)) IGMN=MIN(IGAU,VELCHE(/1)) XCAR(IC)=VELCHE(IGMN,IBMN) endif END DO * *--------------------------------------------------------------------- * * ECOULEMENT SELON LES MODELES * *--------------------------------------------------------------------- * IF(INFIBR.EQ.0)THEN C C MODELE ELASTIQUE LINEAIRE (EXEMPLE) C CALL FIBELA(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C ELSEIF(INFIBR.EQ.1)THEN C C MODELE BETON_UNI C C IF (XMAT(14).LT.0.D0) THEN * write(6,*) ' fcoul2 appel fibeto' CALL FIBETO(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) ELSE * write(6,*) ' fcoul2 appel fibet2' CALL FIBET2(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) ENDIF C IF (EPSUP .LT. VARF(1)) EPSUP=VARF(1) IF (EPINF .GT. VARF(1)) EPINF=VARF(1) C ELSEIF(INFIBR.EQ.2)THEN C C MODELE ACIER_UNI C CALL FIBSTE(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C+PP IST_TOT=IST_TOT+1 * write (6,*) 'fcoul2 apres fibste ',varf(1) IF(INT(VARF(1)).EQ.1)IST_DES=IST_DES+1 C+PP C ELSEIF(INFIBR.EQ.10)THEN C C MODELE ACIER_ANCRAGE C CALL FIBSTA(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C+PP IST_TOT=IST_TOT+1 C write (6,*) 'fcoul2 apres fibsta ',varf(1) IF(INT(VARF(1)).EQ.1)IST_DES=IST_DES+1 C+PP C ELSEIF(INFIBR.EQ.3)THEN C C MODELE MAZARS_FIB C CALL FIBMAZ(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C ELSEIF(INFIBR.EQ.11)THEN C C MODELE CLB_UNI C CALL LABORD(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C C ELSEIF(INFIBR.EQ.4)THEN C C MODELE FRAGILE_UNI C CALL FIBFRA(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C ELSEIF(INFIBR.EQ.5)THEN C C MODELE BETON_BAEL C CALL FIBAEL(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C IF (EPSUP .LT. VARF(2)) EPSUP=VARF(2) IF (EPINF .GT. VARF(2)) EPINF=VARF(2) C ELSEIF(INFIBR.EQ.6)THEN C C MODELE PARFAIT_UNI C CALL FIBPAR(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C ELSEIF(INFIBR.EQ.9)THEN C C MODELE PARFAIT_ANCRAGE C CALL FIBPAA(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C ELSEIF(INFIBR.EQ.12)THEN C C MODELE INTIMP (CHOIX SELON LE TYPE DE CALAGE) C IF (XMAT(18).EQ.0.D0) THEN CALL INTIMP(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) ELSEIF (XMAT(18).EQ.1.D0) THEN CALL INTFIC(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) ELSEIF (XMAT(18).EQ.2.D0) THEN CALL OUGLFI(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) ENDIF C +BR ELSEIF(INFIBR.EQ.13)THEN C C MODELE RICBET_UNI C IF (XMAT(16).EQ.1) THEN CALL RICBETF1(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) ELSEIF (XMAT(16).EQ.2) THEN CALL RICBETF2(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) ENDIF C -BR C +RP ELSEIF(INFIBR.EQ.14)THEN C C OUGLOVA C CALL OUGLOF(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) C -RP ELSEIF(INFIBR.EQ.8)THEN C C MODELE CISAIL_NL C IPOS1=1 CALL COTRAE(WWRK0,WRK2,12,IPOS1,0, NPOINT,KERRE) NTRAP=NPOINT/2 IPOS2=IPOS1+NPOINT CALL COTRAE(WWRK0,WRK2,13,IPOS2,0, NPOINT,KERRE) NTRAN=NPOINT/2 IF(KERRE.EQ.0) THEN * write(6,*) ' fcoul2 appel fibeta' CALL FIBETA(XMAT,XCAR,SIG0,VAR0,SIGF,VARF,DEPS, . WRK2,NTRAP,NTRAN,KERRE) END IF C ELSEIF(INFIBR.EQ.7)THEN C C MODELE STRUT_UNI C C XKEPM = -1. XEULT = XMAT(4+20) C IF (NSTRS2.EQ.3) THEN VARF(30)=VAR0(30)+DEPSI(3) ELSE VARF(30)=VAR0(30)+DEPSI(6) ENDIF C-------------------------------------------------------- C EPSUP - Maximum compression strain C EPINF - Maximum tensile strain C-------------------------------------------------------- C EPS11=(EPSUP+(ABS(XKEPM)-1.D0)*EPINF)/ABS(XKEPM) C EPS22=((ABS(XKEPM)-1.D0)*EPSUP+EPINF)/ABS(XKEPM) C C IF (VARF(30) .GE. 0.D0) THEN C VARF(28)=EPS22 C VARF(29)=EPS11 C ELSE C VARF(28)=EPS11 C VARF(29)=EPS22 C ENDIF C C EPS11 = VAR0(28)+DEPSI(1) C EPS11 = 0.5D0 * (EPSUP +EPINF) EPS22 = EPS11 VARF(28)= EPS11 VARF(29)= EPS11 C VARF(34)=VAR0(34) VARF(35)=VAR0(35) *-------------------------------------------------------- * CHECK IF THE SHEAR DEFORMATION CHANGED SIGN *-------------------------------------------------------- SHEXY=VAR0(25)+DEPS(2) C IF (SHEXY .GE. 0.0D0) THEN KFAC1=28 KFAC2=34 ELSE KFAC1=29 KFAC2=35 ENDIF * *-------------------------------------------------------- * CORRECT THE MAXIMUM ALLOWED AXIAL DEFORMATION *-------------------------------------------------------- * FREEZE THE AVERAGE AXIAL STRAIN WHEN * SHEAR STRAIN CHANGES SIGN *-------------------------------------------------------- IF (((SHEXY*VAR0(25)) .LE. 0.0D0) .AND. . (VAR0(25) .NE. 0.0D0)) THEN IF ((VARF(KFAC1) .GT. 0.0D0) .AND. . (VAR0(KFAC1) .GT. 0.0D0)) THEN FACPR=VAR0(25)/DEPS(2) VARF(KFAC1)=FACPR*(VAR0(KFAC1)-VARF(KFAC1))+ . VAR0(KFAC1) VARF(KFAC2)=VARF(KFAC1) ENDIF ENDIF *-------------------------------------------------------- * CHECK IF THE AXIAL DEFORMATION IS BELOW THE LIMIT *-------------------------------------------------------- IF (VARF(28) .LT. VAR0(28)) THEN IF (VARF(28) .LT. VAR0(34)) THEN IF (VARF(34) .GT. 0.0D0) THEN VARF(34)=VAR0(34) VARF(28)=VAR0(34) ENDIF ENDIF ENDIF IF (VARF(29) .LT. VAR0(29)) THEN IF (VARF(29) .LT. VAR0(35)) THEN IF (VARF(35) .GT. 0.0D0) THEN VARF(35)=VAR0(35) VARF(29)=VAR0(35) ENDIF ENDIF ENDIF C C IF (XEULT.GE. 0.D0) THEN *-------------------------------------------------------- * FREEZE THE AVERAGE AXIAL STRAIN WHEN * THERE ARE CRACKS OPENED *-------------------------------------------------------- C VARF(34)=VAR0(34) C VARF(35)=VAR0(35) C C IF (ETIQE .EQ. 0.D0) THEN C IF (VARF(34) .LT. VAR0(28)) VARF(34)=VAR0(28) C IF (VARF(35) .LT. VAR0(29)) VARF(35)=VAR0(29) C ENDIF C C IF ((EPS22 .LT. VARF(34)) .AND. C . (VARF(34) .GT. 0.D0)) THEN C EPS22=VARF(34) C EPS11=VARF(35) C ENDIF C C VARF(28)=EPS22 C VARF(29)=EPS11 C ENDIF *-------------------------------------------------------- * DAMAG - Maximum compression strain / Ultimate strain *-------------------------------------------------------- DAMAG = 0.D0 IF (ABS(XEULT) .LE. 1.0D0) THEN IF (XEULT.GE.0.D0) THEN C C DAMAG - Position of the neutral axis C IF ((EPINF*EPSUP).LT.0.0D0) THEN DAMGG = EPSUP/( EPSUP - EPINF) ELSE DAMGG =0.D0 ENDIF ELSE C C DAMAG - Maximum compression strain / Ultimate strain C DAMGG=-1.0D0*EPINF/ABS(XEULT) c ENDIF C VARF(32)=VAR0(32) VARF(33)=VAR0(33) C IF (SHEXY .GT. 0.0D0) THEN IF (DAMGG .GE. VARF(32)) VARF(32)=DAMGG ELSE IF (DAMGG .GE. VARF(33)) VARF(33)=DAMGG ENDIF C ELSE *-------------------------------------------------------- * DO NOT CONSIDER DAMAGE IN THE STRUTS *-------------------------------------------------------- VARF(32)=0.0D0 VARF(33)=0.0D0 ENDIF C 2001 CONTINUE IF (NSTRS2.EQ.3) THEN DEPSB(1) = DEPS(1) DEPSB(2) = DEPS(2) DEPSB(3) = 0.D0 SIG0B(1) = SIG0(1) SIG0B(2) = SIG0(2) SIG0B(3) = 0.D0 CALL FIBSTR(XMAT,DEPSB,SIG0B,VAR0,SIGFB,VARF) SIGF(1) = SIGFB(1) SIGF(2) = SIGFB(2) ELSE CALL FIBSTR(XMAT,DEPS,SIG0,VAR0,SIGF,VARF) ENDIF ELSEIF((INFIBR.EQ.15).OR.(INFIBR.EQ.16).OR.(INFIBR.EQ.17).OR. & (INFIBR.EQ.18).OR.(INFIBR.EQ.19)) THEN C C POUR LES MODELES DE FLUAGE C CALL FLUFIB(INFIBR,XMAT,DEPS,SIG0,VAR0,SIGF,VARF,TIME0, & TIMEF,NCXMAT,NVARI) C ENDIF C+PPf C C TRAITEMENT PARTICULIER DES ELEMENTS SEGS(166) ET POJS(167) C C IF(MELE.EQ.167)THEN C+DC DJAC=XCAR(NCARR) IF (NSTRS2.EQ.3) THEN DJAC=XCAR(2) ELSE DJAC=XCAR(3) ENDIF C ELSEIF(MELE.EQ.166)THEN CALL JACOBI(XE,SHP,1,NBNN,DJAC) C+DC DJAC=DJAC*XCAR(NCARR) IF (NSTRS2.EQ.3) THEN DJAC=DJAC*XCAR(2) ELSE DJAC=DJAC*XCAR(3) ENDIF ELSE CALL JACOBI(XE,SHP,2,NBNN,DJAC) ENDIF C+PPf C C CONTRIBUTION A LA CONTRAINTE DE LA SECTION C PGAUSS=POIGAU(IGAU)*ABS(DJAC) IF (NSTRS2.EQ.3) THEN SIGMA(1)=SIGMA(1)+SIGF(1)*PGAUSS SIGMA(2)=SIGMA(2)+XCAR(1)*SIGF(2)*PGAUSS SIGMA(3)=SIGMA(3)-YY*SIGF(1)*PGAUSS ELSE SIGMA(1)=SIGMA(1)+SIGF(1)*PGAUSS SIGMA(2)=SIGMA(2)+XCAR(1)*SIGF(2)*PGAUSS SIGMA(3)=SIGMA(3)+XCAR(2)*SIGF(3)*PGAUSS SIGMA(4)=SIGMA(4)+ $ (-ZZ*XCAR(1)*SIGF(2)+YY*XCAR(2)*SIGF(3))*PGAUSS SIGMA(5)=SIGMA(5)+ZZ*SIGF(1)*PGAUSS SIGMA(6)=SIGMA(6)-YY*SIGF(1)*PGAUSS ENDIF C C REMPLISSAGE DU SEGMENT CONTENANT LES CONTRAINTES A LA FIN C MPTVAL=IVASTF DO 1116 IC=1,NSTRS MELVAL=IVAL(IC) VELCHE(IGAU,IB)=SIGF(IC) 1116 CONTINUE C C ET LES VARIABLES INTERNES FINALES C MPTVAL=IVARIF DO 1117 IC=1,NVARI MELVAL=IVAL(IC) VELCHE(IGAU,IB)=VARF(IC) 1117 CONTINUE C C FIN DE LA BOUCLE SUR LES POINTS DE GAUSS C 1100 CONTINUE C C FIN DE LA BOUCLE SUR LES ELEMENTS C 1000 CONTINUE C+PP IF(IST_DES.NE.0.and.iimpi.eq.19932)THEN WRITE(IOIMP,*)'FCOUL2:',IST_DES,' steel fibres out of ', > IST_TOT,' are destroyed on the current zone' ENDIF C+PP * SEGSUP WRK2 SEGSUP WWRK0,WWRK1,WWRK2 * RETURN END