1 : C MOCU3D SOURCE CHAT 11/03/16 21:28:02 6902 2 : C MOCU3D SOURCE DC 98/06/03 21:45:21 3202 3 : C MOCU3D SOURCE GOUNAND 97/10/13 22:21:06 2872 4 : SUBROUTINE MOCU3D(IPCURY,IPCURZ,IPAXIA,IPMODL,IPCAR, 5 : . IPEPAX,IPMOMY,IPMOMZ,TOL,mtable,zveri) 6 : ********************************************************************** 7 : * 8 : * OPERATEUR MOCU3D (MOment/CUrvature) 9 : * 10 : * IPCURY : Histoire en courbure 11 : * IPCURZ : Histoire en courbure 12 : * IPAXIA : Histoire en effort normal 13 : * IPMODL : Modele 14 : * IPCAR : Materiau 15 : * IPEPAX : Histoire allongement axial 16 : * IPMOMY : Histoire en moment 17 : * IPMOMZ : Histoire en moment 18 : * TOL : Tolerance 19 : ********************************************************************** 20 : * D'APRES MOCUR DC 98 21 : ********************************************************************** 22 : IMPLICIT INTEGER(I-N) 23 : IMPLICIT REAL*8(A-H,O-Z) 24 : C 25 : DIMENSION DEPSI(6),SIGMA(6),CRIGI(12),CMASS(12) 26 : logical bool,boolre,zveri 27 : character*8 char1,charre 28 : C 29 : -INC CCOPTIO 30 : -INC SMLREEL 31 : -INC CCGEOME 32 : -INC CCREEL 33 : -INC SMMODEL 34 : -INC SMCHAML 35 : -INC SMTABLE 36 : POINTEUR IPCURY.MLREEL,IPCURZ.MLREEL,IPAXIA.MLREEL 37 : POINTEUR IPMOMY.MLREEL,IPMOMZ.MLREEL,IPEPAX.MLREEL 38 : C 39 : NMAXIT=50 40 : NSTRS2 = 6 41 : 42 : C 43 : * 44 : ICOUL1=IDCOUL 45 : C 46 : C------- NSTEPS -- 47 : C 48 : SEGACT,IPCURY,IPCURZ,IPAXIA 49 : NSTEPS=IPCURY.PROG(/1) 50 : IF(IPCURZ.PROG(/1).NE.NSTEPS)THEN 51 : SEGDES,IPCURY,IPCURZ,IPAXIA 52 : RETURN 53 : ENDIF 54 : IF(IPAXIA.PROG(/1).NE.NSTEPS)THEN 55 : SEGDES,IPCURY,IPCURZ,IPAXIA 56 : RETURN 57 : ENDIF 58 : c 59 : if (zveri) then 60 : m=3 61 : segini mtable 62 : itab=mtable 63 : m=nsteps+1 64 : segini mtab1 65 : itab1=mtab1 66 : c 66 : call ecctab (itab,'MOT ',ii,xva,'CONTRAINTES',bool,iob, 67 : $ 'TABLE ',ire,xvare,char1,boolre,itab1) 68 : segini mtab2 68 : call ecctab (mtable,'MOT ',ii,xva,'VARIABLES_INTERNES', 69 : $ bool,iob, 70 : $ 'TABLE ',ii,xva,char1,bool,mtab2) 71 : segini mtab3 71 : call ecctab (mtable,'MOT ',ii,xva,'DEFORMATIONS', 72 : $ bool,iob, 73 : $ 'TABLE ',ii,xva,char1,bool,mtab3) 74 : segdes mtable 75 : endif 76 : C 77 : C-------- INITIALIZING --------------- 78 : C 78 : CALL ZEROP(IPMODL,'CONTRAIN',IPSIG0) 78 : CALL ZEROP(IPMODL,'VARINTER',IPVAR0) 79 : DO IE1=1,6 80 : SIGMA(IE1)=0.D0 81 : END DO 82 : SZERO=0.0D0 83 : if (zveri) then 83 : call ecctab(mtab1,'ENTIER ',0,xva,char1,bool,iob, 84 : $ 'MCHAML ',ii,xva,char1,bool,ipsig0) 84 : call ecctab(mtab2,'ENTIER ',0,xva,char1,bool,iob, 85 : $ 'MCHAML ',ii,xva,char1,bool,ipvar0) 85 : call ecctab(mtab3,'ENTIER ',0,xva,char1,bool,iob, 86 : $ 'FLOTTANT',ii,szero,char1,bool,iob) 87 : endif 88 : JG=NSTEPS 89 : SEGINI,IPMOMY,IPMOMZ,IPEPAX 90 : C 91 : C INITIAL STIFNESS 92 : C 92 : CALL FRIGIE(IPMODL,IPCAR,CRIGI,CMASS) 93 : XKFXEX=CRIGI(1) 94 : XKFXCY=CRIGI(3) 95 : XKFXCZ=-CRIGI(2) 96 : C 97 : C NULL DEFORMATION INCREMENT 98 : C 99 : DEPSI(2)=0.D0 100 : DEPSI(3)=0.D0 101 : DEPSI(4)=0.D0 102 : C 103 : YAUX0=0.D0 104 : ZAUX0=0.D0 105 : ZAXI0=0.D0 106 : DEFOA=0.D0 107 : C 108 : C-------- LOOP --------------- 109 : C 110 : DO 2000 I=1,NSTEPS 111 : C 112 : YAUX1=IPCURY.PROG(I) 113 : ZAUX1=IPCURZ.PROG(I) 114 : ZAXI1=IPAXIA.PROG(I) 115 : 116 : DYCURV=YAUX1-YAUX0 117 : DZCURV=ZAUX1-ZAUX0 118 : DFAXIA=ZAXI1-ZAXI0 119 : C 120 : SIGPRE=ZAXI0 121 : XKFXEX=CRIGI(1) 122 : DSITER=(DFAXIA-XKFXCY*DYCURV-XKFXCZ*DZCURV)/XKFXEX 123 : SZERO =DSITER 124 : C 125 : DEPSI(5)=DYCURV 126 : DEPSI(6)=DZCURV 127 : C 128 : DO 5000 NITER=1,NMAXIT 129 : DEPSI(1)=SZERO 130 : C 130 : CALL FCOUL1(DEPSI,IPMODL,IPSIG0,IPVAR0,IPCAR, 131 : $ SIGMA,IPSIGF,IPVARF,IRETO,NSTRS2) 132 : IF(IRETO.EQ.0) THEN 133 : SEGSUP,IPMOMY,IPMOMZ,IPEPAX 134 : RETURN 135 : ENDIF 136 : C 137 : DSITET=(ZAXI1-SIGMA(1))/XKFXEX 138 : C 139 : IF(IIMPI.EQ.1789)THEN 140 : write(IOIMP,999)I,NITER,XKFXEX,DSITET 141 : 999 format(1x,'I=',I4,' NITER=',I2,' XKFXEX=', 142 : $ 1pd12.5,' DSITET=',1pd12.5) 143 : ENDIF 144 : C 145 : IF(ABS(DSITET).LT.TOL)GOTO 5010 146 : C 147 : C on tente une approximation quasi-newton apres la prediction ... 148 : C 149 : IF(NITER.GT.1) XKFXEX=(SIGMA(1)-SIGPRE)/DSITER 150 : SIGPRE=SIGMA(1) 151 : DSITER=(ZAXI1-SIGMA(1))/XKFXEX 152 : SZERO=SZERO+DSITER 153 : C 153 : CALL DTCHAM(IPSIGF) 153 : CALL DTCHAM(IPVARF) 154 : C 155 : 5000 CONTINUE 156 : C 157 : INTERR(1)=NMAXIT 157 : CALL ERREUR(638) 157 : CALL DTCHAM(IPSIG0) 157 : CALL DTCHAM(IPVAR0) 158 : SEGSUP,IPMOMY,IPMOMZ,IPEPAX 159 : SEGDES,IPCURY,IPCURZ,IPAXIA 160 : RETURN 161 : C 162 : 5010 ZAXI0=SIGMA(1) 163 : ZAUX0=ZAUX1 164 : YAUX0=YAUX1 165 : DEFOA=DEFOA+SZERO 166 : C 167 : if (.not.zveri) then 167 : CALL DTCHAM(IPSIG0) 167 : CALL DTCHAM(IPVAR0) 168 : endif 169 : IPSIG0=IPSIGF 170 : IPVAR0=IPVARF 171 : c 172 : if (zveri) then 172 : call ecctab(mtab1,'ENTIER',I,xva,char1,bool,iob, 173 : $ 'MCHAML',ii,xva,char1,bool,ipsig0) 173 : call ecctab(mtab2,'ENTIER',I,xva,char1,bool,iob, 174 : $ 'MCHAML',ii,xva,char1,bool,ipvar0) 174 : call ecctab(mtab3,'ENTIER ',I,xva,char1,bool,iob, 175 : $ 'FLOTTANT ',ii,DEFOA,char1,bool,iob) 176 : endif 177 : C 178 : 179 : 180 : IPMOMY.PROG(I)=SIGMA(5) 181 : IPMOMZ.PROG(I)=SIGMA(6) 182 : IPEPAX.PROG(I)=DEFOA 183 : C 184 : C 185 : 2000 CONTINUE 186 : C 187 : if (zveri) then 188 : segdes,mtab1,mtab2,mtab3 189 : else 189 : CALL DTCHAM(IPSIG0) 189 : CALL DTCHAM(IPVAR0) 190 : endif 191 : SEGDES,IPCURY,IPCURZ,IPAXIA 192 : C 193 : SEGDES,IPMOMY,IPMOMZ,IPEPAX 194 : C 195 : RETURN 196 : END 197 : 198 : 199 : 200 : 201 : 202 : 203 : 204 : 205 : 206 : 207 : 208 : 209 : 210 : 211 : 212 : 213 : 214 : 215 :
© Cast3M 2003 - All rights reserved.
Disclaimer