Shells: Jim reported following issue to bugs on 8/6/12. Bug reported by Stefano on 3/22/07 "it seems that in MAT_22 for shells the output of history variables #1, #2 and #3 from the d3plot file does not agree with what the Keyword manual says. In particular, - history variable #1 is the tensile matrix mode - history variable #2 is the compressive matrix mode - the tensile fiber mode seems lost The problem is present on the latest 971 R3 beta version, but also on older version of ls970 (I tried .5434 and it was there)." A mat_22 shell fails if all the IP have failed by reaching one of the strength values SC, XT, YT, or YC. Nonlinear shear stress parameter, alpha: I checked the mat22 source code for shells and found that ALPH IS used in the evaluation of tau-bar (Eq. 19.22.2 of the 2006 Theory Manual) and thus has units of 1/(stress^3). On the other hand, there is NO effect of ALPH on the shear stress vs. shear strain relationship. That is, the third equation of Eq. 19.22.1 should not include the last term alpha*tau12^3. (Same for mat_54.) ALPH is set to zero in the little examples that I have. __________________________________________________________ Solids: codes on 01/24/2007 at 08:52:02 (version 10136) Fixed rotational instability in material type 22. Solid elements only. ================================================= From the mat22 (solid element) source code, c delamination c el2(i) =sn2* max(0.0,sig3(i))**2+syz2*sig5(i)**2+szx2* 1 sig6(i)**2-1.0 The material delaminates when el2 is greater than zero. In sets2w2, outside of f3dm22, the solid strenth variables are set as ... sn2 = 1/(SN^2) syz2 = 1/(SYZ^2) szx2 = 1/(SZX^2) (proven by print statement added to f3dm22) KFAIL is for solids only and is used to calculate a hydrostatic pressure in the element after tensile fiber failure (for a volumetric state of compression only, i.e., volume less than initial volume). do i=lft,llt prtr =bulk*log(def(i)) << def is relative volume, thus prtr is neg in compression if (prtr.lt.0.0) then prssur =(1.-ef(i))* min(prtr,0.0) << ef fades to zero 100 steps after fiber tensile failure sig1(i) =ef(i)*sig1(i)+prssur sig2(i) =ef(i)*sig2(i)+prssur sig3(i) =ef(i)*sig3(i)+prssur endif enddo The attached one element model (http://ftp.lstc.com/anonymous/outgoing/jday/composites/mat22.pull_then_compress.k) loaded to failure in fiber tension and then compressed shows the effect of KFAIL. Stresses are output in the material coordinate system: x-stress = fiber stress = material a-direction = global z-direction. Here's what I observe: KFAIL.gt.0: When change in volume first becomes negative subsequent to element failure in fiber tension, pressure immediately is reset to zero and then increases in proportion to further change in volume according to the bulk modulus = KFAIL. The physical basis for the compressive stress discontinuity is not obvious to me. KFAIL.eq.0: Compressive stresses do not undergo a discontinuity, i.e., pressure is not reset to zero when volume change becomes negative. Stresses are computed based on the elastic constants and degraded within 100 steps in case of failure criteria has been reached. (last sentence from stefan.hartmann@dynamore.de) ----------------------------------------------------------------------------- The variable "check" is in the 7th spot which gets labeled as "eff. plastic strain" by LS-PrePost. Do not be misled into thinking it is actually effective plastic strain. The 12 extra history variables follow the 7th spot... ef(nlq),em(nlq),ed(nlq), &efs(nlq),ems(nlq),eds(nlq),el(nlq),els(nlq), &q11(nlq),q12(nlq),q13(nlq),q31(nlq), ef = tensile fiber failure flag em = tensile matrix failure flag ed = compressive matrix failure flag efs ems eds el = delamination failure flag (unique to solids; not mentioned in Users Manual or Theory Manual) els where efs, ems, eds, els are stress scaling factors calculated from ef, em, ed, and el. Those scaling factors are applied to the six components of stress thusly efsels =efs(i)*els(i) sig1(i)= min(sig1(i),efs(i)*sig1(i)) sig2(i)= min(sig2(i),efs(i)*sig2(i)) sig3(i)= min(sig3(i),efsels*sig3(i)) sig4(i)=ems(i)*efs(i)*sig4(i) sig5(i)=efsels*sig5(i) sig6(i)=efsels*sig6(i) if (ed(i)*em(i).eq.0.0) then sig2(i)= min(ems(i)*eds(i)*sig2(i),sig2(i)) endif I believe extra history variable 9 to 12... q11 q12 q13 q31 ... are 4 of the 9 transformation matrix terms. Timestep: It looks like material 22 bases the timestep on the max of the diagonal terms of the constitutive matrix. The constitutive matrix terms are evaluated in subroutine sets22 (called during input) and then a value is stored for the time step calculation. The lines below are from sets22. s11=1./prop(1) s22=1./prop(2) if (prop(3).eq.0.0) prop(3)=prop(2) s33=1./prop(3) s12=-prop(4)*s22 s13=-prop(5)*s33 s23=-prop(6)*s33 s1122=s11*s22 s1123=s11*s23 s2213=s22*s13 s1212=s12*s12 s1223=s12*s23 si=1./(s1122*s33-s1123*s23-s2213*s13-s33*s1212+2.*s1223*s13) cm1=si*(s22*s33-s23*s23) cm3=si*(s1223-s2213) cm4=si*(s33*s11-s13*s13) cm5=si*(s12*s13-s1123) cm6=si*(s1122-s1212) if (nipl.gt.0) then cm(21)=max(cm1,cm4) else cm(21)=max(cm1,cm4,cm6) endif In subroutine shl22s, the value of cm(21) is stored in sndspd which, as you noted, is used in the time step calc. soundspeed = sqrt(cm(21)/rho) time step = tssf * L / soundspeed In the last lines above, we use nipl to determine if we use 2 or 3 terms of the diagonal. I think nipl is the number of lines of material angles that are read from the section_shell data. Therefore, we would always use 3 terms for solid elements. However, it seems that we could use 3 terms for shells too if material angles are not input. This seems a little strange, but perhaps it is intended. Of course, I still can't tell why a smaller EB produces a smaller time step. Perhaps you can get all the numbers and plug them into the equations above. The variables are these: EA=prop(1) EB=prop(2) EC=prop(3) PRBA=prop(4) PRCA=prop(5) PRCB=prop(6) Have fun. -Lee qzmt22 reads the keyword data (unit=txt) and loads the prop array: prop1-8 are "Card 3" in structured input, prop9-16 are "Card 4", and so on. matin loads the cm array by reading from the structured input (unit=txts). sets22 redefines cm45-47 The cm array inside f3dm22 stores the following quantities: cm value stored 1 EA 2 EB 3 EC 4 PRBA 5 PRCA 6 PRCB 7 GAB 8 GBC 9 GCA 10 AOPT 11 MACF 12 13 A1 or V1 or XP 14 A2 or V2 or YP 15 A3 or V3 or ZP 16 BETA or D1 17 D2 18 D3 19 20 21 22 23 24 25 26 SC 27 XT 28 YT 29 YC 30 ALPH 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 1./SN^2 46 1./SYZ^2 47 1./SZX^2 48 KFAIL Any cm index left blank above is not used. I believe the product bqs*dd is the contribution of bulk viscosity to internal energy where bqs is a pressure and dd is a volume strain increment. I believe def(nlq) is the relative volume, i.e., current volume divided by initial volume from the start of simulation. There are 4 modes of failure for mat22 solids, the first 3 of which are mentioned in the Users Manual and are stored as extra history variables 1, 2, and 3. The 4th mode of failure is the delamination mode and its history variable el is stored as the 7th extra history variable. Just based on examination of the source code, I would agree with your assessment of the 0.01 factor, i.e., its purpose is to degrade the stress in a controlled fashion over 100 or more time steps rather than all at once. You can confirm this assessment by changing the hardwired factor from 0.01 to 0.001 or from 0.01 to 0.1 and observing the failure behavior of a single element in fiber tension. It would appear that failure of the mat22 material for solids does not trigger element deletion. http://ftp.lstc.com/anonymous/outgoing/jday/composites/typ1sol.mat22.fail.k