NOTE: http://ftp.lstc.com/anonymous/outgoing/support/FAQ/user-material-notes is provided by a colleague and also contains FAQ-type information about user-defined materials. That information is likely to be more comprehensive and/or up-to-date than the notes below. I suggest you refer to the aforementioned file first, and then if you don't find an answer to your question, refer to the notes below. (jd 4/04/18) _____________________________________________________________________ Subject: User-Defined Materials To invoke a user-defined material, you must ... 1. Create a custom executable which includes your material subroutine. The "Makefile" supplied by LSTC should enable you to create this executable by simply typing the word "make". 2. Invoke that subroutine by defining a part in your input deck that uses *mat_user_defined_material_models. To get started, you'll need to download a file from http://ftp.lstc.com/objects . (Contact lstsupport@ansys.com for the required username and password.) The file you need will depend on your hardware and operating system. The file contains object code, source code in a file called dyn21.f, and a Makefile to compile your custom version of LS-DYNA. There is a sample user material routine in Appendix A of the User's Manual. Appendix A includes a description of how to calculate the deformation gradient matrix which would be required for a material model based on total (as opposed to incremental) strain. For an example of a complete user-defined material subroutine, see subroutine umat41 in dyn21.f. This user-defined material behaves like *mat_elastic but is called via *mat_user_defined_material_models with MT set to 41. An example input deck that calls umat41 is http://ftp.lstc.com/anonymous/outgoing/support/FAQ_kw/2x_sphere2plate.k.gz During development and testing of a user-defined subroutine, you should make liberal use of print statements within the subroutine and run several simple test models to confirm that variables are what you think they are. Generally speaking, every new release of LS-DYNA carries with it changes to dyn21.f. A 2-day class on user-defined materials, taught by Dr. Ala Tabiei (atabiei@aol.com), is offered periodically. Click on "Classes" at www.lstc.com for a schedule of upcoming classes. The class is also now offered on-line as well. For information on course notes and the on-line version of the class, go to http://www.lsdyna-online.com/ A short preview of the on-line course is available at either of these sites... http://www.lsdyna-online.com/tutorial-movies.html http://www.youtube.com/user/lsdynatv FREQUENTLY ASKED QUESTIONS _______________________________________________________________________________ Question: Can I add a common block? Answer (7/2009): The user can add new threadprivate common blocks in usercomm of dyn21.F. It is assumed the user understands the difference between global and threadprivate common blocks as it relates to SMP (Shared Memory Parallel) programming. For SMP, one must initialize all threadprivate common blocks in the beginning of the run. See the notes in subroutine usercomm. _______________________________________________________________________________ Question: Is it strain increment or strain rate that's passed into the user-defined material subroutines? Answer: It depends on the version of the code. In version 970, it's strain increment for shells and strain rate for solids. In version 971, it's strain increment for both shells and solids. _______________________________________________________________________________ Question: In what coordinate system are stresses and strains within a user_defined material subroutine? Answer: Non-orthotropic solids: global coordinate system Non-orthotropic shells: local element coordinate system Orthotropic solids and shells: local material coordinate system _______________________________________________________________________________ Question: Is there a difference in UMAT41-UMAT50? Is the routine name of UMAT only connected with the parameter MT in *MAT_USER_DEFINED_MATERIAL_MODEL? Answer: Avoid using umat42 and the vectorized routines umat48v and umat49v; these are special, undocumented subroutines that do not reside in dyn21.f and thus the user cannot modify them. Aside from those exceptions, it doesn't matter which MT (from 41 to 50) you choose. The umat subroutine used is controlled directly by MT and by IVECT. If IVECT is set to 1, the "v" (vectorized) version of the umat subroutine will be called, e.g., subroutine umat41v is called if MT=41 and IVECT=1. _______________________________________________________________________________ Question: What are the differences between UMAT41 and UMAT41V? How are UMAT41 and UMAT41V used properly? Answer: The vectorized umat routines process elements in groups whereas the nonvectorized umat routines process one element at a time. As you might imagine, the vectorized routine is much more efficient. The number of elements processed each time a vectorized umat routine is called the "vector length". This vector length is given in the file nlqparm which must be 'included' in the umatxxv subroutine as follows: include 'nlqparm' In vectorized umat routines (as in umat41v in dyn21.f), a do loop from lft to llt is used to loop over the group of elements being processed. _______________________________________________________________________________ Question: What is the limit for the number of input constants? Answer: The limit is 48 for isotropic materials, 40 for orthotropic materials. To obtain a greater number of constants, you could hardwire the input values but this would require recompiling in order to change the input. A better way is to read input data from an auxilary input file, e.g., ... if (hist(1).lt.0.1) then c ...read inital values for material open(11, file='material.ini') do 1 i=1,cm(1) read (11,*) cm(i+1) 1 continue close (11) c ...initialization done hist(1)=1 endif _______________________________________________________________________________ Question: Is there a standard way to send (error) messages from the material subroutine to the standard output file? Answer: Write error messages to unit 59 (messag) and/or unit 13 (d3hsp). _______________________________________________________________________________ Question: How do I obtain element IDs, element connectivity, nodal coordinates, etc. from within my user-defined material subroutine? Answer: Functions that will return node IDs and element IDs: Internal node or element ID = lqfint(a,b,c) where a = external ID b = data type: 1 = node 2 = brick 3 = beam 4 = shell 5 = thick shell c = returned error flag: 0 = ID found 1 = ID not found Exteral element ID = lqfinvf(internal_element_ID,ityp) c c ityp=2 solid c ityp=3 beam c ityp=4 shell c ityp=5 thick shell External node ID = lqfinv(internal_node_id, 1) ------------------------------- Shortly after the release of 971 R5.0, Chen prepared a subroutine lqfnode which returns nodal coordinates and internal element ID if the external element ID is known. That subroutine, along with a subroutine umat41v that demonstrates its use, are included in dyn21.extIDin.f. The print statements therein are specific to the test case 2x_sphere2plate.k. I made some slight modifications to subroutine lqfnode and the aforementioned umat41v so that nodal coordinates and external element ID will be returned if the internal element ID is known. Those modifications are included in dyn21.intIDin.f. Note that the argument nnm1 was added to subroutine umat41v and to the calls to umat41v in subroutines urmats and urmathn. For reference, dyn21.f from ls971_d_R5_0_intel64_redhat54.tar.gz is the standard source code from which dyn21.extIDin.f and dyn21.intIDin.f sprang. All aforementioned files are in http://ftp.lstc.com/anonymous/outgoing/support/FAQ_kw/lqfnode_example.tar.gz For older releases of LS-DYNA (than 971 R5.0), read on. Subroutine umat41 in dyn21.f of version 970 or subroutine userrrefin in dyn21.f of version 971 illustrates how to access IDs, coords, connectivity. Here's an excerpt from 970 source code... From within umat41-calling subroutine urmathn: include 'mema.inc' integer i_mem common/dynmem/i_mem(1) real r_mem(1) equivalence (i_mem,r_mem) integer*4 dm_x, dm_v, dm_xms, dm_me1 integer*4 dm_x0,dm_v0,dm_xms0,dm_me10 common /dynmem1/ dm_x, dm_v, dm_xms, dm_me1, & dm_x0,dm_v0,dm_xms0,dm_me10 common/bk13/lc0,lc1h,lc1b,lc1s,lc1t,lc2,lc3,lc4,lc5,lc6,lc7,lc9, 1 lc10,lc11,lc12,lc13,lc14,lc15,lc16,lc17,lc18,lb0,lb1,lb2, 2 lc7a,lc7b,lc7c,lc7d,lc7e,lc7f,lc7g,lc7h,lc7i,lc7j,lc7k,lc7l ... 41 call umat41 (cm(mx+1),eps,sig,hsv,dt1,capa,'shell',tt, . temper,nnm1+i,a(lc1s),r_mem(dm_x),5,i) Then in umat41: subroutine umat41 (cm,eps,sig,hisv,dt1,capa,etype,time, . temp,i,ixs,x,k,j) ... c i = internal element ID c iext = external element ID c ip = internal part ID c ipext = external part ID c n1,n2,... = internal node ID for element connectivity c next1, next2,... = external node ID for element connectivity c x1,y1,z1 = coordinates of first node in the connectivity c x2,y2,z2 = coordinates of second node in the connectivity c etc. ... dimension ixs(k,*),x(3,*) ... c 2nd argument of function nelmntid is 0 for solids, 1 for beams, 2 for shells, c 3 for tshells iext=nelmntid(i,0) ip=ixs(1,i) n1=ixs(2,i) n2=ixs(3,i) n3=ixs(4,i) n4=ixs(5,i) n5=ixs(6,i) n6=ixs(7,i) n7=ixs(8,i) n8=ixs(9,i) c external part and node ids ipext=lqfmiv(ip) next1=lqfinv(n1,1) next2=lqfinv(n2,1) next3=lqfinv(n3,1) next4=lqfinv(n4,1) next5=lqfinv(n5,1) next6=lqfinv(n6,1) next7=lqfinv(n7,1) next8=lqfinv(n8,1) c spot check coordinates x8=x(1,n8) y8=x(2,n8) z8=x(3,n8) ... _______________________________________________________________________________ Question: How are history variables which are defined in the user-defined material subroutine accessed when postprocessing? Answer: In pre-970 versions of LS-DYNA, the 1st history variable in the umat subroutine won't be stored as history var#1 in the d3plot database. The storage location is dependent on a number of factors such as whether the subroutine is vectorized or nonvectorized, whethere the elements are shells or solids, etc. More on this subject comes from developer Lee Bindeman: "When using a vectorized subroutine (i.e. umat46v instead of umat46) and a 3D user defined material (for 3D solid elements), there are 6 history variables automatically used for 6 terms of the transformation matrix, whether or not the user defined material is orthotropic. When the material is orthotropic (IORTHO=1), these 6 variables are automatically allocated, however, when the material is not orthotropic (IORTHO=0), the variables are not and it is essential that these be allocated by the user defined material input. Therefore, if your material uses 46 history variables, you need to set NHV=52. In order to write your 46 history variables to the d3plot files, you need to request 52 extra history variables by setting NEIPH=52 on *DATABASE_EXTENT_BINARY. When post processing, history variables 1 to 6 will contain the transformation matrix terms. If the material is isotropic, these will all be zero. History variables 7 to 52 will contain history variables 1 to 46 in your subroutine. The above rules change for 2D materials (for shell elements). In that case, there are only 2 transformation terms stored so only 2 extra history variables need to be allocated and requested. Perhaps some good news is that I made a fix ... to eliminate much of this confusion. The fix is in version 970 revisions 2903 and later. With this fix, it is no longer be necessary to allocate extra history variables with NHV, or to request 6 (or 2) extra history variables in the d3plot file. You simply need to allocate the same number of history variables that you want to use, and to request the number that you want written to the d3plot file. If the material is isotropic (IORTHO=0), then the transformation terms will be omitted and then the history variable numbers in the user subroutine will match those in the d3plot file. However, if the material is orthotropic (IORTHO=1), then the 6 (or 2) transformation terms will be written to the d3plot file so there will be a mismatch in the numbering of history variables in the user subroutine and the d3plot file." _______________________________________________________________________________ Question: Is there a limit to the number of extra history variables? Answer: (7/2009): Chen Tsay modified dynm.F and dyn21.F for ls971.r4, r4.2 and r4.2.1 to allow the user to modify the maximum number of history variables (default is 142). The user can now define the maximum number of history variable by NHISVAR in nhisparm.inc. Answer (4/2006): Looking at subroutines urmats (shells) and urmathn (solids) in dyn21, it appears that the limit for extra history variables is... For 971: 97 for shells. Subtract 2 if orthotropic. Subtract 9 if hyperelastic. 117 for solids. Subtract 6 if orthotropic. Subtract 9 if hyperelastic. For 970: 100 for shells. Presumably the same reductions for orthotropic and hyperelastic apply. 54 for solids. Ditto. _______________________________________________________________________________ Question: How do I flag an element as being 'failed' and delete the element from the simulation? Answer: The following vectorized umat subroutine illustrates how to implement failure for solid and shell elements. Note the added common blocks. subroutine umat41v(cm,d1,d2,d3,d4,d5,d6,sig1,sig2, . sig3,sig4,sig5,sig6,sigma,hist,lft,llt,dt1,capa, . etype,tt,temp) c . etype,tt,temp,nnm1) c****************************************************************** c| livermore software technology corporation (lstc) | c| ------------------------------------------------------------ | c| copyright 1987-1999 | c| all rights reserved | c****************************************************************** c c isotropic elastic material (sample user subroutine) c c variables c c cm(1)=young's modulus c cm(2)=poisson's ratio c c d1(i)=local x strain rate for solids (incremental strain for shells) c d2(i)=local y ... c d3(i)=local z ... c d4(i)=local xy ... c d5(i)=local yz ... c d6(i)=local zx ... c c sig1(i)=total Cauchy stress in local x (shells) c sig2(i)=total Cauchy stress in local y (shells) c sig3(i)=total Cauchy stress in local z (shells) c sig4(i)=total Cauchy stress in local xy (shells) c sig5(i)=total Cauchy stress in local yz (shells) c sig6(i)=total Cauchy stress in local zx (shells) c For solids, replace sig1(i) with sigma(i,1), sig2(i) with sigma(i,2), etc. c c davg is the volumetric strain rate for solids c davg is the incremental volumetric strain for shells c p is the incremental pressure c c hist(i,1)=1st history variable c hist(i,2)=2nd history variable c . c . c hist(i,n)=nth history variable c c dt1=current time step size c capa=reduction factor for transverse shear c etype: c eq."brick" for solid elements c eq."shell" for all shell elements c eq."beam" for all beam elements c c tt=current problem time. c c c all transformations into the element local system are c performed prior to entering this subroutine. transformations c back to the global system are performed after exiting this c routine. c c all history variables are initialized to zero in the input c phase. initialization of history variables to nonzero values c may be done during the first call to this subroutine for each c element. c c energy calculations for the dyna3d energy balance are done c outside this subroutine. c c add failure for bricks and shells include 'nlqparm' cC_TASKCOMMON (subtssloc) C_TASKCOMMON (failcmloc) C_TASKCOMMON (failuloc) c common /subtssloc/ dt1siz(nlq) c logical failur,failgl common/failcm/failur,failgl common/failcmloc/ifail(nlq) common/failuloc/sieu(nlq),fail(nlq) character*(*) etype dimension cm(*) dimension d1(*),d2(*),d3(*),d4(*),d5(*),d6(*) dimension sig1(*),sig2(*),sig3(*),sig4(*),sig5(*),sig6(*) dimension sigma(nlq,*),hist(nlq,*) c c compute shear modulus, g c g2 =cm(1)/(1.+cm(2)) g =.5*g2 c c ====================================== c internal element number ieint c c do i=lft,llt c ieint=nnm1+i cc external element number ieext cc c if (etype.eq.'brick') then c ieext=nelmntid(ieint,0) c elseif (etype.eq.'shell') then c ieext=nelmntid(ieint,2) c elseif (etype.eq.'beam') then c ieext=nelmntid(ieint,1) c endif cc checking cc c if(tt.le.1.e-8) then c print *,'internal elem # = ',ieint,' external elem # =',ieext c endif c enddo cc c ====================================== c do i=lft,llt if (etype.eq.'brick') then davg=(-d1(i)-d2(i)-d3(i))/3. p=-davg*dt1*cm(1)/((1.-2.*cm(2))) sigma(i,1)=sigma(i,1)+p+g2*(d1(i)+davg)*dt1 sigma(i,2)=sigma(i,2)+p+g2*(d2(i)+davg)*dt1 sigma(i,3)=sigma(i,3)+p+g2*(d3(i)+davg)*dt1 sigma(i,4)=sigma(i,4)+g*d4(i)*dt1 sigma(i,5)=sigma(i,5)+g*d5(i)*dt1 sigma(i,6)=sigma(i,6)+g*d6(i)*dt1 c check for brick failure c use hardwired time-based failure criterion just for illustration if (tt.gt..0008) then sigma(i,1)=0. sigma(i,2)=0. sigma(i,3)=0. sigma(i,4)=0. sigma(i,5)=0. sigma(i,6)=0. c ifail necessary to fail and delete a brick; not sure what fail does failur=.true. ifail(i)=1 fail(i) =0. endif c elseif (etype.eq.'shell') then gc =capa*g q1 =cm(1)*cm(2)/((1.0+cm(2))*(1.0-2.0*cm(2))) q3 =1./(q1+g2) d3(i)=-q1*(d1(i)+d2(i))*q3 davg =(-d1(i)-d2(i)-d3(i))/3. p =-davg*cm(1)/((1.-2.*cm(2))) sig1(i)=sig1(i)+p+g2*(d1(i)+davg) sig2(i)=sig2(i)+p+g2*(d2(i)+davg) sig3(i)=0.0 sig4(i)=sig4(i)+g *d4(i) sig5(i)=sig5(i)+gc*d5(i) sig6(i)=sig6(i)+gc*d6(i) c check for shell failure c use hardwired time-based failure criterion just for illustration c Note that the failed shell element will only be deleted if IFAIL in c *MAT_USER... is set to 1. c if (ipt.eq.1) then c fail(i)=0.0 c endif if (tt.gt..0007) then sig1(i)=0. sig2(i)=0. sig4(i)=0. sig5(i)=0. sig6(i)=0. failur=.true. ifail(i)=1 fail(i)=0. c else c fail(i)=1. endif c endif enddo c return end For another response to this question, see "Yahoo Group Yammerings" in the May 2006 FEAInformation newsletter (www.feainformation.com). _______________________________________________________________________________ Question: How do I access data from a curve (defined via *define_curve) in my user-defined material subroutine? Answer: Update: Words were added to the R7 User's Manual (Appendix A) on 5/6/13 addressing use of curves and tables in user-defined material subroutines. Old synopsis: There are two methods. Method 1: Curves are internally rediscretized so that there are 100 points spaced equally along the abscissa. The rediscretized abscissa values are stored in crv(1 to 100, 1, internal_curve_id) and the rediscretized ordinate values are stored in crv(1 to 100, 2, internal_curve_id). crv(101,1,*) is the abscissa increment between each point. * is the internal load curve ID which is equal to lcids(nint(id)) where id is the user-specified (external) curve ID. Through use of print statements, you can confirm the values in the crv array. The subroutine umat43 below is a very simple nonlinear elastic model (nonlinear only in the y-direction). The 5th material constant is a curve defining a scale factor versus y-strain where this scale factor is applied to the linear elastic stress increment to get the nonlinear stress increment. Note that the curve is internally rediscretized by LS-DYNA using this approach. This example material model is ... - nonvectorized - for solid elements - for illustration only and does not represent a real material - an example last tested in v. 950c of LS-DYNA. Some slight modifications may be necessary in later versions of the code but the same basic approach should work. First, the calling subroutine (urmath) is modified to include a pointer to the curve data inside the "a" array. This done by including the following common block: common/slcntr/islcnt(21) ... and adding an argument to the call statement ... 43 call umat43 (cm(mx+1),eps,sig,hsv,dt1,capa,'brick',tt, . a(islcnt(16))) The user-defined material subroutine then looks like this... subroutine umat43 (cm,eps,sig,hisv,dt1,capa,etype,time,crv) c eps is a strain increment c crv contains the curve data c character*(*) etype dimension cm(*),eps(*),sig(*),hisv(*),crv(101,2,*) c c compute shear modulus, g c g2 =cm(1)/(1.+cm(2)) g =.5*g2 c cm(5) = external curve id c isl = internal curve id c is1=lcidi(nint(cm(5))) is1=lcids(nint(cm(5))) c print *,cm(5),is1 c c let 1st history variable be total strain in y direction c hisv(1) is stored in "effective plastic strain" location in d3plot hisv(1)=hisv(1)+eps(2) c show that hisv(2) is stored in "1st additional history variable" location c in d3plot c (give hisv(2) some arbitrary value, say twice the y-strain) hisv(2)=hisv(1)*2. c given x-value (eps(2)), extract y-value (sfact) from curve if(is1.gt.0) then crvinc=crv(101,1,is1) crvbgn=crv(1,1,is1) c do i=1,100 c print *,crv(i,1,is1),crv(i,2,is1) c enddo ni1 =(abs(hisv(1))-crvbgn)/crvinc+1 ni1 =min(ni1,99) ni =ni1+1 slop=(crv(ni,2,is1)-crv(ni1,2,is1))/(crv(ni,1,is1)-crv(ni1,1,is1)) sfact = crv(ni1,2,is1)+ slop*(abs(hisv(1))-crv(ni1,1,is1)) c print *,ni1,ni,slop,hisv(1),sfact endif c davg=(-eps(1)-eps(2)-eps(3))/3. p=-davg*cm(1)/((1.-2.*cm(2))) sig(1)=sig(1)+p+g2*(eps(1)+davg) c sfact applied to sig(2) only (nonlinear only in y) sig(2)=sig(2)+p+g2*(eps(2)+davg)*sfact sig(3)=sig(3)+p+g2*(eps(3)+davg) sig(4)=sig(4)+g*eps(4) sig(5)=sig(5)+g*eps(5) sig(6)=sig(6)+g*eps(6) c c return end ----------------------------------------------- Method 2: To access curve data which has not been rediscretized, refer to the comments in subroutine loadud (v. 970): c p - load curve data pairs (abcissa,ordinate) c npc - pointer into p. (p(npc(lc)) points to the beginning c of load curve ID lc. npoints=npc(lc+1)-npc(lc)= c number of points in the load curve. (Actually, I think npoints is the number of values stored for the curve where each point consists of two values -- an abscissa value and an ordinate value.) The curve data in the p array has not been rediscretized, i.e, it is the original curve data as input by the user. lc is the external load curve ID. Now look at subroutine usrmat. The arguments plc and npc in usrmat are the same as p and npc in subroutine loadud. I hope you can figure it out from here. As a check, print out the values from plc(npc(lc)) to plc(npc(lc) + npoints - 1) to confirm that the curve data is as you expect it to be. ______________________________________________________________________________ Question: How can I get shell thickness from within a umat routine? The subroutine urmats has an example of how to access the array of shell thicknesses - search for "ns05" in that routine and you will see code like this: do i=lft,llt hsvs(i,ind_thick )=a(ns05+9*(nnm1+i-1) ) hsvs(i,ind_thick+1)=a(ns05+9*(nnm1+i-1)+1) hsvs(i,ind_thick+2)=a(ns05+9*(nnm1+i-1)+2) hsvs(i,ind_thick+3)=a(ns05+9*(nnm1+i-1)+3) enddo that is, the shell thicknesses are in the array a(ns05) and are accessed as shown above. The pointer ns05 is included in the subroutine through this statement: #include "bk05.inc" Hughes-Liu shells use all 4 values but Belytschko-Tsay shells use only the first one, set to the average of the four thicknesses. ub 1/6/2016 --------------------------------------------- (The following, older description is perhaps less clear than the one above.) The fibl array is the way to go. Since it's not in dyn21, the easiest way to get the data is to use the memory pointer for fibl which is in a common block, ns05, which is already in dyn21 for subroutine useradap. I would suggest putting bk05 in usrmat along with the include statement for mema.inc. You could then add an argument to the call to urmats which would be a(ns05). In urmats, the variable would be named fibl with dimension fibl(9,*). To get data for the current element, you you need to use nnm1+i as the offset in the fibl array. By the way, although fibl(1,i) to fibl(4,i) are thickness values, I think Belytschko-Tsay and variants use only fibl(1,i) as the thickness. If the user inputs different thickness values, fibl(1,i) gets set to the average. In that case, I don't know what is in the other fibl values that are not used. I think Hughes-Liu uses 4 values and they can be unique. I'm not sure what the type 16 element does. If you like, I can check particular element forms since it may be important for the user to know for their particular element type. lpb 11/2004 ______________________________________________________________________________ Question: How can I access initial density and current volume? Answer (from Yong Guo, 8/2005): The initial density is accessed as rhoa(lft) in common block /aux35loc/ (note that only the first value in the array gives a meaningful value). xm(nlq) in common block /aux43loc/ gives us the inverse of the initial volume. The current volume can be accessed with voln(nlq) in common block /aux9loc/ ... common/aux9loc/vlrho(nlq),voln(nlq) ______________________________________________________________________________ Question: How do I write a umat subroutine for a hyperelastic material? Answer: I'd suggest you use v. 971 for this purpose as things are a bit clearer. The attached subroutine (umat45) is taken directly from dyn21.F of version 971. Appendix A in the DRAFT version 971 User's Manual explains this subroutine. Thomas Borrvall (thomas.borrvall@erab.se) revised this appendix for v. 971. Contact support@lstc.com for the latest pdf of the Users Manual. Note that in writing a umat subroutine for a hyperelastic material, you need to have access to the deformation gradient matrix. Note the parameter IHYPER in the *MAT_USER_... input. ______________________________________________________________________________ Question: How do I write a umat subroutine for implicit analysis? Answer: To start, read the section on implicit in Appendix A of the 971 Users Manual, available at www.lstc.com. You may want to contact support@lstc.com to see if a more recent draft of the Users Manual is available. http://ftp.lstc.com/anonymous/outgoing/support/FAQ_kw/implicit_elastic_umat41.k is an example input deck. This simple, elastic example calls subroutine umat41 and, because the solution is implicit, it also calls subroutine utan41. There are more complex examples of utan** subroutines in dyn21.f as well. _____________________________________________________________________________ Question: What measure of incremental strain is passed into the user-defined subroutine? Answer: The incremental strains passed in the example subroutine umat41 are defined as: eps(1) = du/dx eps(2) = dv/dy eps(3) = dw/dz eps(4) = du/dy+dv/dx eps(5) = dv/dz+dw/dy eps(6) = du/dz+dw/dx where (u,v,w) are the displacements and (x,y,z) are spatial coordinates. (Borrvall, 3/7/07) In the case of the shear strains, i.e., eps4,5,6, these passed values are the engineering shear strains or twice the tensorial shear strains. _____________________________________________________________________________ Question: How do I control the number of integration points in an element that must fail before the element is deleted? Answer: For IFAIL=1 only one integration point must fail before element deletion occurs. In 971 R5 beta, you now have the option to define element failure based on a user-specified number of failed integration points. On *MAT_USER_DEFINED_MATERIAL_MODELS the parameter IFAIL may be set to a negative number, which then means that -IFAIL is a pointer into the material parameters array (just as IBULK and IG). The corresponding material parameter should then specify the number of integration points that must fail before the element is deleted. In the code, the logical variables failel/failels(*) get a new interpretation, they are now input/output variables to the user defined material routines. On input they tell whether this integration point has failed (.true.) or not (.false.) at an earlier point in time, and should be used in a conditional to simply skip the stress update. If .false. on input, the stress update takes place as usual and the user may put failel=.true. as output if a failure condition is met, which means that this integration point is taken out from the remainder of the computations. A counter keeps track of the number of failed integration points and the element is deleted when this number reaches the user specified number. The feature applies currently to brick and shell elements. http://ftp.lstc.com/anonymous/outgoing/support/FAQ_kw/usermat_failipts.k is a rather silly example which at least shows that elements with "numint=1" are eroded before elements with "numint=6". (Borrvall, 8/28/09) Thomas' email to qa: usermat_failipts.k (added 28/8 2009) Illustrates how the user may defined failure in user materials based on number of integration points. In the model, four cantilever beams are present and the same failure criterion is on a small part close to the "wall". The difference being that for two beams erosion occurs when first integration point fails and the other two when six integration points have failed. The result is that the elements are eroded at different points in time. This problem only applies to R4 development rev. 55026 and later. _____________________________________________________________________________ Writing or debugging user-subroutines is not within the normal realm of tech support. It is recommended that users interested in writing a user-defined material model attend the corresponding short course.