1!
2!     CalculiX - A 3-dimensional finite element program
3!              Copyright (C) 1998-2021 Guido Dhondt
4!
5!     This program is free software; you can redistribute it and/or
6!     modify it under the terms of the GNU General Public License as
7!     published by the Free Software Foundation(version 2);
8!
9!
10!     This program is distributed in the hope that it will be useful,
11!     but WITHOUT ANY WARRANTY; without even the implied warranty of
12!     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13!     GNU General Public License for more details.
14!
15!     You should have received a copy of the GNU General Public License
16!     along with this program; if not, write to the Free Software
17!     Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
18!
19      subroutine umat_user(amat,iel,iint,kode,elconloc,emec,emec0,
20     &        beta,xokl,voj,xkl,vj,ithermal,t1l,dtime,time,ttime,
21     &        icmd,ielas,mi,nstate_,xstateini,xstate,stre,stiff,
22     &        iorien,pgauss,orab,pnewdt,ipkon)
23!
24!     calculates stiffness and stresses for a user defined material
25!     law
26!
27!     icmd=3: calculates stress at mechanical strain
28!     else: calculates stress at mechanical strain and the stiffness
29!           matrix
30!
31!     INPUT:
32!
33!     amat               material name
34!     iel                element number
35!     iint               integration point number
36!
37!     kode               material type (-100-#of constants entered
38!                        under *USER MATERIAL): can be used for materials
39!                        with varying number of constants
40!
41!     elconloc(*)        user defined constants defined by the keyword
42!                        card *USER MATERIAL (actual # =
43!                        -kode-100), interpolated for the
44!                        actual temperature t1l
45!
46!     emec(6)            Lagrange mechanical strain tensor (component order:
47!                        11,22,33,12,13,23) at the end of the increment
48!                        (thermal strains are subtracted)
49!     emec0(6)           Lagrange mechanical strain tensor at the start of the
50!                        increment (thermal strains are subtracted)
51!     beta(6)            residual stress tensor (the stress entered under
52!                        the keyword *INITIAL CONDITIONS,TYPE=STRESS)
53!
54!     xokl(3,3)          deformation gradient at the start of the increment
55!     voj                Jacobian at the start of the increment
56!     xkl(3,3)           deformation gradient at the end of the increment
57!     vj                 Jacobian at the end of the increment
58!
59!     ithermal           0: no thermal effects are taken into account
60!                        >0: thermal effects are taken into account (triggered
61!                        by the keyword *INITIAL CONDITIONS,TYPE=TEMPERATURE)
62!     t1l                temperature at the end of the increment
63!     dtime              time length of the increment
64!     time               step time at the end of the current increment
65!     ttime              total time at the start of the current step
66!
67!     icmd               not equal to 3: calculate stress and stiffness
68!                        3: calculate only stress
69!     ielas              0: no elastic iteration: irreversible effects
70!                        are allowed
71!                        1: elastic iteration, i.e. no irreversible
72!                           deformation allowed
73!
74!     mi(1)              max. # of integration points per element in the
75!                        model
76!     nstate_            max. # of state variables in the model
77!
78!     xstateini(nstate_,mi(1),# of elements)
79!                        state variables at the start of the increment
80!     xstate(nstate_,mi(1),# of elements)
81!                        state variables at the end of the increment
82!
83!     stre(6)            Piola-Kirchhoff stress of the second kind
84!                        at the start of the increment
85!
86!     iorien             number of the local coordinate axis system
87!                        in the integration point at stake (takes the value
88!                        0 if no local system applies)
89!     pgauss(3)          global coordinates of the integration point
90!     orab(7,*)          description of all local coordinate systems.
91!                        If a local coordinate system applies the global
92!                        tensors can be obtained by premultiplying the local
93!                        tensors with skl(3,3). skl is  determined by calling
94!                        the subroutine transformatrix:
95!                        call transformatrix(orab(1,iorien),pgauss,skl)
96!
97!
98!     OUTPUT:
99!
100!     xstate(nstate_,mi(1),# of elements)
101!                        updated state variables at the end of the increment
102!     stre(6)            Piola-Kirchhoff stress of the second kind at the
103!                        end of the increment
104!     stiff(21):         consistent tangent stiffness matrix in the material
105!                        frame of reference at the end of the increment. In
106!                        other words: the derivative of the PK2 stress with
107!                        respect to the Lagrangian strain tensor. The matrix
108!                        is supposed to be symmetric, only the upper half is
109!                        to be given in the same order as for a fully
110!                        anisotropic elastic material (*ELASTIC,TYPE=ANISO).
111!                        Notice that the matrix is an integral part of the
112!                        fourth order material tensor, i.e. the Voigt notation
113!                        is not used.
114!     pnewdt             to be specified by the user if the material
115!                        routine is unable to return the stiffness matrix
116!                        and/or the stress due to divergence within the
117!                        routine. pnewdt is the factor by which the time
118!                        increment is to be multiplied in the next
119!                        trial and should exceed zero but be less than 1.
120!                        Default is -1 indicating that the user routine
121!                        has converged.
122!     ipkon(*)           ipkon(iel) points towards the position in field
123!                        kon prior to the first node of the element's
124!                        topology. If ipkon(iel) is smaller than 0, the
125!                        element is not used.
126!
127      implicit none
128!
129      character*80 amat
130!
131      integer ithermal(*),icmd,kode,ielas,iel,iint,nstate_,mi(*),iorien,
132     &  ipkon(*)
133!
134      real*8 elconloc(*),stiff(21),emec(6),emec0(6),beta(6),stre(6),
135     &  vj,t1l,dtime,xkl(3,3),xokl(3,3),voj,pgauss(3),orab(7,*),
136     &  time,ttime,pnewdt
137!
138      real*8 xstate(nstate_,mi(1),*),xstateini(nstate_,mi(1),*)
139!
140!     insert here code to calculate the stresses
141!
142      if(icmd.ne.3) then
143!
144!        insert here code to calculate the stiffness matrix
145!
146      endif
147!
148      return
149      end
150
151
152