OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
c3lke3.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| c3lke3 ../engine/source/elements/sh3n/coque3n/c3lke3.F
25!||--- called by ------------------------------------------------------
26!|| c3ke3 ../engine/source/elements/sh3n/coque3n/c3ke3.f
27!||====================================================================
28 SUBROUTINE c3lke3(JFT,JLT,AREA,THK0,THK2,HM,HF,HC,HZ,
29 1 PX1,PY1,PY2,VOL,
30 2 K11,K12,K13,K22,K23,K33,
31 3 M11,M12,M13,M22,M23,M33,
32 4 MF11,MF12,MF13,MF22,MF23,MF33,
33 5 FM12,FM13,FM23,IKGEO,FOR ,MOM ,
34 6 IORTH,HMOR,HFOR,HMFOR,IDRIL,
35 7 NEL)
36C--------------------------------------------------------------------------------------------------
37C-----------------------------------------------
38C I M P L I C I T T Y P E S
39C-----------------------------------------------
40#include "implicit_f.inc"
41#include "mvsiz_p.inc"
42#include "com04_c.inc"
43C-----------------------------------------------
44C D U M M Y A R G U M E N T S
45C-----------------------------------------------
46C---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8
47 INTEGER JFT,JLT,IKGEO,IORTH,IDRIL,NEL
48 MY_REAL
49 . AREA(*),VOL(*), PX1(*),PY1(*),PY2(*),THK0(*),THK2(*),
50 . HM(MVSIZ,4),HF(MVSIZ,4),HC(MVSIZ,2),HZ(*),HMOR(MVSIZ,2),HFOR(MVSIZ,2),
51 . K11(3,3,*),K12(3,3,*),K13(3,3,*),
52 . K22(3,3,*),K23(3,3,*),K33(3,3,*),
53 . M11(3,3,*),M12(3,3,*),M13(3,3,*),
54 . m22(3,3,*),m23(3,3,*),m33(3,3,*),
55 . mf11(3,3,*),mf12(3,3,*),mf13(3,3,*),
56 . mf22(3,3,*),mf23(3,3,*),mf33(3,3,*),
57 . fm12(3,3,*),fm13(3,3,*),fm23(3,3,*),for(nel,5),mom(nel,3),
58 . hmfor(mvsiz,6)
59C-----------------------------------------------
60C L O C A L V A R I A B L E S
61C-----------------------------------------------
62 INTEGER EP,I,J,L,K,I1,J1,IN(2),NF
63 MY_REAL
64 . DM(2,2,MVSIZ),DF(2,2,MVSIZ),DCX(MVSIZ),DCY(MVSIZ),
65 . PX1PX1(MVSIZ),PX1PY1(MVSIZ),PX1PY2(MVSIZ),PX1PY3(MVSIZ),
66 . PY1PY1(MVSIZ),PY1PY2(MVSIZ),PY1PY3(MVSIZ),PY3(MVSIZ),
67 . PY2PY2(MVSIZ),PY2PY3(MVSIZ),PY3PY3(MVSIZ),
68 . C1,C2,DHZ(MVSIZ),GM(MVSIZ),GF(MVSIZ),C3,
69 . gpx1px1,gpx1py1,gpx1py2,gpx1py3,bcxx1,bcxx2,bcxy3,bcxy1,
70 . bcxy2,bcyy1,bcyy2,bcyy3,bcyx1,bcyx2,bcyx3,
71 . px1dx,py1dy,py2dy,py3dy,fac,fxx,fyy,fxy,fxy2,h33(mvsiz),
72 . h11(mvsiz),h12(mvsiz),h22(mvsiz),h13(mvsiz),h23(mvsiz),
73 . dm5(mvsiz),dm6(mvsiz),df5(mvsiz),df6(mvsiz),dm5_2(mvsiz),
74 . dm6_2(mvsiz),df5_2(mvsiz),df6_2(mvsiz),
75 . dmf(3,3,mvsiz),mfij(2,2,mvsiz)
76C-----------Attention Matrice sym Kii ne calcul que la moitie---------72
77 DO ep=jft,jlt
78 c2=vol(ep)
79 c1=thk2(ep)*c2
80 dm(1,1,ep)=hm(ep,1)*c2
81 dm(2,2,ep)=hm(ep,2)*c2
82 dm(1,2,ep)=hm(ep,3)*c2
83 dm(2,1,ep)=dm(1,2,ep)
84 gm(ep) =hm(ep,4)*c2
85 df(1,1,ep)=hf(ep,1)*c1
86 df(2,2,ep)=hf(ep,2)*c1
87 df(1,2,ep)=hf(ep,3)*c1
88 df(2,1,ep)=df(1,2,ep)
89 gf(ep) =hf(ep,4)*c1
90 dhz(ep)=hz(ep)*c1
91 dcx(ep)=hc(ep,1)*c2
92 dcy(ep)=hc(ep,2)*c2
93 ENDDO
94 DO ep=jft,jlt
95C-------PX2=-PX1 ,PX3=0---------------
96 py3(ep)= -py1(ep)-py2(ep)
97 px1px1(ep)=px1(ep)*px1(ep)
98 px1py1(ep)=px1(ep)*py1(ep)
99 px1py2(ep)=px1(ep)*py2(ep)
100 px1py3(ep)=-px1py1(ep)-px1py2(ep)
101 py1py1(ep)=py1(ep)*py1(ep)
102 py1py2(ep)=py1(ep)*py2(ep)
103 py1py3(ep)=-py1py1(ep)-py1py2(ep)
104 py2py2(ep)=py2(ep)*py2(ep)
105 py2py3(ep)=py2(ep)*py3(ep)
106 py3py3(ep)=py3(ep)*py3(ep)
107 ENDDO
108C------Membrainaires-----------
109 DO ep=jft,jlt
110 k11(1,1,ep)=px1px1(ep)*dm(1,1,ep)
111 k11(1,2,ep)=px1py1(ep)*dm(1,2,ep)
112 k11(2,2,ep)=py1py1(ep)*dm(2,2,ep)
113 k22(1,1,ep)=k11(1,1,ep)
114 k22(1,2,ep)=-px1py2(ep)*dm(1,2,ep)
115 k22(2,2,ep)=py2py2(ep)*dm(2,2,ep)
116 k33(2,2,ep)=py3py3(ep)*dm(2,2,ep)
117C
118 k12(1,1,ep)=-k11(1,1,ep)
119 k12(1,2,ep)=-k22(1,2,ep)
120 k12(2,1,ep)=-k11(1,2,ep)
121 k12(2,2,ep)=py1py2(ep)*dm(2,2,ep)
122 k13(1,2,ep)=px1py3(ep)*dm(1,2,ep)
123 k13(2,2,ep)=py1py3(ep)*dm(2,2,ep)
124 k23(1,2,ep)=-k13(1,2,ep)
125 k23(2,2,ep)=py2py3(ep)*dm(2,2,ep)
126 ENDDO
127C------Flexion-----------
128 DO ep=jft,jlt
129 m11(2,2,ep)=px1px1(ep)*df(1,1,ep)
130 m11(1,2,ep)=-px1py1(ep)*df(1,2,ep)
131 m11(1,1,ep)=py1py1(ep)*df(2,2,ep)
132 m22(2,2,ep)=m11(2,2,ep)
133 m22(1,2,ep)=px1py2(ep)*df(1,2,ep)
134 m22(1,1,ep)=py2py2(ep)*df(2,2,ep)
135 m33(1,1,ep)=py3py3(ep)*df(2,2,ep)
136C
137 m12(2,2,ep)=-m11(2,2,ep)
138 m12(1,2,ep)=-m11(1,2,ep)
139 m12(2,1,ep)=-m22(1,2,ep)
140 m12(1,1,ep)=py1py2(ep)*df(2,2,ep)
141 m13(2,1,ep)=-px1py3(ep)*df(1,2,ep)
142 m13(1,1,ep)=py1py3(ep)*df(2,2,ep)
143 m23(2,1,ep)=-m13(2,1,ep)
144 m23(1,1,ep)=py2py3(ep)*df(2,2,ep)
145 ENDDO
146C------G terms-----------
147 DO ep=jft,jlt
148 gpx1px1=px1px1(ep)*gm(ep)
149 gpx1py1=px1py1(ep)*gm(ep)
150 gpx1py2=px1py2(ep)*gm(ep)
151 gpx1py3=-gpx1py1-gpx1py2
152 k11(1,1,ep)=k11(1,1,ep)+py1py1(ep)*gm(ep)
153 k11(1,2,ep)=k11(1,2,ep)+gpx1py1
154 k11(2,2,ep)=k11(2,2,ep)+gpx1px1
155 k22(1,1,ep)=k22(1,1,ep)+py2py2(ep)*gm(ep)
156 k22(1,2,ep)=k22(1,2,ep)-gpx1py2
157 k22(2,2,ep)=k22(2,2,ep)+gpx1px1
158 k33(1,1,ep)=k33(1,1,ep)+py3py3(ep)*gm(ep)
159C
160 k12(1,1,ep)=k12(1,1,ep)+py1py2(ep)*gm(ep)
161 k12(1,2,ep)=k12(1,2,ep)-gpx1py1
162 k12(2,1,ep)=k12(2,1,ep)+gpx1py2
163 k12(2,2,ep)=k12(2,2,ep)-gpx1px1
164 k13(1,1,ep)=k13(1,1,ep)+py1py3(ep)*gm(ep)
165 k13(2,1,ep)=k13(2,1,ep)+gpx1py3
166 k23(1,1,ep)=k23(1,1,ep)+py2py3(ep)*gm(ep)
167 k23(2,1,ep)=k23(2,1,ep)-gpx1py3
168 ENDDO
169C
170 DO ep=jft,jlt
171 gpx1px1=px1px1(ep)*gf(ep)
172 gpx1py1=px1py1(ep)*gf(ep)
173 gpx1py2=px1py2(ep)*gf(ep)
174 gpx1py3=-gpx1py1-gpx1py2
175 m11(2,2,ep)=m11(2,2,ep)+py1py1(ep)*gf(ep)
176 m11(1,2,ep)=m11(1,2,ep)-gpx1py1
177 m11(1,1,ep)=m11(1,1,ep)+gpx1px1
178 m22(2,2,ep)=m22(2,2,ep)+py2py2(ep)*gf(ep)
179 m22(1,2,ep)=m22(1,2,ep)+gpx1py2
180 m22(1,1,ep)=m22(1,1,ep)+gpx1px1
181 m33(2,2,ep)=m33(2,2,ep)+py3py3(ep)*gf(ep)
182C
183 m12(2,2,ep)=m12(2,2,ep)+py1py2(ep)*gf(ep)
184 m12(1,2,ep)=m12(1,2,ep)-gpx1py2
185 m12(2,1,ep)=m12(2,1,ep)+gpx1py1
186 m12(1,1,ep)=m12(1,1,ep)-gpx1px1
187 m13(1,2,ep)=m13(1,2,ep)-gpx1py3
188 m13(2,2,ep)=m13(2,2,ep)+py1py3(ep)*gf(ep)
189 m23(1,2,ep)=m23(1,2,ep)+gpx1py3
190 m23(2,2,ep)=m23(2,2,ep)+py2py3(ep)*gf(ep)
191 ENDDO
192 IF (iorth==1) THEN
193#include "vectorize.inc"
194 DO ep=jft,jlt
195 c2=vol(ep)
196 c1=thk2(ep)*c2
197 dm5(ep)=hmor(ep,1)*c2
198 dm6(ep)=hmor(ep,2)*c2
199 df5(ep)=hfor(ep,1)*c1
200 df6(ep)=hfor(ep,2)*c1
201 ENDDO
202 DO ep=jft,jlt
203 dm5_2(ep)=two*dm5(ep)
204 dm6_2(ep)=two*dm6(ep)
205 df5_2(ep)=two*df5(ep)
206 df6_2(ep)=two*df6(ep)
207 ENDDO
208C----Px2=-Px1,--Px3=0---
209 DO ep=jft,jlt
210 k11(1,1,ep)=k11(1,1,ep)+px1py1(ep)*dm5_2(ep)
211 k11(1,2,ep)=k11(1,2,ep)+
212 1 px1px1(ep)*dm5(ep)+py1py1(ep)*dm6(ep)
213 k11(2,2,ep)=k11(2,2,ep)+px1py1(ep)*dm6_2(ep)
214 k12(1,1,ep)=k12(1,1,ep)+(px1py2(ep)-px1py1(ep))*dm5(ep)
215 c1 = -px1px1(ep)*dm5(ep)+py1py2(ep)*dm6(ep)
216 k12(1,2,ep)=k12(1,2,ep)+c1
217 k12(2,1,ep)=k12(2,1,ep)+c1
218 k12(2,2,ep)=k12(2,2,ep)+(px1py2(ep)-px1py1(ep))*dm6(ep)
219 k13(1,1,ep)=k13(1,1,ep)+px1py3(ep)*dm5(ep)
220 c1 = py1py3(ep)*dm6(ep)
221 k13(1,2,ep)=k13(1,2,ep)+c1
222 k13(2,1,ep)=k13(2,1,ep)+c1
223 k13(2,2,ep)=k13(2,2,ep)+px1py3(ep)*dm6(ep)
224C
225 k22(1,1,ep)=k22(1,1,ep)-px1py2(ep)*dm5_2(ep)
226 k22(1,2,ep)=k22(1,2,ep)+
227 1 px1px1(ep)*dm5(ep)+py2py2(ep)*dm6(ep)
228 k22(2,2,ep)=k22(2,2,ep)-px1py2(ep)*dm6_2(ep)
229 k23(1,1,ep)=k23(1,1,ep)-px1py3(ep)*dm5(ep)
230C
231 c1 = py2py3(ep)*dm6(ep)
232 k23(1,2,ep)=k23(1,2,ep)+c1
233 k23(2,1,ep)=k23(2,1,ep)+c1
234 k23(2,2,ep)=k23(2,2,ep)-px1py3(ep)*dm6(ep)
235C
236 k33(1,2,ep)=k33(1,2,ep)+py3py3(ep)*dm6(ep)
237C
238 m11(1,1,ep)=m11(1,1,ep)+px1py1(ep)*df6_2(ep)
239 m11(1,2,ep)=m11(1,2,ep)-
240 1 px1px1(ep)*df5(ep)-py1py1(ep)*df6(ep)
241 m11(2,2,ep)=m11(2,2,ep)+px1py1(ep)*df5_2(ep)
242 m12(1,1,ep)=m12(1,1,ep)+(px1py2(ep)-px1py1(ep))*df6(ep)
243 c2 = -px1px1(ep)*df5(ep)+py1py2(ep)*df6(ep)
244 m12(1,2,ep)=m12(1,2,ep)-c2
245 m12(2,1,ep)=m12(2,1,ep)-c2
246 m12(2,2,ep)=m12(2,2,ep)+(px1py2(ep)-px1py1(ep))*df5(ep)
247 m13(1,1,ep)=m13(1,1,ep)+px1py3(ep)*df6(ep)
248 m13(1,2,ep)=m13(1,2,ep)-py1py3(ep)*df6(ep)
249 m13(2,1,ep)=m13(2,1,ep)-py1py3(ep)*df6(ep)
250 m13(2,2,ep)=m13(2,2,ep)+px1py3(ep)*df5(ep)
251 m22(1,1,ep)=m22(1,1,ep)-px1py2(ep)*df6_2(ep)
252 m22(1,2,ep)=m22(1,2,ep)-
253 1 px1px1(ep)*df5(ep)-py2py2(ep)*df6(ep)
254 m22(2,2,ep)=m22(2,2,ep)-px1py2(ep)*df5_2(ep)
255 m23(1,1,ep)=m23(1,1,ep)-px1py3(ep)*df6(ep)
256 m23(1,2,ep)=m23(1,2,ep)-py2py3(ep)*df6(ep)
257 m23(2,1,ep)=m23(2,1,ep)-py2py3(ep)*df6(ep)
258 m23(2,2,ep)=m23(2,2,ep)-px1py3(ep)*df5(ep)
259 m33(1,2,ep)=m33(1,2,ep)-py3py3(ep)*df6(ep)
260 ENDDO
261 DO ep=jft,jlt
262 c2=vol(ep)*thk0(ep)
263 dmf(1,1,ep)=hmfor(ep,1)*c2
264 dmf(2,2,ep)=hmfor(ep,2)*c2
265 dmf(1,2,ep)=hmfor(ep,3)*c2
266 dmf(1,3,ep)=hmfor(ep,5)*c2
267 dmf(2,3,ep)=hmfor(ep,6)*c2
268 dmf(2,1,ep)=dmf(1,2,ep)
269 dmf(3,1,ep)=dmf(1,3,ep)
270 dmf(3,2,ep)=dmf(2,3,ep)
271 dmf(3,3,ep)=hmfor(ep,4)*c2
272 ENDDO
273C ----add mem-bending coupling terms of orthotrope--coplanar first---
274 DO ep=jft,jlt
275 mf11(1,1,ep)=
276 1 -dmf(1,2,ep)*px1py1(ep)-dmf(1,3,ep)*px1px1(ep)
277 2 -dmf(2,3,ep)*py1py1(ep)-dmf(3,3,ep)*px1py1(ep)
278 mf11(1,2,ep)=
279 1 dmf(1,1,ep)*px1px1(ep)+dmf(1,3,ep)*px1py1(ep)
280 2 +dmf(1,3,ep)*px1py1(ep)+dmf(3,3,ep)*py1py1(ep)
281 mf11(2,1,ep)=
282 1 -dmf(2,2,ep)*py1py1(ep)-dmf(2,3,ep)*px1py1(ep)
283 2 -dmf(2,3,ep)*px1py1(ep)-dmf(3,3,ep)*px1px1(ep)
284 mf11(2,2,ep)=
285 1 dmf(1,2,ep)*px1py1(ep)+dmf(2,3,ep)*py1py1(ep)
286 2 +dmf(1,3,ep)*px1px1(ep)+dmf(3,3,ep)*px1py1(ep)
287 ENDDO
288 DO ep=jft,jlt
289 mf12(1,1,ep)=
290 1 -dmf(1,2,ep)*px1py2(ep)+dmf(1,3,ep)*px1px1(ep)
291C--------------------------------------------------------YIXJ-----
292 2 -dmf(2,3,ep)*py1py2(ep)+dmf(3,3,ep)*px1py1(ep)
293 mf12(1,2,ep)=
294 1 -dmf(1,1,ep)*px1px1(ep)+dmf(1,3,ep)*px1py2(ep)
295C-------------------------YIXJ-----
296 2 -dmf(1,3,ep)*px1py1(ep)+dmf(3,3,ep)*py1py2(ep)
297 mf12(2,1,ep)=
298C---------------------------------------------------------YIXJ-----
299 1 -dmf(2,2,ep)*py1py2(ep)+dmf(2,3,ep)*px1py1(ep)
300 2 -dmf(2,3,ep)*px1py2(ep)+dmf(3,3,ep)*px1px1(ep)
301 mf12(2,2,ep)=
302C---------------------------------YIXJ-----
303 1 -dmf(1,2,ep)*px1py1(ep)+dmf(2,3,ep)*py1py2(ep)
304 2 -dmf(1,3,ep)*px1px1(ep)+dmf(3,3,ep)*px1py2(ep)
305 ENDDO
306 DO ep=jft,jlt
307 mf22(1,1,ep)=
308 1 dmf(1,2,ep)*px1py2(ep)-dmf(1,3,ep)*px1px1(ep)
309C--------------------------------------------------------YIXJ-----
310 2 -dmf(2,3,ep)*py2py2(ep)+dmf(3,3,ep)*px1py2(ep)
311 mf22(1,2,ep)=
312 1 dmf(1,1,ep)*px1px1(ep)-dmf(1,3,ep)*px1py2(ep)
313C-------------------------YIXJ-----
314 2 -dmf(1,3,ep)*px1py2(ep)+dmf(3,3,ep)*py2py2(ep)
315 mf22(2,1,ep)=
316C---------------------------------------------------------YIXJ-----
317 1 -dmf(2,2,ep)*py2py2(ep)+dmf(2,3,ep)*px1py2(ep)
318 2 +dmf(2,3,ep)*px1py2(ep)-dmf(3,3,ep)*px1px1(ep)
319 mf22(2,2,ep)=
320C---------------------------------YIXJ-----
321 1 -dmf(1,2,ep)*px1py2(ep)+dmf(2,3,ep)*py2py2(ep)
322 2 +dmf(1,3,ep)*px1px1(ep)-dmf(3,3,ep)*px1py2(ep)
323 ENDDO
324 DO ep=jft,jlt
325 mf23(1,1,ep)=
326 1 dmf(1,2,ep)*px1py3(ep)
327C--------------------------------------------------------YIXJ-----
328 2 -dmf(2,3,ep)*py2py3(ep)
329 mf23(1,2,ep)=
330 1 -dmf(1,3,ep)*px1py3(ep)
331C-------------------------YIXJ-----
332 2 +dmf(3,3,ep)*py2py3(ep)
333 mf23(2,1,ep)=
334C---------------------------------------------------------YIXJ-----
335 1 -dmf(2,2,ep)*py2py3(ep)
336 2 +dmf(2,3,ep)*px1py3(ep)
337 mf23(2,2,ep)=
338C---------------------------------YIXJ-----
339 1 dmf(2,3,ep)*py2py3(ep)
340 2 -dmf(3,3,ep)*px1py3(ep)
341 ENDDO
342 DO ep=jft,jlt
343 mf13(1,1,ep)=
344 1 -dmf(1,2,ep)*px1py3(ep)
345C--------------------------------------------------------YIXJ-----
346 2 -dmf(2,3,ep)*py1py3(ep)
347 mf13(1,2,ep)=
348 1 dmf(1,3,ep)*px1py3(ep)
349C-------------------------YIXJ-----
350 2 +dmf(3,3,ep)*py1py3(ep)
351 mf13(2,1,ep)=
352C---------------------------------------------------------YIXJ-----
353 1 -dmf(2,2,ep)*py1py3(ep)
354 2 -dmf(2,3,ep)*px1py3(ep)
355 mf13(2,2,ep)=
356C---------------------------------YIXJ-----
357 1 dmf(2,3,ep)*py1py3(ep)
358 2 +dmf(3,3,ep)*px1py3(ep)
359 ENDDO
360 DO ep=jft,jlt
361 mf33(1,1,ep)=
362 2 -dmf(2,3,ep)*py3py3(ep)
363 mf33(1,2,ep)=
364 2 dmf(3,3,ep)*py3py3(ep)
365 mf33(2,1,ep)=
366 1 -dmf(2,2,ep)*py3py3(ep)
367 mf33(2,2,ep)=
368 1 dmf(2,3,ep)*py3py3(ep)
369 ENDDO
370C---------FMIJ(i,j)=MFJI(j,i)-----------
371 DO ep=jft,jlt
372 fm12(1,1,ep)=
373 1 dmf(1,2,ep)*px1py1(ep)+dmf(1,3,ep)*px1px1(ep)
374C--------------------------------------------------------YIXJ-----
375 2 -dmf(2,3,ep)*py1py2(ep)-dmf(3,3,ep)*px1py2(ep)
376 fm12(2,1,ep)=
377 1 -dmf(1,1,ep)*px1px1(ep)-dmf(1,3,ep)*px1py1(ep)
378c 1 -DMF(1,1,EP)*PX1PX1(EP)+DMF(1,3,EP)*PX1PY2(EP)
379C-------------------------YIXJ-----
380 2 +dmf(1,3,ep)*px1py2(ep)+dmf(3,3,ep)*py1py2(ep)
381 fm12(1,2,ep)=
382C---------------------------------------------------------YIXJ-----
383 1 -dmf(2,2,ep)*py1py2(ep)-dmf(2,3,ep)*px1py2(ep)
384 2 +dmf(2,3,ep)*px1py1(ep)+dmf(3,3,ep)*px1px1(ep)
385 fm12(2,2,ep)=
386C---------------------------------YIXJ-----
387 1 dmf(1,2,ep)*px1py2(ep)+dmf(2,3,ep)*py1py2(ep)
388 2 -dmf(1,3,ep)*px1px1(ep)-dmf(3,3,ep)*px1py1(ep)
389 ENDDO
390 DO ep=jft,jlt
391 fm13(1,1,ep)=
392c 1 -DMF(1,2,EP)*PX3PY1(EP)-DMF(1,3,EP)*PX3PX1(EP)
393C--------------------------------------------------------YIXJ-----
394 2 -dmf(2,3,ep)*py1py3(ep)-dmf(3,3,ep)*px1py3(ep)
395 fm13(2,1,ep)=
396c 1 +DMF(1,1,EP)*PX3PX1(EP)+DMF(1,3,EP)*PX3PY1(EP)
397C-------------------------YIXJ-----
398 2 dmf(1,3,ep)*px1py3(ep)+dmf(3,3,ep)*py1py3(ep)
399 fm13(1,2,ep)=
400C---------------------------------------------------------YIXJ-----
401 1 -dmf(2,2,ep)*py1py3(ep)-dmf(2,3,ep)*px1py3(ep)
402c 2 -DMF(2,3,EP)*PX3PY1(EP)-DMF(3,3,EP)*PX3PX1(EP)
403 fm13(2,2,ep)=
404C---------------------------------YIXJ-----
405 1 dmf(1,2,ep)*px1py3(ep)+dmf(2,3,ep)*py1py3(ep)
406c 2 +DMF(1,3,EP)*PX3PX1(EP)+DMF(3,3,EP)*PX3PY1(EP)
407 ENDDO
408 DO ep=jft,jlt
409 fm23(1,1,ep)=
410c 1 -DMF(1,2,EP)*PX3PY2(EP)-DMF(1,3,EP)*PX2PX3(EP)
411C--------------------------------------------------------YIXJ-----
412 2 -dmf(2,3,ep)*py2py3(ep)+dmf(3,3,ep)*px1py3(ep)
413 fm23(2,1,ep)=
414c 1 DMF(1,1,EP)*PX2PX3(EP)+DMF(1,3,EP)*PX3PY2(EP)
415C-------------------------YIXJ-----
416 2 -dmf(1,3,ep)*px1py3(ep)+dmf(3,3,ep)*py2py3(ep)
417 fm23(1,2,ep)=
418C---------------------------------------------------------YIXJ-----
419 1 -dmf(2,2,ep)*py2py3(ep)+dmf(2,3,ep)*px1py3(ep)
420c 2 -DMF(2,3,EP)*PX3PY2(EP)-DMF(3,3,EP)*PX2PX3(EP)
421 fm23(2,2,ep)=
422C---------------------------------YIXJ-----
423 1 -dmf(1,2,ep)*px1py3(ep)+dmf(2,3,ep)*py2py3(ep)
424c 2 +DMF(1,3,EP)*PX2PX3(EP)+DMF(3,3,EP)*PX3PY2(EP)
425 ENDDO
426 ENDIF
427C-------------C.T.---------------
428 DO ep=jft,jlt
429 gpx1px1=px1px1(ep)*dcx(ep)
430 k11(3,3,ep)=gpx1px1+py1py1(ep)*dcy(ep)
431 k22(3,3,ep)=gpx1px1+py2py2(ep)*dcy(ep)
432 k33(3,3,ep)=py3py3(ep)*dcy(ep)
433 k12(3,3,ep)=-gpx1px1+py1py2(ep)*dcy(ep)
434 k13(3,3,ep)=py1py3(ep)*dcy(ep)
435 k23(3,3,ep)=py2py3(ep)*dcy(ep)
436 ENDDO
437 DO ep=jft,jlt
438 fac = third*area(ep)
439 bcxx1=-fac*px1px1(ep)
440 bcxx2= -bcxx1
441 bcxy3= fac*(px1py1(ep)+px1py2(ep))
442 bcxy1= bcxy3+bcxy3+fac*px1py2(ep)
443 bcxy2= bcxy3+bcxy3+fac*px1py1(ep)
444 bcyy1= fac*(py1py1(ep)+py1py2(ep)+py1py2(ep))
445 bcyy2= -fac*(py2py2(ep)+py1py2(ep)+py1py2(ep))
446 bcyy3= -bcyy1-bcyy2
447 bcyx1= -bcxy3-fac*px1py1(ep)
448 bcyx2= -bcxy3-fac*px1py2(ep)
449 bcyx3= bcyx1+bcyx2
450 m11(1,1,ep)=m11(1,1,ep)+bcxx1*bcxx1*dcx(ep)+bcyx1*bcyx1*dcy(ep)
451 m11(1,2,ep)=m11(1,2,ep)+bcxx1*bcxy1*dcx(ep)+bcyx1*bcyy1*dcy(ep)
452 m11(2,2,ep)=m11(2,2,ep)+bcxy1*bcxy1*dcx(ep)+bcyy1*bcyy1*dcy(ep)
453 m22(1,1,ep)=m22(1,1,ep)+bcxx2*bcxx2*dcx(ep)+bcyx2*bcyx2*dcy(ep)
454 m22(1,2,ep)=m22(1,2,ep)+bcxx2*bcxy2*dcx(ep)+bcyx2*bcyy2*dcy(ep)
455 m22(2,2,ep)=m22(2,2,ep)+bcxy2*bcxy2*dcx(ep)+bcyy2*bcyy2*dcy(ep)
456 m33(1,1,ep)=m33(1,1,ep)+bcyx3*bcyx3*dcy(ep)
457 m33(1,2,ep)=m33(1,2,ep)+bcyx3*bcyy3*dcy(ep)
458 m33(2,2,ep)=m33(2,2,ep)+bcxy3*bcxy3*dcx(ep)+bcyy3*bcyy3*dcy(ep)
459 m12(1,1,ep)=m12(1,1,ep)+bcxx1*bcxx2*dcx(ep)+bcyx1*bcyx2*dcy(ep)
460 m12(1,2,ep)=m12(1,2,ep)+bcxx1*bcxy2*dcx(ep)+bcyx1*bcyy2*dcy(ep)
461 m12(2,1,ep)=m12(2,1,ep)+bcxx2*bcxy1*dcx(ep)+bcyx2*bcyy1*dcy(ep)
462 m12(2,2,ep)=m12(2,2,ep)+bcxy1*bcxy2*dcx(ep)+bcyy1*bcyy2*dcy(ep)
463 m13(1,1,ep)=m13(1,1,ep)+bcyx1*bcyx3*dcy(ep)
464 m13(1,2,ep)=m13(1,2,ep)+bcxx1*bcxy3*dcx(ep)+bcyx1*bcyy3*dcy(ep)
465 m13(2,1,ep)=m13(2,1,ep)+bcyx3*bcyy1*dcy(ep)
466 m13(2,2,ep)=m13(2,2,ep)+bcxy1*bcxy3*dcx(ep)+bcyy1*bcyy3*dcy(ep)
467 m23(1,1,ep)=m23(1,1,ep)+bcyx2*bcyx3*dcy(ep)
468 m23(1,2,ep)=m23(1,2,ep)+bcxx2*bcxy3*dcx(ep)+bcyx2*bcyy3*dcy(ep)
469 m23(2,1,ep)=m23(2,1,ep)+bcyx3*bcyy2*dcy(ep)
470 m23(2,2,ep)=m23(2,2,ep)+bcxy2*bcxy3*dcx(ep)+bcyy2*bcyy3*dcy(ep)
471C---
472 px1dx=px1(ep)*dcx(ep)
473 py1dy=py1(ep)*dcy(ep)
474 py2dy=py2(ep)*dcy(ep)
475 py3dy=py3(ep)*dcy(ep)
476 mf11(3,1,ep)=bcxx1*px1dx+bcyx1*py1dy
477 mf11(3,2,ep)=bcxy1*px1dx+bcyy1*py1dy
478 mf22(3,1,ep)=-bcxx2*px1dx+bcyx2*py2dy
479 mf22(3,2,ep)=-bcxy2*px1dx+bcyy2*py2dy
480 mf33(3,1,ep)=bcyx3*py3dy
481 mf33(3,2,ep)=bcyy3*py3dy
482 mf12(3,1,ep)=bcxx2*px1dx+bcyx2*py1dy
483 mf12(3,2,ep)=bcxy2*px1dx+bcyy2*py1dy
484 mf13(3,1,ep)=bcyx3*py1dy
485 mf13(3,2,ep)=bcxy3*px1dx+bcyy3*py1dy
486 mf23(3,1,ep)=bcyx3*py2dy
487 mf23(3,2,ep)=-bcxy3*px1dx+bcyy3*py2dy
488 fm12(1,3,ep)=-bcxx1*px1dx+bcyx1*py2dy
489 fm12(2,3,ep)=-bcxy1*px1dx+bcyy1*py2dy
490 fm13(1,3,ep)=bcyx1*py3dy
491 fm13(2,3,ep)=bcyy1*py3dy
492 fm23(1,3,ep)=bcyx2*py3dy
493 fm23(2,3,ep)=bcyy2*py3dy
494 ENDDO
495C
496C------ Mzz pour tous-------
497 DO ep=jft,jlt
498 m11(3,3,ep)=m11(3,3,ep)+dhz(ep)*(px1px1(ep)+py1py1(ep))
499 m12(3,3,ep)=m12(3,3,ep)+dhz(ep)*(py1py2(ep)-px1px1(ep))
500 m13(3,3,ep)=m13(3,3,ep)+dhz(ep)*py1py3(ep)
501 m22(3,3,ep)=m22(3,3,ep)+dhz(ep)*(px1px1(ep)+py2py2(ep))
502 m23(3,3,ep)=m23(3,3,ep)+dhz(ep)*py2py3(ep)
503 m33(3,3,ep)=m33(3,3,ep)+dhz(ep)*py3py3(ep)
504 ENDDO
505C------ renforce positive ----------
506 IF (neig==0) THEN
507 DO ep=jft,jlt
508 c1 = em8*m11(3,3,ep)
509 c2 = em8*m22(3,3,ep)
510 c3 = em8*m33(3,3,ep)
511 m11(3,3,ep)=m11(3,3,ep)+c1
512 m22(3,3,ep)=m22(3,3,ep)+c2
513 m33(3,3,ep)=m33(3,3,ep)+c3
514 ENDDO
515 ENDIF
516C
517 IF (ikgeo == 1 ) THEN
518C---------membrane--only-----
519 DO ep=jft,jlt
520 c2=vol(ep)
521 fxx=for(ep,1)*c2
522 fyy=for(ep,2)*c2
523 fxy=for(ep,3)*c2
524 fxy2=two*fxy
525 h11(ep)=fxx*px1px1(ep)+fyy*py1py1(ep)+fxy2*px1py1(ep)
526 h12(ep)=-fxx*px1px1(ep)+fyy*py1py2(ep)
527 . +fxy*(px1py2(ep)-px1py1(ep))
528 h22(ep)=fxx*px1px1(ep)+fyy*py2py2(ep)-fxy2*px1py2(ep)
529 h13(ep)=fyy*py1py3(ep)+fxy*px1py3(ep)
530 h23(ep)=fyy*py2py3(ep)-fxy*px1py3(ep)
531 h33(ep)=fyy*py3py3(ep)
532 ENDDO
533 DO i=1,3
534 DO ep=jft,jlt
535 k11(i,i,ep) = k11(i,i,ep)+h11(ep)
536 k12(i,i,ep) = k12(i,i,ep)+h12(ep)
537 k22(i,i,ep) = k22(i,i,ep)+h22(ep)
538 k13(i,i,ep) = k13(i,i,ep)+h13(ep)
539 k23(i,i,ep) = k23(i,i,ep)+h23(ep)
540 k33(i,i,ep) = k33(i,i,ep)+h33(ep)
541 ENDDO
542 ENDDO
543C------ renforce positive ----------
544 IF (neig==0.AND.idril==0.AND.iorth>0) THEN
545 DO ep=jft,jlt
546 c1 = min(h11(ep),h22(ep),h33(ep))
547 IF (c1 < zero) THEN
548 c2 =min(m11(3,3,ep),m22(3,3,ep),m33(3,3,ep))
549 c1 = min(-c1,ten*c2)
550 m11(3,3,ep)=m11(3,3,ep) + c1
551 m22(3,3,ep)=m22(3,3,ep) + c1
552 m33(3,3,ep)=m33(3,3,ep) + c1
553 END IF
554 ENDDO
555 ENDIF
556 END IF !IF (IKGEO ==1)
557C
558 RETURN
559 END
560!||====================================================================
561!|| c3lkerz3 ../engine/source/elements/sh3n/coque3n/c3lke3.F
562!||--- called by ------------------------------------------------------
563!|| c3ke3 ../engine/source/elements/sh3n/coque3n/c3ke3.F
564!||====================================================================
565 SUBROUTINE c3lkerz3(JFT,JLT,HM,HZ,
566 1 PX1,PY1,PY2,VOL,AREA,
567 2 K11,K12,K13,K22,K23,K33,
568 3 M11,M12,M13,M22,M23,M33,
569 4 MF11,MF12,MF13,MF22,MF23,MF33,
570 5 FM12,FM13,FM23,IORTH,HMOR,
571 6 BM0RZ,B0RZ,BKRZ,BERZ,THK0,HMFOR )
572C--------------------------------------------------------------------------------------------------
573C-----------------------------------------------
574C I M P L I C I T T Y P E S
575C-----------------------------------------------
576#include "implicit_f.inc"
577#include "mvsiz_p.inc"
578C-----------------------------------------------
579C D U M M Y A R G U M E N T S
580C-----------------------------------------------
581C---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8
582 INTEGER JFT,JLT,IORTH
583 MY_REAL
584 . VOL(*), PX1(*),PY1(*),PY2(*),AREA(*),THK0(*),
585 . HM(MVSIZ,4),HZ(*),HMOR(MVSIZ,2),
586 . K11(3,3,*),K12(3,3,*),K13(3,3,*),
587 . K22(3,3,*),K23(3,3,*),K33(3,3,*),
588 . M11(3,3,*),M12(3,3,*),M13(3,3,*),
589 . M22(3,3,*),M23(3,3,*),M33(3,3,*),
590 . MF11(3,3,*),MF12(3,3,*),MF13(3,3,*),
591 . MF22(3,3,*),MF23(3,3,*),MF33(3,3,*),
592 . FM12(3,3,*),FM13(3,3,*),FM23(3,3,*),
593 . BM0RZ(MVSIZ,3,2),B0RZ(MVSIZ,3),BKRZ(MVSIZ,2),BERZ(MVSIZ,2),HMFOR(MVSIZ,6)
594C-----------------------------------------------
595c FUNCTION: elementary sub-stiffness matrix relative to the drilling dof for Tria
596c
597c Note:
598c ARGUMENTS: (I: input, O: output, IO: input * output, W: workspace)
599c
600c TYPE NAME FUNCTION
601c I JFT,JLT - element id limit
602c I HM(4,NEL) - membrane stiffness modulus (plane stress)
603c I HZ(NEL) -drilling dof modulus
604C I PX1,PY1,PX2=-PX1,PY2(NEL): standard [B] of Tria
605c I VOL,AREA - element volume,AREA
606c IO Kij,FMij,MFij,Mij(3,3,NEL) sub-stiffness matrix
607C---------------|[KIJ][MFIJ]|----
608C-----KE(6x6)= | |
609C---------------|[FMIJ]{MIJ]|----
610c I IORTH - flag for orthotropic mat
611c I HMOR(2,NEL) - supplementary membrane stiffness modulus for orth
612c O BM0RZ(I,J,NEL) - constant terms of derivations for membrane
613C I=1:A*Nx,x;I=2:A*Ny,y;I=3:A*(Nx,y+Ny,x); J=1,2(node)
614C only store J=1,2 as f(j=3)=-f(j=1)-f(j=2)
615C O B0RZ(J,NEL) A*(-Nx,y+Ny,x -2Ni) for asymmetric rotation
616c O BKRZ(J,NEL) - Ksi terms of derivation : A*(-Nx,y+Ny,x -2Ni)
617c O BERZ(J,NEL) - Eta terms of derivation : A*(-Nx,y+Ny,x -2Ni)
618C-----------------------------------------------
619C L O C A L V A R I A B L E S
620C-----------------------------------------------
621 INTEGER EP,I,J,NG,NPG
622 MY_REAL
623 . DM(3,3,MVSIZ),
624 . PX1PX1(MVSIZ),PX1PY1(MVSIZ),PX1PY2(MVSIZ),PX1PY3(MVSIZ),
625 . PY1PY1(MVSIZ),PY1PY2(MVSIZ),PY1PY3(MVSIZ),PY3(MVSIZ),
626 . PY2PY2(MVSIZ),PY2PY3(MVSIZ),PY3PY3(MVSIZ),
627 . C1,C2,DHZ(MVSIZ),C3,
628 . KZ11(2,2,MVSIZ),KZ12(2,2,MVSIZ),KZ22(2,2,MVSIZ),
629 . KZ13(2,2,MVSIZ),KZ23(2,2,MVSIZ),KZ33(2,2,MVSIZ)
630 PARAMETER (NPG = 3)
631 my_real
632 . dprz(3,mvsiz),a_hammer(npg,2),bn1rz,bn2rz,bn3rz,dz12(mvsiz),
633 . prz(3,mvsiz),prx(3,mvsiz),pry(3,mvsiz),prxy(3,mvsiz),
634 . cbrx(3,mvsiz),cbry(3,mvsiz),cbrz(3,mvsiz),a_i(mvsiz),
635 . dmf(3,3,mvsiz)
636 DATA a_hammer /
637 1 0.166666666666667,0.666666666666667,0.166666666666667,
638 2 0.166666666666667,0.166666666666667,0.666666666666667/
639C---------------|[KIJ][MFIJ]|----
640C-----KE(6x6)= | |
641C---------------|[FMIJ]{MIJ]|----
642C-----------Attention Matrice sym Kii ne calcul que la moitie---------72
643C--------Bi: Bi*A--- DM: DM/A
644 DO ep=jft,jlt
645 a_i(ep)=one/max(area(ep),em20)
646 ENDDO
647 DO ep=jft,jlt
648 c2=vol(ep)
649 dm(1,1,ep)=hm(ep,1)*c2
650 dm(2,2,ep)=hm(ep,2)*c2
651 dm(1,2,ep)=hm(ep,3)*c2
652 dm(2,1,ep)=dm(1,2,ep)
653 dm(3,3,ep)=hm(ep,4)*c2
654 dhz(ep)=hz(ep)*fourth*c2
655 dz12(ep)=dhz(ep)*third
656 dm(1,3,ep)=zero
657 dm(2,3,ep)=zero
658 ENDDO
659 IF (iorth==1) THEN
660 DO ep=jft,jlt
661 c2=vol(ep)
662 dm(1,3,ep)=hmor(ep,1)*c2
663 dm(2,3,ep)=hmor(ep,2)*c2
664 ENDDO
665 ENDIF
666 DO ep=jft,jlt
667 dm(3,1,ep)=dm(1,3,ep)
668 dm(3,2,ep)=dm(2,3,ep)
669 ENDDO
670C------ PX2=-PX1 , PX3=0-------------
671 DO ep=jft,jlt
672 py3(ep)= -py1(ep)-py2(ep)
673 px1px1(ep)=px1(ep)*px1(ep)
674 px1py1(ep)=px1(ep)*py1(ep)
675 px1py2(ep)=px1(ep)*py2(ep)
676 px1py3(ep)=-px1py1(ep)-px1py2(ep)
677 py1py1(ep)=py1(ep)*py1(ep)
678 py1py2(ep)=py1(ep)*py2(ep)
679 py1py3(ep)=-py1py1(ep)-py1py2(ep)
680 py2py2(ep)=py2(ep)*py2(ep)
681 py2py3(ep)=py2(ep)*py3(ep)
682 py3py3(ep)=py3(ep)*py3(ep)
683 ENDDO
684C-------0.5*[-By Bx 0 0 0 BRZ]^tKG[-By Bx 0 0 0 BRZ]*0.5
685C------ BRZ = (-BRXY+BRYX)-0.5
686 DO ep=jft,jlt
687 kz11(1,1,ep)=py1py1(ep)*dhz(ep)
688 kz11(1,2,ep)=-px1py1(ep)*dhz(ep)
689 kz11(2,2,ep)=px1px1(ep)*dhz(ep)
690C
691 kz22(1,1,ep)=py2py2(ep)*dhz(ep)
692 kz22(1,2,ep)=px1py2(ep)*dhz(ep)
693 kz22(2,2,ep)=kz11(2,2,ep)
694C
695 kz33(1,1,ep)=py3py3(ep)*dhz(ep)
696 kz33(1,2,ep)=zero
697 kz33(2,2,ep)=zero
698C
699 kz12(1,1,ep)=py1py2(ep)*dhz(ep)
700C KZ12(1,2,EP)=PX1PY1(EP)*DHZ(EP)
701 kz12(1,2,ep)=-kz11(1,2,ep)
702 kz12(2,2,ep)=-kz11(2,2,ep)
703C KZ12(2,1,EP)=-PX1PY2(EP)*DHZ(EP)
704 kz12(2,1,ep)=-kz22(1,2,ep)
705C
706 kz13(1,1,ep)=py1py3(ep)*dhz(ep)
707 kz13(1,2,ep)=zero
708 kz13(2,2,ep)=zero
709 kz13(2,1,ep)=-px1py3(ep)*dhz(ep)
710C
711 kz23(1,1,ep)=py2py3(ep)*dhz(ep)
712 kz23(1,2,ep)=zero
713 kz23(2,2,ep)=zero
714C KZ23(2,1,EP)=PX1PY3(EP)*DHZ(EP)
715 kz23(2,1,ep)=-kz13(2,1,ep)
716 ENDDO
717C
718 DO i=1,2
719 DO j=i,2
720 DO ep=jft,jlt
721 k11(i,j,ep)=k11(i,j,ep)+kz11(i,j,ep)
722 k12(i,j,ep)=k12(i,j,ep)+kz12(i,j,ep)
723 k22(i,j,ep)=k22(i,j,ep)+kz22(i,j,ep)
724 ENDDO
725 ENDDO
726 ENDDO
727C
728 DO ep=jft,jlt
729 k12(2,1,ep)=k12(2,1,ep)+kz12(2,1,ep)
730 k13(1,1,ep)=k13(1,1,ep)+kz13(1,1,ep)
731 k13(2,1,ep)=k13(2,1,ep)+kz13(2,1,ep)
732 k23(1,1,ep)=k23(1,1,ep)+kz23(1,1,ep)
733 k23(2,1,ep)=k23(2,1,ep)+kz23(2,1,ep)
734 k33(1,1,ep)=k33(1,1,ep)+kz33(1,1,ep)
735 ENDDO
736C
737C------ Mzz reset zero----MF,FM are not initialized to zero------
738 DO ep=jft,jlt
739 m11(3,3,ep)=zero
740 m12(3,3,ep)=zero
741 m13(3,3,ep)=zero
742 m22(3,3,ep)=zero
743 m23(3,3,ep)=zero
744 m33(3,3,ep)=zero
745C
746 mf11(1,3,ep)=zero
747 mf11(2,3,ep)=zero
748 mf22(1,3,ep)=zero
749 mf22(2,3,ep)=zero
750 mf33(1,3,ep)=zero
751 mf33(2,3,ep)=zero
752 mf12(1,3,ep)=zero
753 mf12(2,3,ep)=zero
754 mf13(1,3,ep)=zero
755 mf13(2,3,ep)=zero
756 mf23(1,3,ep)=zero
757 mf23(2,3,ep)=zero
758C
759 fm12(3,1,ep)=zero
760 fm12(3,2,ep)=zero
761 fm13(3,1,ep)=zero
762 fm13(3,2,ep)=zero
763 fm23(3,1,ep)=zero
764 fm23(3,2,ep)=zero
765 ENDDO
766C
767 DO ng =1,npg
768 DO ep=jft,jlt
769 bn1rz=bkrz(ep,1)*a_hammer(ng,1)+berz(ep,1)*a_hammer(ng,2)
770 bn2rz=bkrz(ep,2)*a_hammer(ng,1)+berz(ep,2)*a_hammer(ng,2)
771 bn3rz=-bn1rz-bn2rz
772 prz(1,ep)=(b0rz(ep,1)+bn1rz)*a_i(ep)
773 prz(2,ep)=(b0rz(ep,2)+bn2rz)*a_i(ep)
774 prz(3,ep)=(b0rz(ep,3)+bn3rz)*a_i(ep)
775 dprz(1,ep)=prz(1,ep)*dz12(ep)
776 dprz(2,ep)=prz(2,ep)*dz12(ep)
777 dprz(3,ep)=prz(3,ep)*dz12(ep)
778 ENDDO
779C
780 DO ep=jft,jlt
781C------MFIJ(1,3)=-PYI*DPRZ(J); MFIJ(2,3)=PXI*DPRZ(J)
782C------FMIJ(3,1)=-PYJ*DPRZ(I); FMIJ(3,2)=PXJ*DPRZ(I)
783 mf11(1,3,ep)=mf11(1,3,ep)-py1(ep)*dprz(1,ep)
784 mf11(2,3,ep)=mf11(2,3,ep)+px1(ep)*dprz(1,ep)
785 mf22(1,3,ep)=mf22(1,3,ep)-py2(ep)*dprz(2,ep)
786 mf22(2,3,ep)=mf22(2,3,ep)-px1(ep)*dprz(2,ep)
787 mf33(1,3,ep)=mf33(1,3,ep)-py3(ep)*dprz(3,ep)
788C MF33(2,3,EP)=MF33(2,3,EP)+ PX3(EP)*DPRZ(3,EP)
789 mf12(1,3,ep)=mf12(1,3,ep)-py1(ep)*dprz(2,ep)
790 mf12(2,3,ep)=mf12(2,3,ep)+px1(ep)*dprz(2,ep)
791 mf13(1,3,ep)=mf13(1,3,ep)-py1(ep)*dprz(3,ep)
792 mf13(2,3,ep)=mf13(2,3,ep)+px1(ep)*dprz(3,ep)
793 mf23(1,3,ep)=mf23(1,3,ep)-py2(ep)*dprz(3,ep)
794 mf23(2,3,ep)=mf23(2,3,ep)-px1(ep)*dprz(3,ep)
795C
796 fm12(3,1,ep)=fm12(3,1,ep)-py2(ep)*dprz(1,ep)
797 fm12(3,2,ep)=fm12(3,2,ep)-px1(ep)*dprz(1,ep)
798 fm13(3,1,ep)=fm13(3,1,ep)-py3(ep)*dprz(1,ep)
799C FM13(3,2,EP)=FM13(3,2,EP)+PX3(EP)*DPRZ(1,EP)
800 fm23(3,1,ep)=fm23(3,1,ep)-py3(ep)*dprz(2,ep)
801C FM23(3,2,EP)=FM23(3,2,EP)+PX3(EP)*DPRZ(2,EP)
802 ENDDO
803C------ Mzz ----------
804 DO ep=jft,jlt
805 m11(3,3,ep)=m11(3,3,ep)+prz(1,ep)*dprz(1,ep)
806 m12(3,3,ep)=m12(3,3,ep)+prz(1,ep)*dprz(2,ep)
807 m13(3,3,ep)=m13(3,3,ep)+prz(1,ep)*dprz(3,ep)
808 m22(3,3,ep)=m22(3,3,ep)+prz(2,ep)*dprz(2,ep)
809 m23(3,3,ep)=m23(3,3,ep)+prz(2,ep)*dprz(3,ep)
810 m33(3,3,ep)=m33(3,3,ep)+prz(3,ep)*dprz(3,ep)
811 ENDDO
812 END DO !NG =1,NPG
813C-------[MFIJ]=[Bm]^t[C][BRm]; [MIJ]=[BRm]^t[C][BRm];----
814C---------------|0 0 BRX |----
815C-----BR(3x3)= |0 0 BRY |
816C---------------|0 0 BRXY+BRYX |----
817 DO ep=jft,jlt
818 DO j=1,2
819 prx(j,ep) =bm0rz(ep,1,j)*a_i(ep)
820 pry(j,ep) =bm0rz(ep,2,j)*a_i(ep)
821 prxy(j,ep)=bm0rz(ep,3,j)*a_i(ep)
822 ENDDO
823 prx(3,ep) =-prx(1,ep)-prx(2,ep)
824 pry(3,ep) =-pry(1,ep)-pry(2,ep)
825 prxy(3,ep)=-prxy(1,ep)-prxy(2,ep)
826 END DO
827C
828 DO ep=jft,jlt
829 cbrx(1,ep) =dm(1,1,ep)*prx(1,ep)+dm(1,2,ep)*pry(1,ep)+
830 . dm(1,3,ep)*prxy(1,ep)
831 cbry(1,ep) =dm(2,1,ep)*prx(1,ep)+dm(2,2,ep)*pry(1,ep)+
832 . dm(2,3,ep)*prxy(1,ep)
833 cbrz(1,ep) =dm(3,1,ep)*prx(1,ep)+dm(3,2,ep)*pry(1,ep)+
834 . dm(3,3,ep)*prxy(1,ep)
835 cbrx(2,ep) =dm(1,1,ep)*prx(2,ep)+dm(1,2,ep)*pry(2,ep)+
836 . dm(1,3,ep)*prxy(2,ep)
837 cbry(2,ep) =dm(2,1,ep)*prx(2,ep)+dm(2,2,ep)*pry(2,ep)+
838 . dm(2,3,ep)*prxy(2,ep)
839 cbrz(2,ep) =dm(3,1,ep)*prx(2,ep)+dm(3,2,ep)*pry(2,ep)+
840 . dm(3,3,ep)*prxy(2,ep)
841 cbrx(3,ep) =-cbrx(1,ep)-cbrx(2,ep)
842 cbry(3,ep) =-cbry(1,ep)-cbry(2,ep)
843 cbrz(3,ep) =-cbrz(1,ep)-cbrz(2,ep)
844 ENDDO
845C
846C-- [MFIJ(,3)]=|Bxi 0 Byi|{CPRXj} [FMIJ(3,)]={CPRXi CPRYi CPRZi}|Bxj 0 0|
847C-- |0 Byi Bxi|{CPRYj} |0 Byj 0|
848C-- |0 0 0 |{CPRZj} |ByjBxj 0|
849 DO ep=jft,jlt
850 mf11(1,3,ep)=mf11(1,3,ep)+ px1(ep)*cbrx(1,ep)+py1(ep)*cbrz(1,ep)
851 mf11(2,3,ep)=mf11(2,3,ep)+ py1(ep)*cbry(1,ep)+px1(ep)*cbrz(1,ep)
852 mf12(1,3,ep)=mf12(1,3,ep)+ px1(ep)*cbrx(2,ep)+py1(ep)*cbrz(2,ep)
853 mf12(2,3,ep)=mf12(2,3,ep)+ py1(ep)*cbry(2,ep)+px1(ep)*cbrz(2,ep)
854 mf13(1,3,ep)=mf13(1,3,ep)+ px1(ep)*cbrx(3,ep)+py1(ep)*cbrz(3,ep)
855 mf13(2,3,ep)=mf13(2,3,ep)+ py1(ep)*cbry(3,ep)+px1(ep)*cbrz(3,ep)
856 mf22(1,3,ep)=mf22(1,3,ep)- px1(ep)*cbrx(2,ep)+py2(ep)*cbrz(2,ep)
857 mf22(2,3,ep)=mf22(2,3,ep)+ py2(ep)*cbry(2,ep)-px1(ep)*cbrz(2,ep)
858 mf23(1,3,ep)=mf23(1,3,ep)- px1(ep)*cbrx(3,ep)+py2(ep)*cbrz(3,ep)
859 mf23(2,3,ep)=mf23(2,3,ep)+ py2(ep)*cbry(3,ep)-px1(ep)*cbrz(3,ep)
860 mf33(1,3,ep)=mf33(1,3,ep) +py3(ep)*cbrz(3,ep)
861 mf33(2,3,ep)=mf33(2,3,ep)+ py3(ep)*cbry(3,ep)
862C
863 fm12(3,1,ep)=fm12(3,1,ep)- px1(ep)*cbrx(1,ep)+py2(ep)*cbrz(1,ep)
864 fm12(3,2,ep)=fm12(3,2,ep)+ py2(ep)*cbry(1,ep)-px1(ep)*cbrz(1,ep)
865 fm13(3,1,ep)=fm13(3,1,ep) +py3(ep)*cbrz(1,ep)
866 fm13(3,2,ep)=fm13(3,2,ep)+ py3(ep)*cbry(1,ep)
867 fm23(3,1,ep)=fm23(3,1,ep) +py3(ep)*cbrz(2,ep)
868 fm23(3,2,ep)=fm23(3,2,ep)+ py3(ep)*cbry(2,ep)
869 ENDDO
870C------ Mzz ----------
871 DO ep=jft,jlt
872 m11(3,3,ep)=m11(3,3,ep)+prx(1,ep)*cbrx(1,ep)+
873 . pry(1,ep)*cbry(1,ep)+prxy(1,ep)*cbrz(1,ep)
874 m12(3,3,ep)=m12(3,3,ep)+prx(1,ep)*cbrx(2,ep)+
875 . pry(1,ep)*cbry(2,ep)+prxy(1,ep)*cbrz(2,ep)
876 m13(3,3,ep)=m13(3,3,ep)+prx(1,ep)*cbrx(3,ep)+
877 . pry(1,ep)*cbry(3,ep)+prxy(1,ep)*cbrz(3,ep)
878 m22(3,3,ep)=m22(3,3,ep)+prx(2,ep)*cbrx(2,ep)+
879 . pry(2,ep)*cbry(2,ep)+prxy(2,ep)*cbrz(2,ep)
880 m23(3,3,ep)=m23(3,3,ep)+prx(2,ep)*cbrx(3,ep)+
881 . pry(2,ep)*cbry(3,ep)+prxy(2,ep)*cbrz(3,ep)
882 m33(3,3,ep)=m33(3,3,ep)+prx(3,ep)*cbrx(3,ep)+
883 . pry(3,ep)*cbry(3,ep)+prxy(3,ep)*cbrz(3,ep)
884 ENDDO
885 IF (iorth>0) THEN
886C ----add mem-bending coupling terms of orthotrope--coplanar first---
887#include "vectorize.inc"
888 DO ep=jft,jlt
889 c2=vol(ep)*thk0(ep)
890 dmf(1,1,ep)=hmfor(ep,1)*c2
891 dmf(2,2,ep)=hmfor(ep,2)*c2
892 dmf(1,2,ep)=hmfor(ep,3)*c2
893 dmf(1,3,ep)=hmfor(ep,5)*c2
894 dmf(2,3,ep)=hmfor(ep,6)*c2
895 dmf(2,1,ep)=dmf(1,2,ep)
896 dmf(3,1,ep)=dmf(1,3,ep)
897 dmf(3,2,ep)=dmf(2,3,ep)
898 dmf(3,3,ep)=hmfor(ep,4)*c2
899 ENDDO
900C-------[C][BRm];----
901 DO ep=jft,jlt
902 cbrx(1,ep) =dmf(1,1,ep)*prx(1,ep)+dmf(1,2,ep)*pry(1,ep)+
903 . dmf(1,3,ep)*prxy(1,ep)
904 cbry(1,ep) =dmf(2,1,ep)*prx(1,ep)+dmf(2,2,ep)*pry(1,ep)+
905 . dmf(2,3,ep)*prxy(1,ep)
906 cbrz(1,ep) =dmf(3,1,ep)*prx(1,ep)+dmf(3,2,ep)*pry(1,ep)+
907 . dmf(3,3,ep)*prxy(1,ep)
908 cbrx(2,ep) =dmf(1,1,ep)*prx(2,ep)+dmf(1,2,ep)*pry(2,ep)+
909 . dmf(1,3,ep)*prxy(2,ep)
910 cbry(2,ep) =dmf(2,1,ep)*prx(2,ep)+dmf(2,2,ep)*pry(2,ep)+
911 . dmf(2,3,ep)*prxy(2,ep)
912 cbrz(2,ep) =dmf(3,1,ep)*prx(2,ep)+dmf(3,2,ep)*pry(2,ep)+
913 . dmf(3,3,ep)*prxy(2,ep)
914 cbrx(3,ep) =-cbrx(1,ep)-cbrx(2,ep)
915 cbry(3,ep) =-cbry(1,ep)-cbry(2,ep)
916 cbrz(3,ep) =-cbrz(1,ep)-cbrz(2,ep)
917 ENDDO
918C ----add Rz-bending coupling terms of orthotrope--coplanar first---
919C MIJ(1,3)=-BY(I)*CBR2(J)-BX(I)*CBR3(J)
920C MIJ(2,3)=BX(I)*CBR1(J)+BY(I)*CBR3(J)
921C MIJ(3,1)=-BY(J)*CBR2(I)-BX(J)*CBR3(I)
922C MIJ(3,2)=BX(J)*CBR1(I)+BY(J)*CBR3(I)
923C
924C------ PX2=-PX1 , PX3=0-------------
925 DO ep=jft,jlt
926 m11(1,3,ep)=-py1(ep)*cbry(1,ep)-px1(ep)*cbrz(1,ep)
927 m11(2,3,ep)= px1(ep)*cbrx(1,ep)+py1(ep)*cbrz(1,ep)
928 m12(1,3,ep)=-py1(ep)*cbry(2,ep)-px1(ep)*cbrz(2,ep)
929 m12(2,3,ep)= px1(ep)*cbrx(2,ep)+py1(ep)*cbrz(2,ep)
930 m13(1,3,ep)=-py1(ep)*cbry(3,ep)-px1(ep)*cbrz(3,ep)
931 m13(2,3,ep)= px1(ep)*cbrx(3,ep)+py1(ep)*cbrz(3,ep)
932 m22(1,3,ep)=-py2(ep)*cbry(2,ep)+px1(ep)*cbrz(2,ep)
933 m22(2,3,ep)=-px1(ep)*cbrx(2,ep)+py2(ep)*cbrz(2,ep)
934 m23(1,3,ep)=-py2(ep)*cbry(3,ep)+px1(ep)*cbrz(3,ep)
935 m23(2,3,ep)=-px1(ep)*cbrx(3,ep)+py2(ep)*cbrz(3,ep)
936 m33(1,3,ep)= py3(ep)*cbry(3,ep)
937 m33(2,3,ep)=-py3(ep)*cbrz(3,ep)
938C
939 m12(3,1,ep)=-py2(ep)*cbry(1,ep)+px1(ep)*cbrz(1,ep)
940 m12(3,2,ep)=-px1(ep)*cbrx(1,ep)+py2(ep)*cbrz(1,ep)
941 m13(3,1,ep)= py3(ep)*cbry(1,ep)
942 m13(3,2,ep)= -py3(ep)*cbrz(1,ep)
943 m23(3,1,ep)=-py3(ep)*cbry(2,ep)
944 m23(3,2,ep)= py3(ep)*cbrz(2,ep)
945 ENDDO
946 ENDIF !IF (IORTH==1)
947C
948 RETURN
949 END
950
c3ke3
subroutine c3ke3(jft, jlt, nft, npt, mtn, ithk, ncycle, istrain, ipla, pm, geo, ixtg, elbuf_str, bufmat, offset, indxof, etag, iddl, ndof, k_diag, k_lt, iadk, jdik, ihbe, thke, ismstr, x, ikgeo, ipm, igeo, iexpan, iparg, isubstack, stack, drape_sh3n, indx_drape, sedrape, numel_drape)
Definition c3ke3.F:54
c3lke3
subroutine c3lke3(jft, jlt, area, thk0, thk2, hm, hf, hc, hz, px1, py1, py2, vol, k11, k12, k13, k22, k23, k33, m11, m12, m13, m22, m23, m33, mf11, mf12, mf13, mf22, mf23, mf33, fm12, fm13, fm23, ikgeo, for, mom, iorth, hmor, hfor, hmfor, idril, nel)
Definition c3lke3.F:36
c3lkerz3
subroutine c3lkerz3(jft, jlt, hm, hz, px1, py1, py2, vol, area, k11, k12, k13, k22, k23, k33, m11, m12, m13, m22, m23, m33, mf11, mf12, mf13, mf22, mf23, mf33, fm12, fm13, fm23, iorth, hmor, bm0rz, b0rz, bkrz, berz, thk0, hmfor)
Definition c3lke3.F:572
area
subroutine area(d1, x, x2, y, y2, eint, stif0)
Definition hm_read_prop32.F:567
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
for
for(i8=*sizetab-1;i8 >=0;i8--)
Definition mumps_common.c:91