FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
static_LIB_beam.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
7
8 use hecmw
9 use mmechgauss
10 use m_utilities
11 use elementinfo
12
13 implicit none
14
15contains
16
17
18 subroutine framtr(refx, xl, le, t)
19 real(kind=kreal), intent(in) :: refx(3)
20 real(kind=kreal), intent(in) :: xl(3,2)
21 real(kind=kreal), intent(out) :: le
22 real(kind=kreal), intent(out) :: t(3,3)
23
24 real(kind=kreal) :: dl
25 real(kind=kreal), parameter :: tol = 1.d-08
26
27 t(1,1) = xl(1,2) - xl(1,1)
28 t(1,2) = xl(2,2) - xl(2,1)
29 t(1,3) = xl(3,2) - xl(3,1)
30 le = sqrt(t(1,1)*t(1,1)+t(1,2)*t(1,2)+t(1,3)*t(1,3))
31 dl = 1.0d0/le
32 t(1,1) = t(1,1)*dl
33 t(1,2) = t(1,2)*dl
34 t(1,3) = t(1,3)*dl
35
36 t(3,1) = refx(1)
37 t(3,2) = refx(2)
38 t(3,3) = refx(3)
39
40
41 t(2,1) = (t(3,2)*t(1,3) - t(3,3)*t(1,2))
42 t(2,2) = (t(3,3)*t(1,1) - t(3,1)*t(1,3))
43 t(2,3) = (t(3,1)*t(1,2) - t(3,2)*t(1,1))
44 dl = sqrt(t(2,1)*t(2,1)+t(2,2)*t(2,2)+t(2,3)*t(2,3))
45 if(dl<tol*le) then
46 stop "Bad reference for beam element!"
47 else
48 t(2,1) = t(2,1)/dl
49 t(2,2) = t(2,2)/dl
50 t(2,3) = t(2,3)/dl
51 t(3,1) = t(1,2)*t(2,3) - t(1,3)*t(2,2)
52 t(3,2) = t(1,3)*t(2,1) - t(1,1)*t(2,3)
53 t(3,3) = t(1,1)*t(2,2) - t(1,2)*t(2,1)
54 endif
55
56 end subroutine framtr
57
59 subroutine stf_beam(etype,nn,ecoord,section,E,P,STIFF)
60 integer, intent(in) :: etype
61 integer, intent(in) :: nn
62 real(kind=kreal), intent(in) :: ecoord(3,nn)
63 real(kind=kreal), intent(in) :: section(:)
64 real(kind=kreal), intent(in) :: e,p
65 real(kind=kreal), intent(out) :: stiff(nn*6,nn*6)
66
67 real(kind=kreal) :: le, trans(3,3), refv(3), transt(3,3)
68 real(kind=kreal) :: g
69 real(kind=kreal) :: l2, l3, a, iy, iz, jx, ea, twoe, foure, twelvee, sixe
70
71 refv = section(1:3)
72 call framtr(refv, ecoord, le, trans)
73 transt= transpose(trans)
74 l2 = le*le
75 l3 = l2*le
76
77 g = e/(2.d0*(1.d0 + p))
78
79 a = section(4); iy=section(5); iz=section(6); jx=section(7)
80
81 ea = e*a/le
82 twoe = 2.d0*e/le
83 foure = 4.d0*e/le
84 twelvee = 12.d0*e/l3
85 sixe = 6.d0*e/l2
86
87 stiff = 0.d0
88 stiff(1,1) = ea;
89 stiff(7,1) = -ea;
90
91 stiff(2,2) = twelvee*iz;
92 stiff(6,2) = sixe*iz;
93 stiff(8,2) = -twelvee*iz;
94 stiff(12,2) = sixe*iz;
95
96 stiff(3,3) = twelvee*iy;
97 stiff(5,3) = -sixe*iy;
98 stiff(9,3) = -twelvee*iy;
99 stiff(11,3) = -sixe*iy;
100
101 stiff(4,4) = g*jx/le;
102 stiff(10,4) = -g*jx/le;
103
104 stiff(3,5) = -sixe*iy;
105 stiff(5,5) = foure*iy;
106 stiff(9,5) = sixe*iy;
107 stiff(11,5) = twoe*iy;
108
109 stiff(2,6) = sixe*iz;
110 stiff(6,6) = foure*iz;
111 stiff(8,6) = -sixe*iz;
112 stiff(12,6) = twoe*iz;
113
114 stiff(1,7) = -ea;
115 stiff(7,7) = ea;
116
117 stiff(2,8) = -twelvee*iz;
118 stiff(6,8) = -sixe*iz;
119 stiff(8,8) = twelvee*iz;
120 stiff(12,8) = -sixe*iz;
121
122 stiff(3,9) = -twelvee*iy;
123 stiff(5,9) = sixe*iy;
124 stiff(9,9) = twelvee*iy;
125 stiff(11,9) = sixe*iy;
126
127 stiff(4,10) = -g*jx/le;
128 stiff(10,10) = g*jx/le;
129
130 stiff(3,11) = -sixe*iy;
131 stiff(5,11) = twoe*iy;
132 stiff(9,11) = sixe*iy;
133 stiff(11,11) = foure*iy;
134
135 stiff(2,12) = sixe*iz;
136 stiff(6,12) = twoe*iz;
137 stiff(8,12) = -sixe*iz;
138 stiff(12,12) = foure*iz;
139
140 stiff(1:3,:) = matmul( transt, stiff(1:3,:) )
141 stiff(4:6,:) = matmul( transt, stiff(4:6,:) )
142 stiff(7:9,:) = matmul( transt, stiff(7:9,:) )
143 stiff(10:12,:) = matmul( transt, stiff(10:12,:) )
144
145 stiff(:,1:3) = matmul( stiff(:,1:3), trans )
146 stiff(:,4:6) = matmul( stiff(:,4:6), trans )
147 stiff(:,7:9) = matmul( stiff(:,7:9), trans )
148 stiff(:,10:12) = matmul( stiff(:,10:12), trans )
149
150 end subroutine stf_beam
151
152 !####################################################################
153 subroutine updatest_beam(etype,nn,ecoord,u,du,section,gausses,QF)
154 integer, intent(in) :: etype
155 integer, intent(in) :: nn
156 real(kind=kreal), intent(in) :: ecoord(3,nn)
157 real(kind=kreal), intent(in) :: u(6,nn)
158 real(kind=kreal), intent(in) :: du(6,nn)
159 real(kind=kreal), intent(in) :: section(:)
160 type(tgaussstatus), intent(in) :: gausses(:)
161 real(kind=kreal), intent(out) :: qf(nn*6)
162
163 real(kind=kreal) :: stiff(nn*6, nn*6), totaldisp(nn*6)
164 integer(kind=kint) :: i, j
165 real(kind=kreal) :: e,p
166
167 e = gausses(1)%pMaterial%variables(m_youngs)
168 p = gausses(1)%pMaterial%variables(m_poisson)
169
170 call stf_beam(etype,nn,ecoord,section,e,p,stiff)
171
172 do i=1,nn
173 do j=1,6
174 totaldisp(6*(i-1)+j) = u(j,i) + du(j,i)
175 end do
176 end do
177
178 qf = matmul(stiff,totaldisp)
179
180 end subroutine updatest_beam
181
182 ! (Gaku Hashimoto, The University of Tokyo, 2014/02/06) <
184 !####################################################################
185 subroutine stf_beam_641 &
186 (etype, nn, ecoord, gausses, section, stiff, tt, t0)
187 !####################################################################
188
189 use mmechgauss
190
191 !--------------------------------------------------------------------
192
193 integer, intent(in) :: etype
194 integer, intent(in) :: nn
195 real(kind=kreal), intent(in) :: ecoord(3, nn)
196 type(tgaussstatus), intent(in) :: gausses(:)
197 real(kind=kreal), intent(in) :: section(:)
198 real(kind=kreal), intent(out) :: stiff(nn*3, nn*3)
199 real(kind=kreal), intent(in), optional :: tt(nn), t0(nn)
200
201 !--------------------------------------------------------------------
202
203 real(kind = kreal) :: refv(3)
204 real(kind = kreal) :: trans(3, 3), transt(3, 3)
205 real(kind = kreal) :: ec(3, 2)
206 real(kind = kreal) :: tempc
207 real(kind = kreal) :: ina1(1), outa1(2)
208 real(kind = kreal) :: ee, pp
209 real(kind = kreal) :: le
210 real(kind = kreal) :: l2, l3, g, a, iy, iz, jx
211 real(kind = kreal) :: ea, twoe, foure, twelvee, sixe
212
213 logical :: ierr
214
215 !--------------------------------------------------------------------
216
217 refv(1) = section(1)
218 refv(2) = section(2)
219 refv(3) = section(3)
220
221 ec(1, 1) = ecoord(1, 1)
222 ec(2, 1) = ecoord(2, 1)
223 ec(3, 1) = ecoord(3, 1)
224 ec(1, 2) = ecoord(1, 2)
225 ec(2, 2) = ecoord(2, 2)
226 ec(3, 2) = ecoord(3, 2)
227
228 call framtr(refv, ec, le, trans)
229
230 transt= transpose( trans )
231
232 l2 = le*le
233 l3 = l2*le
234
235 !--------------------------------------------------------------------
236
237 if( present( tt ) ) then
238
239 tempc = 0.5d0*( tt(1)+tt(2) )
240
241 end if
242
243 !--------------------------------------------------------------------
244
245 if( present( tt ) ) then
246
247 ina1(1) = tempc
248
249 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa1, ierr, ina1 )
250
251 else
252
253 ierr = .true.
254
255 end if
256
257 !--------------------------------------------------------------
258
259 if( ierr ) then
260
261 ee = gausses(1)%pMaterial%variables(m_youngs)
262 pp = gausses(1)%pMaterial%variables(m_poisson)
263
264 else
265
266 ee = outa1(1)
267 pp = outa1(2)
268
269 end if
270
271 !--------------------------------------------------------------------
272
273 g = ee/( 2.0d0*( 1.0d0+pp ) )
274
275 a = section(4)
276
277 iy = section(5)
278 iz = section(6)
279 jx = section(7)
280
281 !--------------------------------------------------------------------
282
283 ea = ee*a/le
284
285 twoe = 2.0d0*ee/le
286 foure = 4.0d0*ee/le
287 twelvee = 12.0d0*ee/l3
288 sixe = 6.0d0*ee/l2
289
290 !--------------------------------------------------------------------
291
292 stiff = 0.0d0
293
294 stiff(1, 1) = ea
295 !stiff(7, 1) = -ea
296 stiff(4, 1) = -ea
297
298 stiff(2, 2) = twelvee*iz
299 !stiff(6, 2) = sixe*iz
300 stiff(9, 2) = sixe*iz
301 !stiff(8, 2) = -twelvee*iz
302 stiff(5, 2) = -twelvee*iz
303 stiff(12, 2) = sixe*iz
304
305 stiff(3, 3) = twelvee*iy
306 !stiff(5, 3) = -sixe*iy
307 stiff(8, 3) = -sixe*iy
308 !stiff(9, 3) = -twelvee*iy
309 stiff(6, 3) = -twelvee*iy
310 stiff(11, 3) = -sixe*iy
311
312 !stiff(4, 4) = g*jx/le
313 stiff(7, 7) = g*jx/le
314 !stiff(10, 4) = -g*jx/le
315 stiff(10, 7) = -g*jx/le
316
317 !stiff(3, 5) = -sixe*iy
318 stiff(3, 8) = -sixe*iy
319 !stiff(5, 5) = foure*iy
320 stiff(8, 8) = foure*iy
321 !stiff(9, 5) = sixe*iy
322 stiff(6, 8) = sixe*iy
323 !stiff(11, 5) = twoe*iy
324 stiff(11, 8) = twoe*iy
325
326 !stiff(2, 6) = sixe*iz
327 stiff(2, 9) = sixe*iz
328 !stiff(6, 6) = foure*iz
329 stiff(9, 9) = foure*iz
330 !stiff(8, 6) = -sixe*iz
331 stiff(5, 9) = -sixe*iz
332 !stiff(12, 6) = twoe*iz
333 stiff(12, 9) = twoe*iz
334
335 !stiff(1, 7) = -ea
336 stiff(1, 4) = -ea
337 !stiff(7, 7) = ea
338 stiff(4, 4) = ea
339
340 !stiff(2, 8) = -twelvee*iz
341 stiff(2, 5) = -twelvee*iz
342 !stiff(6, 8) = -sixe*iz
343 stiff(9, 5) = -sixe*iz
344 !stiff(8, 8) = twelvee*iz
345 stiff(5, 5) = twelvee*iz
346 !stiff(12, 8) = -sixe*iz
347 stiff(12, 5) = -sixe*iz
348
349 !stiff(3, 9) = -twelvee*iy
350 stiff(3, 6) = -twelvee*iy
351 !stiff(5, 9) = sixe*iy
352 stiff(8, 6) = sixe*iy
353 !stiff(9, 9) = twelvee*iy
354 stiff(6, 6) = twelvee*iy
355 !stiff(11, 9) = sixe*iy
356 stiff(11, 6) = sixe*iy
357
358 !stiff(4, 10) = -g*jx/le
359 stiff(7, 10) = -g*jx/le
360 stiff(10, 10) = g*jx/le
361
362 stiff(3, 11) = -sixe*iy
363 !stiff(5, 11) = twoe*iy
364 stiff(8, 11) = twoe*iy
365 !stiff(9, 11) = sixe*iy
366 stiff(6, 11) = sixe*iy
367 stiff(11, 11) = foure*iy
368
369 stiff(2, 12) = sixe*iz
370 !stiff(6, 12) = twoe*iz
371 stiff(9, 12) = twoe*iz
372 !stiff(8, 12) = -sixe*iz
373 stiff(5, 12) = -sixe*iz
374 stiff(12, 12) = foure*iz
375
376 !--------------------------------------------------------------------
377
378 stiff( 1:3, :) = matmul( transt, stiff( 1:3, :) )
379 stiff( 4:6, :) = matmul( transt, stiff( 4:6, :) )
380 stiff( 7:9, :) = matmul( transt, stiff( 7:9, :) )
381 stiff(10:12, :) = matmul( transt, stiff(10:12, :) )
382
383 stiff(:, 1:3) = matmul( stiff(:, 1:3), trans )
384 stiff(:, 4:6) = matmul( stiff(:, 4:6), trans )
385 stiff(:, 7:9) = matmul( stiff(:, 7:9), trans )
386 stiff(:, 10:12) = matmul( stiff(:, 10:12), trans )
387
388 !--------------------------------------------------------------------
389
390 return
391
392 !####################################################################
393 end subroutine stf_beam_641
394 !####################################################################
395 ! > (Gaku Hashimoto, The University of Tokyo, 2014/02/06)
396
398!####################################################################
399 SUBROUTINE nqm_beam_641 &
400 (etype, nn, ecoord, gausses, section, stiff, tt, t0, tdisp, rnqm)
401!####################################################################
402
403 USE mmechgauss
404
405!--------------------------------------------------------------------
406
407 INTEGER, INTENT(IN) :: etype
408 INTEGER, INTENT(IN) :: nn
409 REAL(kind=kreal), INTENT(IN) :: ecoord(3, nn)
410 TYPE(tgaussstatus), INTENT(IN) :: gausses(:)
411 REAL(kind=kreal), INTENT(IN) :: section(:)
412 REAL(kind=kreal), INTENT(OUT) :: stiff(nn*3, nn*3)
413 REAL(kind=kreal), INTENT(IN), OPTIONAL :: tt(nn), t0(nn)
414
415 REAL(kind=kreal), INTENT(INOUT) :: tdisp(nn*3)
416 REAL(kind=kreal), INTENT(OUT) :: rnqm(nn*3)
417
418 REAL(kind=kreal) :: tdisp1(nn*3)
419
420!--------------------------------------------------------------------
421
422 REAL(kind = kreal) :: refv(3)
423 REAL(kind = kreal) :: trans(3, 3), transt(3, 3)
424 REAL(kind = kreal) :: ec(3, 2)
425 REAL(kind = kreal) :: tempc
426 REAL(kind = kreal) :: ina1(1), outa1(2)
427 REAL(kind = kreal) :: ee, pp
428 REAL(kind = kreal) :: le
429 REAL(kind = kreal) :: l2, l3, g, a, iy, iz, jx
430 REAL(kind = kreal) :: ea, twoe, foure, twelvee, sixe
431
432 LOGICAL :: ierr
433
434!--------------------------------------------------------------------
435
436 refv(1) = section(1)
437 refv(2) = section(2)
438 refv(3) = section(3)
439
440 ec(1, 1) = ecoord(1, 1)
441 ec(2, 1) = ecoord(2, 1)
442 ec(3, 1) = ecoord(3, 1)
443 ec(1, 2) = ecoord(1, 2)
444 ec(2, 2) = ecoord(2, 2)
445 ec(3, 2) = ecoord(3, 2)
446
447 CALL framtr(refv, ec, le, trans)
448
449 transt= transpose( trans )
450
451 l2 = le*le
452 l3 = l2*le
453
454!--------------------------------------------------------------------
455
456 IF( PRESENT( tt ) ) THEN
457
458 tempc = 0.5d0*( tt(1)+tt(2) )
459
460 END IF
461
462!--------------------------------------------------------------------
463
464 IF( PRESENT( tt ) ) THEN
465
466 ina1(1) = tempc
467
468 CALL fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa1, ierr, ina1 )
469
470 ELSE
471
472 ierr = .true.
473
474 END IF
475
476 !--------------------------------------------------------------
477
478 IF( ierr ) THEN
479
480 ee = gausses(1)%pMaterial%variables(m_youngs)
481 pp = gausses(1)%pMaterial%variables(m_poisson)
482
483 ELSE
484
485 ee = outa1(1)
486 pp = outa1(2)
487
488 END IF
489
490!--------------------------------------------------------------------
491
492 g = ee/( 2.0d0*( 1.0d0+pp ) )
493
494 a = section(4)
495
496 iy = section(5)
497 iz = section(6)
498 jx = section(7)
499
500! write (6,'(a,4e15.5)') 'a,iy,iz,jx',a,iy,iz,jx
501
502!--------------------------------------------------------------------
503
504 ea = ee*a/le
505
506 twoe = 2.0d0*ee/le
507 foure = 4.0d0*ee/le
508 twelvee = 12.0d0*ee/l3
509 sixe = 6.0d0*ee/l2
510
511!--------------------------------------------------------------------
512
513 stiff = 0.0d0
514
515 stiff(1, 1) = ea
516 !stiff(7, 1) = -ea
517 stiff(4, 1) = -ea
518
519 stiff(2, 2) = twelvee*iz
520 !stiff(6, 2) = sixe*iz
521 stiff(9, 2) = sixe*iz
522 !stiff(8, 2) = -twelvee*iz
523 stiff(5, 2) = -twelvee*iz
524 stiff(12, 2) = sixe*iz
525
526 stiff(3, 3) = twelvee*iy
527 !stiff(5, 3) = -sixe*iy
528 stiff(8, 3) = -sixe*iy
529 !stiff(9, 3) = -twelvee*iy
530 stiff(6, 3) = -twelvee*iy
531 stiff(11, 3) = -sixe*iy
532
533 !stiff(4, 4) = g*jx/le
534 stiff(7, 7) = g*jx/le
535 !stiff(10, 4) = -g*jx/le
536 stiff(10, 7) = -g*jx/le
537
538 !stiff(3, 5) = -sixe*iy
539 stiff(3, 8) = -sixe*iy
540 !stiff(5, 5) = foure*iy
541 stiff(8, 8) = foure*iy
542 !stiff(9, 5) = sixe*iy
543 stiff(6, 8) = sixe*iy
544 !stiff(11, 5) = twoe*iy
545 stiff(11, 8) = twoe*iy
546
547 !stiff(2, 6) = sixe*iz
548 stiff(2, 9) = sixe*iz
549 !stiff(6, 6) = foure*iz
550 stiff(9, 9) = foure*iz
551 !stiff(8, 6) = -sixe*iz
552 stiff(5, 9) = -sixe*iz
553 !stiff(12, 6) = twoe*iz
554 stiff(12, 9) = twoe*iz
555
556 !stiff(1, 7) = -ea
557 stiff(1, 4) = -ea
558 !stiff(7, 7) = ea
559 stiff(4, 4) = ea
560
561 !stiff(2, 8) = -twelvee*iz
562 stiff(2, 5) = -twelvee*iz
563 !stiff(6, 8) = -sixe*iz
564 stiff(9, 5) = -sixe*iz
565 !stiff(8, 8) = twelvee*iz
566 stiff(5, 5) = twelvee*iz
567 !stiff(12, 8) = -sixe*iz
568 stiff(12, 5) = -sixe*iz
569
570 !stiff(3, 9) = -twelvee*iy
571 stiff(3, 6) = -twelvee*iy
572 !stiff(5, 9) = sixe*iy
573 stiff(8, 6) = sixe*iy
574 !stiff(9, 9) = twelvee*iy
575 stiff(6, 6) = twelvee*iy
576 !stiff(11, 9) = sixe*iy
577 stiff(11, 6) = sixe*iy
578
579 !stiff(4, 10) = -g*jx/le
580 stiff(7, 10) = -g*jx/le
581 stiff(10, 10) = g*jx/le
582
583 stiff(3, 11) = -sixe*iy
584 !stiff(5, 11) = twoe*iy
585 stiff(8, 11) = twoe*iy
586 !stiff(9, 11) = sixe*iy
587 stiff(6, 11) = sixe*iy
588 stiff(11, 11) = foure*iy
589
590 stiff(2, 12) = sixe*iz
591 !stiff(6, 12) = twoe*iz
592 stiff(9, 12) = twoe*iz
593 !stiff(8, 12) = -sixe*iz
594 stiff(5, 12) = -sixe*iz
595 stiff(12, 12) = foure*iz
596
597!--------------------------------------------------------------------
598 tdisp1( 1: 3 ) = matmul( trans, tdisp( 1: 3 ) )
599 tdisp1( 4: 6 ) = matmul( trans, tdisp( 4: 6 ) )
600 tdisp1( 7: 9 ) = matmul( trans, tdisp( 7: 9 ) )
601 tdisp1( 10:12 ) = matmul( trans, tdisp( 10:12 ) )
602!--------------------------------------------------------------------
603 rnqm( 1:12 ) = matmul( stiff, tdisp1 )
604
605 RETURN
606
607!####################################################################
608 END SUBROUTINE nqm_beam_641
609!####################################################################
610
611 ! (Gaku Hashimoto, The University of Tokyo, 2014/02/06) <
612 !####################################################################
613 subroutine dl_beam_641(etype, nn, xx, yy, zz, rho, ltype, params, &
614 section, vect, nsize)
615 !####################################################################
616 !**
617 !** SET DLOAD
618 !**
619 ! BX LTYPE=1 :BODY FORCE IN X-DIRECTION
620 ! BY LTYPE=2 :BODY FORCE IN Y-DIRECTION
621 ! BZ LTYPE=3 :BODY FORCE IN Z-DIRECTION
622 ! GRAV LTYPE=4 :GRAVITY FORCE
623 ! CENT LTYPE=5 :CENTRIFUGAL LOAD
624 ! P1 LTYPE=10 :TRACTION IN NORMAL-DIRECTION FOR FACE-1
625 ! P2 LTYPE=20 :TRACTION IN NORMAL-DIRECTION FOR FACE-2
626 ! P3 LTYPE=30 :TRACTION IN NORMAL-DIRECTION FOR FACE-3
627 ! P4 LTYPE=40 :TRACTION IN NORMAL-DIRECTION FOR FACE-4
628 ! P5 LTYPE=50 :TRACTION IN NORMAL-DIRECTION FOR FACE-5
629 ! P6 LTYPE=60 :TRACTION IN NORMAL-DIRECTION FOR FACE-6
630 ! I/F VARIABLES
631 integer(kind = kint), intent(in) :: etype, nn
632 real(kind = kreal), intent(in) :: xx(:), yy(:), zz(:)
633 real(kind = kreal), intent(in) :: params(0:6)
634 real(kind = kreal), intent(in) :: section(:)
635 real(kind = kreal), intent(inout) :: vect(:)
636 real(kind = kreal) :: rho
637 integer(kind = kint) :: ltype, nsize
638 ! LOCAL VARIABLES
639 integer(kind = kint) :: ndof
640 parameter(ndof = 3)
641 integer(kind = kint) :: ivol, isuf, nod(nn)
642 integer(kind = kint) :: i ,surtype, nsur
643 real(kind = kreal) :: vx, vy, vz, val, a, aa
644
645 !--------------------------------------------------------------------
646
647 val = params(0)
648
649 !--------------------------------------------------------------
650
651 ivol = 0
652 isuf = 0
653
654 if( ltype .LT. 10 ) then
655
656 ivol = 1
657
658 else if( ltype .GE. 10 ) then
659
660 isuf = 1
661
662 call getsubface(etype, ltype/10, surtype, nod)
663
664 nsur = getnumberofnodes(surtype)
665
666 end if
667
668 !--------------------------------------------------------------------
669
670 nsize = nn*ndof
671
672 !--------------------------------------------------------------------
673
674 vect(1:nsize) = 0.0d0
675
676 !--------------------------------------------------------------
677
678 ! Volume force
679
680 if( ivol .EQ. 1 ) then
681
682 if( ltype .EQ. 4 ) then
683
684 aa = dsqrt( ( xx(2)-xx(1) )*( xx(2)-xx(1) ) &
685 +( yy(2)-yy(1) )*( yy(2)-yy(1) ) &
686 +( zz(2)-zz(1) )*( zz(2)-zz(1) ) )
687
688 a = section(4)
689
690 vx = params(1)
691 vy = params(2)
692 vz = params(3)
693 vx = vx/dsqrt( params(1)**2+params(2)**2+params(3)**2 )
694 vy = vy/dsqrt( params(1)**2+params(2)**2+params(3)**2 )
695 vz = vz/dsqrt( params(1)**2+params(2)**2+params(3)**2 )
696
697 do i = 1, 2
698
699 vect(3*i-2) = val*rho*a*0.5d0*aa*vx
700 vect(3*i-1) = val*rho*a*0.5d0*aa*vy
701 vect(3*i ) = val*rho*a*0.5d0*aa*vz
702
703 end do
704
705 do i = 3, 4
706
707 vect(3*i-2) = 0.0d0
708 vect(3*i-1) = 0.0d0
709 vect(3*i ) = 0.0d0
710
711 end do
712
713 end if
714
715 end if
716
717 !--------------------------------------------------------------------
718
719 return
720
721 !####################################################################
722 end subroutine dl_beam_641
723 !####################################################################
724 ! > (Gaku Hashimoto, The University of Tokyo, 2014/02/06)
725
726
727 ! (Gaku Hashimoto, The University of Tokyo, 2014/02/06) <
728 !####################################################################
729 subroutine tload_beam_641 &
730 (etype, nn, ndof, xx, yy, zz, tt, t0, &
731 gausses, section, vect)
732 !####################################################################
733
734 use hecmw
735 use m_fstr
736 use m_utilities
737 use mmechgauss
739
740 !--------------------------------------------------------------------
741
742 integer(kind = kint), intent(in) :: etype
743 integer(kind = kint), intent(in) :: nn
744 integer(kind = kint), intent(in) :: ndof
745 type(tgaussstatus), intent(in) :: gausses(:)
746 real(kind = kreal), intent(in) :: section(:)
747 real(kind = kreal), intent(in) :: xx(nn), yy(nn), zz(nn)
748 real(kind = kreal), intent(in) :: tt(nn), t0(nn)
749 real(kind = kreal), intent(out) :: vect(nn*ndof)
750
751 !--------------------------------------------------------------------
752
753 real(kind = kreal) :: tempc, temp0
754 real(kind = kreal) :: ecoord(3, nn)
755 real(kind = kreal) :: ec(3, 2)
756 real(kind = kreal) :: ina1(1), outa1(2)
757 real(kind = kreal) :: ina2(1), outa2(1)
758 real(kind = kreal) :: alp, alp0
759 real(kind = kreal) :: ee, pp
760 real(kind = kreal) :: a
761 real(kind = kreal) :: refv(3)
762 real(kind = kreal) :: g
763 real(kind = kreal) :: le
764 real(kind = kreal) :: trans(3, 3), transt(3, 3)
765
766 logical :: ierr
767
768 !--------------------------------------------------------------------
769
770 ecoord(1, 1:nn) = xx(1:nn)
771 ecoord(2, 1:nn) = yy(1:nn)
772 ecoord(3, 1:nn) = zz(1:nn)
773
774 !--------------------------------------------------------------------
775
776 tempc = 0.5d0*( tt(1)+tt(2) )
777 temp0 = 0.5d0*( t0(1)+t0(2) )
778
779 !--------------------------------------------------------------
780
781 ina1(1) = tempc
782
783 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa1, ierr, ina1 )
784
785 if( ierr ) then
786
787 ee = gausses(1)%pMaterial%variables(m_youngs)
788 pp = gausses(1)%pMaterial%variables(m_poisson)
789
790 else
791
792 ee = outa1(1)
793 pp = outa1(2)
794
795 end if
796
797 !--------------------------------------------------------------
798
799 ina2(1) = tempc
800
801 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa2(:), ierr, ina2 )
802
803 if( ierr ) stop "Fails in fetching expansion coefficient!"
804
805 alp = outa2(1)
806
807 !--------------------------------------------------------------
808
809 ina2(1) = temp0
810
811 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa2(:), ierr, ina2 )
812
813 if( ierr ) stop "Fails in fetching expansion coefficient!"
814
815 alp0 = outa2(1)
816
817 !--------------------------------------------------------------------
818
819 refv(1) = section(1)
820 refv(2) = section(2)
821 refv(3) = section(3)
822
823 ec(1, 1) = ecoord(1, 1)
824 ec(2, 1) = ecoord(2, 1)
825 ec(3, 1) = ecoord(3, 1)
826 ec(1, 2) = ecoord(1, 2)
827 ec(2, 2) = ecoord(2, 2)
828 ec(3, 2) = ecoord(3, 2)
829
830 call framtr(refv, ec, le, trans)
831
832 transt= transpose( trans )
833
834 !--------------------------------------------------------------------
835
836 a = section(4)
837
838 g = ee/( 2.0d0*( 1.0d0+pp ))
839
840 !--------------------------------------------------------------------
841
842 vect( 1) = -a*ee*( alp*( tempc-ref_temp )-alp0*( temp0-ref_temp ) )
843 vect( 2) = 0.0d0
844 vect( 3) = 0.0d0
845
846 vect( 4) = a*ee*( alp*( tempc-ref_temp )-alp0*( temp0-ref_temp ) )
847 vect( 5) = 0.0d0
848 vect( 6) = 0.0d0
849
850 vect( 7) = 0.0d0
851 vect( 8) = 0.0d0
852 vect( 9) = 0.0d0
853
854 vect(10) = 0.0d0
855 vect(11) = 0.0d0
856 vect(12) = 0.0d0
857
858 !--------------------------------------------------------------------
859
860 vect( 1:3) = matmul( transt, vect(1:3) )
861 vect( 4:6) = matmul( transt, vect(4:6) )
862 vect( 7:9) = matmul( transt, vect(7:9) )
863 vect(10:12) = matmul( transt, vect(10:12) )
864
865 !--------------------------------------------------------------------
866
867 return
868
869 !####################################################################
870 end subroutine tload_beam_641
871 !####################################################################
872 ! > (Gaku Hashimoto, The University of Tokyo, 2013/09/13)
873
874
875 ! (Gaku Hashimoto, The University of Tokyo, 2014/02/06) <
876 !####################################################################
878 (etype, nn, ecoord, gausses, section, edisp, &
879 ndstrain, ndstress, tt, t0, ntemp)
880 !####################################################################
881
882 use m_fstr
883 use mmechgauss
884
885 !--------------------------------------------------------------------
886
887 integer(kind = kint), intent(in) :: etype
888 integer(kind = kint), intent(in) :: nn
889 real(kind = kreal), intent(in) :: ecoord(3, nn)
890 type(tgaussstatus), intent(inout) :: gausses(:)
891 real(kind = kreal), intent(in) :: section(:)
892 real(kind = kreal), intent(in) :: edisp(3, nn)
893 real(kind = kreal), intent(out) :: ndstrain(nn, 6)
894 real(kind = kreal), intent(out) :: ndstress(nn, 6)
895 real(kind=kreal), intent(in), optional :: tt(nn), t0(nn)
896 integer(kind = kint), intent(in) :: ntemp
897
898 !--------------------------------------------------------------------
899
900 real(kind=kreal) :: stiffx(12, 12)
901 real(kind=kreal) :: tdisp(12)
902 real(kind=kreal) :: rnqm(12)
903
904 !--------------------------------------------------------------------
905
906 integer(kind = kint) :: i, j, k, jj
907
908 real(kind = kreal) :: tempc, temp0
909 real(kind = kreal) :: ina1(1), outa1(2)
910 real(kind = kreal) :: ina2(1), outa2(1)
911 real(kind = kreal) :: alp, alp0
912 real(kind = kreal) :: ee, pp
913 real(kind = kreal) :: a, radius, angle(6)
914 real(kind = kreal) :: refv(3)
915 real(kind = kreal) :: le, l2, l3
916 real(kind = kreal) :: trans(3, 3), transt(3, 3)
917 real(kind = kreal) :: edisp_hat(3, nn)
918 real(kind = kreal) :: ec(3, 2)
919 real(kind = kreal) :: t(3, 3), t_hat(3, 3)
920 real(kind = kreal) :: t_hat_tmp(3, 3)
921 real(kind = kreal) :: e(3, 3), e_hat(3, 3)
922 real(kind = kreal) :: e_hat_tmp(3, 3)
923 real(kind = kreal) :: x1_hat, x2_hat, x3_hat
924 real(kind = kreal) :: pi
925
926 logical :: ierr
927
928 alp = 0.0d0; alp0 = 0.0d0
929 tempc = 0.0d0; temp0 = 0.0d0
930
931 !--------------------------------------------------------------------
932
933 pi = 4.0d0*datan( 1.0d0 )
934
935 !--------------------------------------------------------------------
936
937 if( present( tt ) .AND. present( t0 ) ) then
938
939 tempc = 0.5d0*( tt(1)+tt(2) )
940 temp0 = 0.5d0*( t0(1)+t0(2) )
941
942 end if
943
944 !--------------------------------------------------------------------
945
946 if( ntemp .EQ. 1 ) then
947
948 ina1(1) = tempc
949
950 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa1, ierr, ina1 )
951
952 else
953
954 ierr = .true.
955
956 end if
957
958 !--------------------------------------------------------------
959
960 if( ierr ) then
961
962 ee = gausses(1)%pMaterial%variables(m_youngs)
963 pp = gausses(1)%pMaterial%variables(m_poisson)
964
965 else
966
967 ee = outa1(1)
968 pp = outa1(2)
969
970 end if
971
972 !--------------------------------------------------------------------
973
974 if( ntemp .EQ. 1 ) then
975
976 ina2(1) = tempc
977
978 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa2(:), ierr, ina2 )
979
980 if( ierr ) stop "Fails in fetching expansion coefficient!"
981
982 alp = outa2(1)
983
984 end if
985
986 !--------------------------------------------------------------
987
988 if( ntemp .EQ. 1 ) then
989
990 ina2(1) = temp0
991
992 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa2(:), ierr, ina2 )
993
994 if( ierr ) stop "Fails in fetching expansion coefficient!"
995
996 alp0 = outa2(1)
997
998 end if
999
1000 !--------------------------------------------------------------------
1001
1002 refv(1) = section(1)
1003 refv(2) = section(2)
1004 refv(3) = section(3)
1005
1006 ec(1, 1) = ecoord(1, 1)
1007 ec(2, 1) = ecoord(2, 1)
1008 ec(3, 1) = ecoord(3, 1)
1009 ec(1, 2) = ecoord(1, 2)
1010 ec(2, 2) = ecoord(2, 2)
1011 ec(3, 2) = ecoord(3, 2)
1012
1013 call framtr(refv, ec, le, trans)
1014
1015 transt= transpose( trans )
1016
1017 l2 = le*le
1018 l3 = l2*le
1019
1020 !--------------------------------------------------------------------
1021
1022 a = section(4)
1023
1024 radius = gausses(1)%pMaterial%variables(m_beam_radius)
1025
1026 angle(1) = gausses(1)%pMaterial%variables(m_beam_angle1)
1027 angle(2) = gausses(1)%pMaterial%variables(m_beam_angle2)
1028 angle(3) = gausses(1)%pMaterial%variables(m_beam_angle3)
1029 angle(4) = gausses(1)%pMaterial%variables(m_beam_angle4)
1030 angle(5) = gausses(1)%pMaterial%variables(m_beam_angle5)
1031 angle(6) = gausses(1)%pMaterial%variables(m_beam_angle6)
1032
1033 !--------------------------------------------------------------------
1034
1035 do k = 1, 6
1036
1037 !--------------------------------------------------------
1038
1039 angle(k) = angle(k)/180.0d0*pi
1040
1041 x2_hat = radius*dcos( angle(k) )
1042 x3_hat = radius*dsin( angle(k) )
1043
1044 !--------------------------------------------------------
1045
1046 jj = 0
1047 do j = 1, nn
1048
1049 do i = 1, 3
1050
1051 edisp_hat(i, j) = trans(i, 1)*edisp(1, j) &
1052 +trans(i, 2)*edisp(2, j) &
1053 +trans(i, 3)*edisp(3, j)
1054
1055 jj = jj + 1
1056 tdisp(jj) = edisp(i,j)
1057
1058 end do
1059
1060 end do
1061
1062 !--------------------------------------------------------
1063
1064 x1_hat = 0.5d0*le
1065
1066 e_hat = 0.0d0
1067 e_hat(1, 1) = ( edisp_hat(1, 2)-edisp_hat(1, 1) )/le
1068
1069 t_hat = 0.0d0
1070 t_hat(1, 1) = ee*( edisp_hat(1, 2)-edisp_hat(1, 1) )/le &
1071 -ee*x2_hat*( ( -6.0d0/l2+12.0d0*x1_hat/l3 )*edisp_hat(2, 1) &
1072 +( -4.0d0/le+6.0d0*x1_hat/l2 )*edisp_hat(3, 3) &
1073 +( 6.0d0/l2-12.0d0*x1_hat/l3 )*edisp_hat(2, 2) &
1074 +( -2.0d0/le+6.0d0*x1_hat/l2 )*edisp_hat(3, 4) ) &
1075 -ee*x3_hat*( ( -6.0d0/l2+12.0d0*x1_hat/l3 )*edisp_hat(3, 1) &
1076 +( 4.0d0/le-6.0d0*x1_hat/l2 )*edisp_hat(2, 3) &
1077 +( 6.0d0/l2-12.0d0*x1_hat/l3 )*edisp_hat(3, 2) &
1078 +( 2.0d0/le-6.0d0*x1_hat/l2 )*edisp_hat(2, 4) )
1079
1080 if( ntemp .EQ. 1 ) then
1081
1082 t_hat(1, 1) &
1083 = t_hat(1, 1) &
1084 -ee*( alp*( tempc-ref_temp )-alp0*( temp0-ref_temp ) )
1085
1086 end if
1087
1088 e_hat_tmp(1:3,:) = matmul( trans, e_hat(1:3,:) )
1089 t_hat_tmp(1:3,:) = matmul( trans, t_hat(1:3,:) )
1090
1091 e(:, 1:3) = matmul( e_hat_tmp(:,1:3), transt )
1092 t(:, 1:3) = matmul( t_hat_tmp(:,1:3), transt )
1093
1094 gausses(1)%strain(k) = e_hat(1, 1)
1095 gausses(1)%stress(k) = t_hat(1, 1)
1096
1097 !set stress and strain for output
1098 gausses(1)%strain_out(k) = gausses(1)%strain(k)
1099 gausses(1)%stress_out(k) = gausses(1)%stress(k)
1100
1101 !--------------------------------------------------------
1102
1103 ndstrain(1, k) = 0.0d0
1104 ndstrain(2, k) = 0.0d0
1105 ndstrain(3, k) = 0.0d0
1106 ndstrain(4, k) = 0.0d0
1107
1108 ndstress(1, k) = 0.0d0
1109 ndstress(2, k) = 0.0d0
1110 ndstress(3, k) = 0.0d0
1111 ndstress(4, k) = 0.0d0
1112
1113 !--------------------------------------------------------
1114
1115 x1_hat = 0.0d0
1116
1117 e_hat = 0.0d0
1118 e_hat(1, 1) = ( edisp_hat(1, 2)-edisp_hat(1, 1) )/le
1119
1120 t_hat = 0.0d0
1121 t_hat(1, 1) = ee*( edisp_hat(1, 2)-edisp_hat(1, 1) )/le &
1122 -ee*x2_hat*( ( -6.0d0/l2+12.0d0*x1_hat/l3 )*edisp_hat(2, 1) &
1123 +( -4.0d0/le+6.0d0*x1_hat/l2 )*edisp_hat(3, 3) &
1124 +( 6.0d0/l2-12.0d0*x1_hat/l3 )*edisp_hat(2, 2) &
1125 +( -2.0d0/le+6.0d0*x1_hat/l2 )*edisp_hat(3, 4) ) &
1126 -ee*x3_hat*( ( -6.0d0/l2+12.0d0*x1_hat/l3 )*edisp_hat(3, 1) &
1127 +( 4.0d0/le-6.0d0*x1_hat/l2 )*edisp_hat(2, 3) &
1128 +( 6.0d0/l2-12.0d0*x1_hat/l3 )*edisp_hat(3, 2) &
1129 +( 2.0d0/le-6.0d0*x1_hat/l2 )*edisp_hat(2, 4) )
1130
1131 if( ntemp .EQ. 1 ) then
1132
1133 t_hat(1, 1) &
1134 = t_hat(1, 1) &
1135 -ee*( alp*( tempc-ref_temp )-alp0*( temp0-ref_temp ) )
1136
1137 end if
1138
1139 e_hat_tmp(1:3, :) = matmul( trans, e_hat(1:3, :) )
1140 t_hat_tmp(1:3, :) = matmul( trans, t_hat(1:3, :) )
1141
1142 e(:, 1:3) = matmul( e_hat_tmp(:, 1:3), transt )
1143 t(:, 1:3) = matmul( t_hat_tmp(:, 1:3), transt )
1144
1145 ndstrain(1, k) = e_hat(1, 1)
1146 ndstress(1, k) = t_hat(1, 1)
1147
1148 !--------------------------------------------------------
1149
1150 x1_hat = le
1151
1152 e_hat = 0.0d0
1153 e_hat(1, 1) = ( edisp_hat(1, 2)-edisp_hat(1, 1) )/le
1154
1155 t_hat = 0.0d0
1156 t_hat(1, 1) = ee*( edisp_hat(1, 2)-edisp_hat(1, 1) )/le &
1157 -ee*x2_hat*( ( -6.0d0/l2+12.0d0*x1_hat/l3 )*edisp_hat(2, 1) &
1158 +( -4.0d0/le+6.0d0*x1_hat/l2 )*edisp_hat(3, 3) &
1159 +( 6.0d0/l2-12.0d0*x1_hat/l3 )*edisp_hat(2, 2) &
1160 +( -2.0d0/le+6.0d0*x1_hat/l2 )*edisp_hat(3, 4) ) &
1161 -ee*x3_hat*( ( -6.0d0/l2+12.0d0*x1_hat/l3 )*edisp_hat(3, 1) &
1162 +( 4.0d0/le-6.0d0*x1_hat/l2 )*edisp_hat(2, 3) &
1163 +( 6.0d0/l2-12.0d0*x1_hat/l3 )*edisp_hat(3, 2) &
1164 +( 2.0d0/le-6.0d0*x1_hat/l2 )*edisp_hat(2, 4) )
1165
1166 if( ntemp .EQ. 1 ) then
1167
1168 t_hat(1, 1) &
1169 = t_hat(1, 1) &
1170 -ee*( alp*( tempc-ref_temp )-alp0*( temp0-ref_temp ) )
1171
1172 end if
1173
1174 e_hat_tmp(1:3, :) = matmul( trans, e_hat(1:3, :) )
1175 t_hat_tmp(1:3, :) = matmul( trans, t_hat(1:3, :) )
1176
1177 e(:, 1:3) = matmul( e_hat_tmp(:, 1:3), transt )
1178 t(:, 1:3) = matmul( t_hat_tmp(:, 1:3), transt )
1179
1180 ndstrain(2, k) = e_hat(1, 1)
1181 ndstress(2, k) = t_hat(1, 1)
1182
1183 !--------------------------------------------------------
1184
1185 end do
1186
1187 !--------------------------------------------------------------------
1188 stiffx = 0.0
1189
1190 call nqm_beam_641 &
1191 (etype, nn, ecoord, gausses, section, stiffx, tt, t0, tdisp, rnqm )
1192
1193 gausses(1)%nqm(1:12) = rnqm(1:12)
1194
1195! write (6,'(a5,6a15)') 'dis-ij','x','y','z','theta-x','theta-y','theta-z'
1196! write (6,'(a,1p,6e15.5,0p)') 'dis-i',(tdisp(j),j= 1, 3),(tdisp(j),j= 7, 9)
1197! write (6,'(a,1p,6e15.5,0p)') 'dis-j',(tdisp(j),j= 4, 6),(tdisp(j),j=10,12)
1198! write (6,'(a5,6a15)') 'nqm-ij','N','Qy','QZ','Mx','My','Mz'
1199! write (6,'(a,1p,6e15.5,0p)') 'nqm-i',(rnqm(j),j= 1, 3),(rnqm(j),j= 7, 9)
1200! write (6,'(a,1p,6e15.5,0p)') 'nqm-j',(rnqm(j),j= 4, 6),(rnqm(j),j=10,12)
1201! write (6,'(a)') ''
1202
1203 !--------------------------------------------------------------------
1204
1205 return
1206
1207 !####################################################################
1208 end subroutine nodalstress_beam_641
1209 !####################################################################
1210 ! > (Gaku Hashimoto, The University of Tokyo, 2013/09/13)
1211
1212 !####################################################################
1214 ( gausses, estrain, estress, enqm )
1215 !####################################################################
1216 use m_fstr
1217 use mmechgauss
1218 implicit none
1219
1220 !--------------------------------------------------------------------
1221
1222 type(tgaussstatus), intent(inout) :: gausses(:)
1223 real(kind = kreal), intent(out) :: estrain(6)
1224 real(kind = kreal), intent(out) :: estress(6)
1225 real(kind = kreal), intent(out) :: enqm(12)
1226
1227 !--------------------------------------------------------------------
1228
1229 estrain(1:6) = gausses(1)%strain_out(1:6)
1230 estress(1:6) = gausses(1)%stress_out(1:6)
1231 enqm(1:12) = gausses(1)%nqm(1:12)
1232
1233 end subroutine elementalstress_beam_641
1234
1235 !####################################################################
1237 (etype, nn, ecoord, u, du, gausses, section, qf, tt, t0)
1238 !####################################################################
1239
1240 use mmechgauss
1241
1242 !--------------------------------------------------------------------
1243
1244 integer, intent(in) :: etype
1245 integer, intent(in) :: nn
1246 real(kind=kreal), intent(in) :: ecoord(3, nn)
1247 real(kind=kreal), intent(in) :: u(3, nn)
1248 real(kind=kreal), intent(in) :: du(3, nn)
1249 type(tgaussstatus), intent(in) :: gausses(:)
1250 real(kind=kreal), intent(in) :: section(:)
1251 real(kind=kreal), intent(out) :: qf(nn*3)
1252 real(kind=kreal), intent(in), optional :: tt(nn), t0(nn)
1253
1254 !--------------------------------------------------------------------
1255 real(kind = kreal) :: stiff(nn*3, nn*3), totaldisp(nn*3)
1256 integer(kind = kint) :: i
1257
1258 call stf_beam_641(etype, nn, ecoord, gausses, section, stiff, tt, t0)
1259
1260 totaldisp = 0.d0
1261 do i=1,nn
1262 totaldisp(3*i-2:3*i) = u(1:3,i) + du(1:3,i)
1263 end do
1264
1265 qf = matmul(stiff,totaldisp)
1266
1267 end subroutine updatest_beam_641
1268
1269end module
This module encapsulate the basic functions of all elements provide by this software.
Definition: element.f90:43
integer(kind=kind(2)) function getnumberofnodes(etype)
Obtain number of nodes of the element.
Definition: element.f90:126
subroutine getsubface(intype, innumber, outtype, nodes)
Find the definition of surface of the element.
Definition: element.f90:188
This module provides data for gauss quadrature.
Definition: GaussM.f90:6
Definition: hecmw.f90:6
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
real(kind=kreal), pointer ref_temp
REFTEMP.
Definition: m_fstr.f90:120
This module provide common functions of beam elements.
subroutine nqm_beam_641(etype, nn, ecoord, gausses, section, stiff, tt, t0, tdisp, rnqm)
Calculate N,Q,M vector of BEAM elements.
subroutine updatest_beam_641(etype, nn, ecoord, u, du, gausses, section, qf, tt, t0)
subroutine framtr(refx, xl, le, t)
subroutine dl_beam_641(etype, nn, xx, yy, zz, rho, ltype, params, section, vect, nsize)
subroutine elementalstress_beam_641(gausses, estrain, estress, enqm)
subroutine stf_beam(etype, nn, ecoord, section, e, p, stiff)
Calculate stiff matrix of BEAM elements.
subroutine nodalstress_beam_641(etype, nn, ecoord, gausses, section, edisp, ndstrain, ndstress, tt, t0, ntemp)
subroutine tload_beam_641(etype, nn, ndof, xx, yy, zz, tt, t0, gausses, section, vect)
subroutine stf_beam_641(etype, nn, ecoord, gausses, section, stiff, tt, t0)
Calculate stiff matrix of BEAM elements.
subroutine updatest_beam(etype, nn, ecoord, u, du, section, gausses, qf)
This module provides aux functions.
Definition: utilities.f90:6
This modules defines a structure to record history dependent parameter in static analysis.
Definition: mechgauss.f90:6
All data should be recorded in every quadrature points.
Definition: mechgauss.f90:13