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 calcmechstrain(vkl,vokl,emec,eth,iperturb)
20!
21!     calculates the mechanical strain from the displacement gradients
22!     and the thermal stretches
23!
24      implicit none
25!
26      integer iperturb(*)
27!
28      real*8 elineng(6),vkl(0:3,3),vokl(3,3),emec(6),eth(6),
29     &     wkl(3,3),wokl(3,3)
30!
31!
32!
33c      write(*,*) 'calcmechstrain',iperturb(1),iperturb(2)
34!
35!     subtracting the thermal stretch from the deformation gradients
36!     at the end of the increment
37!
38c      write(*,*) vkl(1,1),eth(1)
39      wkl(1,1)=vkl(1,1)-eth(1)
40      wkl(2,2)=vkl(2,2)-eth(2)
41      wkl(3,3)=vkl(3,3)-eth(3)
42      wkl(1,2)=vkl(1,2)-eth(4)
43      wkl(1,3)=vkl(1,3)-eth(5)
44      wkl(2,3)=vkl(2,3)-eth(6)
45      wkl(2,1)=vkl(2,1)-eth(4)
46      wkl(3,1)=vkl(3,1)-eth(5)
47      wkl(3,2)=vkl(3,2)-eth(6)
48!
49!     attention! elineng(4),elineng(5) and elineng(6) are engineering strains!
50!
51      elineng(1)=wkl(1,1)
52      elineng(2)=wkl(2,2)
53      elineng(3)=wkl(3,3)
54      elineng(4)=wkl(1,2)+wkl(2,1)
55      elineng(5)=wkl(1,3)+wkl(3,1)
56      elineng(6)=wkl(2,3)+wkl(3,2)
57c      write(*,*) 'elineng,wkl ',elineng(1),wkl(1,1)
58!
59      if(iperturb(2).eq.1) then
60!
61!     Lagrangian strain
62!
63         elineng(1)=elineng(1)+
64     &        (wkl(1,1)**2+wkl(2,1)**2+wkl(3,1)**2)/2.d0
65         elineng(2)=elineng(2)+
66     &        (wkl(1,2)**2+wkl(2,2)**2+wkl(3,2)**2)/2.d0
67         elineng(3)=elineng(3)+
68     &        (wkl(1,3)**2+wkl(2,3)**2+wkl(3,3)**2)/2.d0
69         elineng(4)=elineng(4)+wkl(1,1)*wkl(1,2)+wkl(2,1)*wkl(2,2)+
70     &        wkl(3,1)*wkl(3,2)
71         elineng(5)=elineng(5)+wkl(1,1)*wkl(1,3)+wkl(2,1)*wkl(2,3)+
72     &        wkl(3,1)*wkl(3,3)
73         elineng(6)=elineng(6)+wkl(1,2)*wkl(1,3)+wkl(2,2)*wkl(2,3)+
74     &        wkl(3,2)*wkl(3,3)
75!
76!     for frequency analysis or buckling with preload the
77!     strains are calculated with respect to the deformed
78!     configuration
79!
80      elseif(iperturb(1).eq.1) then
81!
82!     subtracting the thermal stretch from the deformation gradients
83!     at the start of the increment
84!
85         wokl(1,1)=vokl(1,1)-eth(1)
86         wokl(2,2)=vokl(2,2)-eth(2)
87         wokl(3,3)=vokl(3,3)-eth(3)
88         wokl(1,2)=vokl(1,2)-eth(4)
89         wokl(1,3)=vokl(1,3)-eth(5)
90         wokl(2,3)=vokl(2,3)-eth(6)
91         wokl(2,1)=vokl(2,1)-eth(4)
92         wokl(3,1)=vokl(3,1)-eth(5)
93         wokl(3,2)=vokl(3,2)-eth(6)
94!
95         elineng(1)=elineng(1)+wokl(1,1)*wkl(1,1)+wokl(2,1)*wkl(2,1)+
96     &        wokl(3,1)*wkl(3,1)
97         elineng(2)=elineng(2)+wokl(1,2)*wkl(1,2)+wokl(2,2)*wkl(2,2)+
98     &        wokl(3,2)*wkl(3,2)
99         elineng(3)=elineng(3)+wokl(1,3)*wkl(1,3)+wokl(2,3)*wkl(2,3)+
100     &        wokl(3,3)*wkl(3,3)
101         elineng(4)=elineng(4)+wokl(1,1)*wkl(1,2)+wokl(1,2)*wkl(1,1)+
102     &        wokl(2,1)*wkl(2,2)+wokl(2,2)*wkl(2,1)+
103     &        wokl(3,1)*wkl(3,2)+wokl(3,2)*wkl(3,1)
104         elineng(5)=elineng(5)+wokl(1,1)*wkl(1,3)+wokl(1,3)*wkl(1,1)+
105     &        wokl(2,1)*wkl(2,3)+wokl(2,3)*wkl(2,1)+
106     &        wokl(3,1)*wkl(3,3)+wokl(3,3)*wkl(3,1)
107         elineng(6)=elineng(6)+wokl(1,2)*wkl(1,3)+wokl(1,3)*wkl(1,2)+
108     &        wokl(2,2)*wkl(2,3)+wokl(2,3)*wkl(2,2)+
109     &        wokl(3,2)*wkl(3,3)+wokl(3,3)*wkl(3,2)
110      endif
111!
112!     storing the local strains
113!
114      if(iperturb(1).ne.-1) then
115         emec(1)=elineng(1)
116         emec(2)=elineng(2)
117         emec(3)=elineng(3)
118         emec(4)=elineng(4)/2.d0
119         emec(5)=elineng(5)/2.d0
120         emec(6)=elineng(6)/2.d0
121      else
122!
123!        linear iteration within a nonlinear increment:
124!
125!        subtracting the thermal stretch from the deformation gradients
126!        at the start of the increment
127!
128         wokl(1,1)=vokl(1,1)-eth(1)
129         wokl(2,2)=vokl(2,2)-eth(2)
130         wokl(3,3)=vokl(3,3)-eth(3)
131         wokl(1,2)=vokl(1,2)-eth(4)
132         wokl(1,3)=vokl(1,3)-eth(5)
133         wokl(2,3)=vokl(2,3)-eth(6)
134         wokl(2,1)=vokl(2,1)-eth(4)
135         wokl(3,1)=vokl(3,1)-eth(5)
136         wokl(3,2)=vokl(3,2)-eth(6)
137!
138         emec(1)=wokl(1,1)+
139     &        (wokl(1,1)**2+wokl(2,1)**2+wokl(3,1)**2)/2.d0
140         emec(2)=wokl(2,2)+
141     &        (wokl(1,2)**2+wokl(2,2)**2+wokl(3,2)**2)/2.d0
142         emec(3)=wokl(3,3)+
143     &        (wokl(1,3)**2+wokl(2,3)**2+wokl(3,3)**2)/2.d0
144         emec(4)=(wokl(1,2)+wokl(2,1)+wokl(1,1)*wokl(1,2)+
145     &        wokl(2,1)*wokl(2,2)+wokl(3,1)*wokl(3,2))/2.d0
146         emec(5)=(wokl(1,3)+wokl(3,1)+wokl(1,1)*wokl(1,3)+
147     &        wokl(2,1)*wokl(2,3)+wokl(3,1)*wokl(3,3))/2.d0
148         emec(6)=(wokl(2,3)+wokl(3,2)+wokl(1,2)*wokl(1,3)+
149     &        wokl(2,2)*wokl(2,3)+wokl(3,2)*wokl(3,3))/2.d0
150      endif
151c      write(*,*) 'emec ',emec(1)
152!
153      return
154      end
155