1 : C SBSTEP SOURCE CHAT 05/01/13 03:12:07 5004 2 : CCC 3 : C ********************************************************************** 4 : CCC 5 : SUBROUTINE SUBSTEP 6 : . (sigini,varini,dsigtr,sigfin,varfin, 7 : . ddefpl,Kmat,ndims,ndimv,ndimk, 8 : . xmat,kerre,precis,nitmax,nescri,ues, 9 : . nmodel,nnumer,deltax,nsteps, 10 : . isub,ntotiter,iincre) 11 : c 12 : c IN: 13 : c sigini Initial stresses 14 : c varini Initial internal variables 15 : c dsigtr Increment of trial stresses 16 : c ndims Stress dimension 17 : c ndimv Internal variable dimension 18 : c ndimk Tangent stiffness matrix dimension 19 : c xmat Material parameters 20 : c precis Tolerance for local problem iterations 21 : c nitmax Maximun iterations for local problem 22 : c nescri Write intermediate results (0:no, 1: yes) 23 : c ues Output unit 24 : c nmodel Plasticity model 25 : c nnumer Generalized flow vector computation 26 : c 0:analytical 27 : c 1:1st forward differences 28 : c 2:1st centered differences 29 : c 3:1st complex approximation 30 : c deltax Relative stepsize value 31 : c nsteps Number of substeps 32 : c OUT: 33 : c sigfin Final stresses 34 : c varfin Final internal variables 35 : c ddefpl Increment of plastic strains 36 : c Kmat Consistent tangent stiffness matrix 37 : c kerre Error control flag 38 : IMPLICIT INTEGER(I-N) 39 : 40 : integer ndims,ndimv,ndimk,kerre,nitmax,nescri,iterlocal, 41 : . nnumer,ues,nmodel,nsteps,isub,ntotiter,iincre 42 : real*8 sigini(*),varini(*),ddefpl(*), 43 : . dsigtr(*),sigfin(*),varfin(*), 44 : . Kmat(*),xmat(*),precis,deltax 45 : integer ips,i,k,ivoid,ndimx,flagini,iaux,ic,jf 46 : real*8 sig1(6),sig2(6),ddsigtr(6),dddefpl(6), 47 : . var1(4),var2(4),Kmat0(64),Kmat1(64),Kmat2(64), 48 : . alphai 49 : 50 : do i=1,2 51 : var1(i)=0.D0 52 : var2(i)=0.D0 53 : enddo 54 : do i=1,6 55 : sig1(i)=0.D0 56 : sig2(i)=0.D0 57 : ddsigtr(i)=0.D0 58 : dddefpl(i)=0.D0 59 : enddo 60 : do i=1,64 61 : Kmat0(i)=0.D0 62 : Kmat1(i)=0.D0 63 : Kmat2(i)=0.D0 64 : enddo 65 : c nescri=1 66 : c variable de ecoj2 que no se usa 67 : ivoid=0 68 : c dimension de sigma ampliado 68 : call fijarndimx(nmodel,xmat,ndims,ndimx) 69 : c posicion en el step (valor de i inicial) 70 : ips =0 71 : c contador de substeps 72 : isub =0 73 : c contador de iteraciones 74 : ntotiter=0 75 : c flag de primer substep para la mtc 76 : flagini=1 77 : c poner sigma y varint iniciales 78 : do i=1,ndims 79 : sig1(i)=sigini(i) 80 : enddo 81 : do i=1,ndimv 82 : var1(i)=varini(i) 83 : enddo 84 : do while (ips.lt.nsteps) 85 : 100 continue 86 : c proporcion de carga para el substep 87 : alphai=float(iincre)/float(nsteps) 88 : c dividir incremento de sigma 89 : do k=1,ndims 90 : ddsigtr(k)=dsigtr(k)*alphai 91 : enddo 92 : c integra un substep 93 : if (nmodel.eq.1) then 93 : call eco_j2(sig1,var1,ddsigtr,sig2,var2, 94 : . dddefpl,ivoid,ndims,xmat,kerre, 95 : . precis,nitmax,nescri,ues,iterlocal) 96 : else if (nmodel.eq.2) then 96 : call eco_rhmc(sig1,var1,ddsigtr,sig2,var2, 97 : . dddefpl,ivoid,ndims,xmat,kerre, 98 : . precis,nitmax,nescri,ues,iterlocal) 99 : else 99 : call eco_MRSMAC(sig1,var1,ddsigtr,sig2,var2, 100 : . dddefpl,ivoid,ndims,xmat,kerre, 101 : . precis,nitmax,nescri,ues,nnumer,deltax, 102 : . iterlocal) 103 : endif 104 : c criterio de aceptacion 105 : if (kerre.eq.1) then 106 : if (iincre.eq.1) then 107 : write(*,*) ' Numero de substeps max. insuficiente' 108 : return 109 : else 110 : iincre=int(iincre/2) 111 : iaux=nsteps-(ips+iincre) 112 : if (iaux.lt.0) iincre=nsteps-ips 113 : c repite el step dividido por la mitad 114 : kerre=0 115 : goto 100 116 : endif 117 : endif 118 : c añadir defor plastica del step 119 : do i=1,ndims 120 : ddefpl(i)=dddefpl(i) 121 : enddo 122 : do i=1,6 123 : dddefpl(i)=0.D0 124 : enddo 125 : c poner como inicial las sigma y varint finales del step 126 : do i=1,ndims 127 : sig1(i)=sig2(i) 128 : enddo 129 : do i=1,ndimv 130 : var1(i)=var2(i) 131 : enddo 132 : c actualizar la i de posicion en el step 133 : ips=ips+iincre 134 : c añadir uno al contador de steps 135 : isub=isub+1 136 : c añadir el numero de iteraciones al total 137 : ntotiter=ntotiter+(iterlocal) 138 : c calcula la nx.nx de la inversa del jacobiano (sin E^-1 en la diagonal) 139 : if (nmodel.eq.1) then 139 : call jac_J2(Kmat1,ndimx,sig2,ndims,var2(1),var2(2),xmat, 140 : . nescri,ues,kerre) 141 : else if (nmodel.eq.2) then 141 : call jac_rhmc(Kmat1,ndimx,sig2,ndims,var2(2),xmat, 142 : . nescri,ues,kerre) 143 : else 143 : call jac_mrsmac(Kmat1,ndimx,sig2,ndims, 144 : . var2(1),var2(2),var2(3),var2(4), 145 : . xmat,nescri,ues,nnumer,deltax,kerre) 146 : endif 147 : if (flagini.eq.1) then 148 : c se guardan ns columnas de Kmat1 en Kmat0 multiplicadas por alphai 149 : do ic=1,ndims 150 : do jf=1,ndimx 151 : Kmat0(ndimx*(ic-1)+jf)=Kmat1(ndimx*(ic-1)+jf)*alphai 152 : enddo 153 : enddo 154 : flagini=0 155 : do i=1,ndimx*ndimx 156 : Kmat1(i)=0.D0 157 : enddo 158 : else 159 : c suma (Anterior + P*alphai) [nx.ns+nx.ns] 160 : do ic=1,ndims 161 : Kmat0(ndimx*(ic-1)+ic)=Kmat0(ndimx*(ic-1)+ic)+alphai 162 : enddo 163 : c multiplica Kmat2=Kmat1*Kmat0 [nx.ns=nx.nx*nx.ns]; 164 : c corresponde a JSS^-1*() 165 : do ic=1,ndims 166 : do jf=1,ndimx 167 : do k=1,ndimx 168 : Kmat2(ndimx*(ic-1)+jf)=Kmat2(ndimx*(ic-1)+jf)+ 169 : . Kmat1(ndimx*(k-1)+jf)*Kmat0(ndimx*(ic-1)+k) 170 : enddo 171 : enddo 172 : enddo 173 : c guarda Kmat2 en Kmat0 para el siguiente paso 174 : do i=1,ndimx*ndims 175 : Kmat0(i)=Kmat2(i) 176 : enddo 177 : do i=1,ndimx*ndimx 178 : Kmat1(i)=0.D0 179 : enddo 180 : do i=1,ndimx*ndims 181 : Kmat2(i)=0.D0 182 : enddo 183 : c write(ues,999)((Kmat0(ndims*(jf-1)+i),jf=1,ndims),i=1,ndims) 184 : endif 185 : enddo 186 : c poner las sigma y varint finales 187 : do i=1,ndims 188 : sigfin(i)=sig2(i) 189 : enddo 190 : do i=1,ndimv 191 : varfin(i)=var2(i) 192 : enddo 193 : c Kmat2 = Matriz de hook 194 : ivoid=1 194 : call MatHok(Kmat2,ndimx,ndims,ivoid) 195 : c calcular (Pt*Kmat0)*E [(ns.nx*nx.ns)*ns.ns] 196 : do ic=1,ndims 197 : do jf=1,ndims 198 : Kmat(ndimk*(ic-1)+jf)=0.D0 199 : do k=1,ndims 200 : Kmat(ndimk*(ic-1)+jf)=Kmat(ndimk*(ic-1)+jf)+ 201 : . Kmat0(ndimx*(k-1)+jf)*Kmat2(ndimx*(ic-1)+k) 202 : enddo 203 : enddo 204 : enddo 205 : end 206 : 207 : 208 :
© Cast3M 2003 - All rights reserved.
Disclaimer