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)
24 real(kind=kreal) :: dl
25 real(kind=kreal),
parameter :: tol = 1.d-08
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))
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))
46 stop
"Bad reference for beam element!"
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)
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)
67 real(kind=kreal) :: le, trans(3,3), refv(3), transt(3,3)
69 real(kind=kreal) :: l2, l3, a, iy, iz, jx, ea, twoe, foure, twelvee, sixe
72 call framtr(refv, ecoord, le, trans)
73 transt= transpose(trans)
77 g = e/(2.d0*(1.d0 + p))
79 a = section(4); iy=section(5); iz=section(6); jx=section(7)
91 stiff(2,2) = twelvee*iz;
93 stiff(8,2) = -twelvee*iz;
94 stiff(12,2) = sixe*iz;
96 stiff(3,3) = twelvee*iy;
97 stiff(5,3) = -sixe*iy;
98 stiff(9,3) = -twelvee*iy;
99 stiff(11,3) = -sixe*iy;
101 stiff(4,4) = g*jx/le;
102 stiff(10,4) = -g*jx/le;
104 stiff(3,5) = -sixe*iy;
105 stiff(5,5) = foure*iy;
106 stiff(9,5) = sixe*iy;
107 stiff(11,5) = twoe*iy;
109 stiff(2,6) = sixe*iz;
110 stiff(6,6) = foure*iz;
111 stiff(8,6) = -sixe*iz;
112 stiff(12,6) = twoe*iz;
117 stiff(2,8) = -twelvee*iz;
118 stiff(6,8) = -sixe*iz;
119 stiff(8,8) = twelvee*iz;
120 stiff(12,8) = -sixe*iz;
122 stiff(3,9) = -twelvee*iy;
123 stiff(5,9) = sixe*iy;
124 stiff(9,9) = twelvee*iy;
125 stiff(11,9) = sixe*iy;
127 stiff(4,10) = -g*jx/le;
128 stiff(10,10) = g*jx/le;
130 stiff(3,11) = -sixe*iy;
131 stiff(5,11) = twoe*iy;
132 stiff(9,11) = sixe*iy;
133 stiff(11,11) = foure*iy;
135 stiff(2,12) = sixe*iz;
136 stiff(6,12) = twoe*iz;
137 stiff(8,12) = -sixe*iz;
138 stiff(12,12) = foure*iz;
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,:) )
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 )
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(:)
161 real(kind=kreal),
intent(out) :: qf(nn*6)
163 real(kind=kreal) :: stiff(nn*6, nn*6), totaldisp(nn*6)
164 integer(kind=kint) :: i, j
165 real(kind=kreal) :: e,p
167 e = gausses(1)%pMaterial%variables(m_youngs)
168 p = gausses(1)%pMaterial%variables(m_poisson)
170 call stf_beam(etype,nn,ecoord,section,e,p,stiff)
174 totaldisp(6*(i-1)+j) = u(j,i) + du(j,i)
178 qf = matmul(stiff,totaldisp)
186 (etype, nn, ecoord, gausses, section, stiff, tt, t0)
193 integer,
intent(in) :: etype
194 integer,
intent(in) :: nn
195 real(kind=kreal),
intent(in) :: ecoord(3, nn)
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)
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
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)
228 call framtr(refv, ec, le, trans)
230 transt= transpose( trans )
237 if(
present( tt ) )
then
239 tempc = 0.5d0*( tt(1)+tt(2) )
245 if(
present( tt ) )
then
249 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa1, ierr, ina1 )
261 ee = gausses(1)%pMaterial%variables(m_youngs)
262 pp = gausses(1)%pMaterial%variables(m_poisson)
273 g = ee/( 2.0d0*( 1.0d0+pp ) )
287 twelvee = 12.0d0*ee/l3
298 stiff(2, 2) = twelvee*iz
300 stiff(9, 2) = sixe*iz
302 stiff(5, 2) = -twelvee*iz
303 stiff(12, 2) = sixe*iz
305 stiff(3, 3) = twelvee*iy
307 stiff(8, 3) = -sixe*iy
309 stiff(6, 3) = -twelvee*iy
310 stiff(11, 3) = -sixe*iy
313 stiff(7, 7) = g*jx/le
315 stiff(10, 7) = -g*jx/le
318 stiff(3, 8) = -sixe*iy
320 stiff(8, 8) = foure*iy
322 stiff(6, 8) = sixe*iy
324 stiff(11, 8) = twoe*iy
327 stiff(2, 9) = sixe*iz
329 stiff(9, 9) = foure*iz
331 stiff(5, 9) = -sixe*iz
333 stiff(12, 9) = twoe*iz
341 stiff(2, 5) = -twelvee*iz
343 stiff(9, 5) = -sixe*iz
345 stiff(5, 5) = twelvee*iz
347 stiff(12, 5) = -sixe*iz
350 stiff(3, 6) = -twelvee*iy
352 stiff(8, 6) = sixe*iy
354 stiff(6, 6) = twelvee*iy
356 stiff(11, 6) = sixe*iy
359 stiff(7, 10) = -g*jx/le
360 stiff(10, 10) = g*jx/le
362 stiff(3, 11) = -sixe*iy
364 stiff(8, 11) = twoe*iy
366 stiff(6, 11) = sixe*iy
367 stiff(11, 11) = foure*iy
369 stiff(2, 12) = sixe*iz
371 stiff(9, 12) = twoe*iz
373 stiff(5, 12) = -sixe*iz
374 stiff(12, 12) = foure*iz
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, :) )
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 )
400 (etype, nn, ecoord, gausses, section, stiff, tt, t0, tdisp, rnqm)
407 INTEGER,
INTENT(IN) :: etype
408 INTEGER,
INTENT(IN) :: nn
409 REAL(kind=kreal),
INTENT(IN) :: ecoord(3, nn)
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)
415 REAL(kind=kreal),
INTENT(INOUT) :: tdisp(nn*3)
416 REAL(kind=kreal),
INTENT(OUT) :: rnqm(nn*3)
418 REAL(kind=kreal) :: tdisp1(nn*3)
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
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)
447 CALL framtr(refv, ec, le, trans)
449 transt= transpose( trans )
456 IF(
PRESENT( tt ) )
THEN
458 tempc = 0.5d0*( tt(1)+tt(2) )
464 IF(
PRESENT( tt ) )
THEN
468 CALL fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa1, ierr, ina1 )
480 ee = gausses(1)%pMaterial%variables(m_youngs)
481 pp = gausses(1)%pMaterial%variables(m_poisson)
492 g = ee/( 2.0d0*( 1.0d0+pp ) )
508 twelvee = 12.0d0*ee/l3
519 stiff(2, 2) = twelvee*iz
521 stiff(9, 2) = sixe*iz
523 stiff(5, 2) = -twelvee*iz
524 stiff(12, 2) = sixe*iz
526 stiff(3, 3) = twelvee*iy
528 stiff(8, 3) = -sixe*iy
530 stiff(6, 3) = -twelvee*iy
531 stiff(11, 3) = -sixe*iy
534 stiff(7, 7) = g*jx/le
536 stiff(10, 7) = -g*jx/le
539 stiff(3, 8) = -sixe*iy
541 stiff(8, 8) = foure*iy
543 stiff(6, 8) = sixe*iy
545 stiff(11, 8) = twoe*iy
548 stiff(2, 9) = sixe*iz
550 stiff(9, 9) = foure*iz
552 stiff(5, 9) = -sixe*iz
554 stiff(12, 9) = twoe*iz
562 stiff(2, 5) = -twelvee*iz
564 stiff(9, 5) = -sixe*iz
566 stiff(5, 5) = twelvee*iz
568 stiff(12, 5) = -sixe*iz
571 stiff(3, 6) = -twelvee*iy
573 stiff(8, 6) = sixe*iy
575 stiff(6, 6) = twelvee*iy
577 stiff(11, 6) = sixe*iy
580 stiff(7, 10) = -g*jx/le
581 stiff(10, 10) = g*jx/le
583 stiff(3, 11) = -sixe*iy
585 stiff(8, 11) = twoe*iy
587 stiff(6, 11) = sixe*iy
588 stiff(11, 11) = foure*iy
590 stiff(2, 12) = sixe*iz
592 stiff(9, 12) = twoe*iz
594 stiff(5, 12) = -sixe*iz
595 stiff(12, 12) = foure*iz
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 ) )
603 rnqm( 1:12 ) = matmul( stiff, tdisp1 )
613 subroutine dl_beam_641(etype, nn, xx, yy, zz, rho, ltype, params, &
614 section, vect, nsize)
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
639 integer(kind = kint) :: ndof
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
654 if( ltype .LT. 10 )
then
658 else if( ltype .GE. 10 )
then
662 call getsubface(etype, ltype/10, surtype, nod)
674 vect(1:nsize) = 0.0d0
680 if( ivol .EQ. 1 )
then
682 if( ltype .EQ. 4 )
then
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) ) )
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 )
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
730 (etype, nn, ndof, xx, yy, zz, tt, t0, &
731 gausses, section, vect)
742 integer(kind = kint),
intent(in) :: etype
743 integer(kind = kint),
intent(in) :: nn
744 integer(kind = kint),
intent(in) :: ndof
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)
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)
770 ecoord(1, 1:nn) = xx(1:nn)
771 ecoord(2, 1:nn) = yy(1:nn)
772 ecoord(3, 1:nn) = zz(1:nn)
776 tempc = 0.5d0*( tt(1)+tt(2) )
777 temp0 = 0.5d0*( t0(1)+t0(2) )
783 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa1, ierr, ina1 )
787 ee = gausses(1)%pMaterial%variables(m_youngs)
788 pp = gausses(1)%pMaterial%variables(m_poisson)
801 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa2(:), ierr, ina2 )
803 if( ierr ) stop
"Fails in fetching expansion coefficient!"
811 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa2(:), ierr, ina2 )
813 if( ierr ) stop
"Fails in fetching expansion coefficient!"
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)
830 call framtr(refv, ec, le, trans)
832 transt= transpose( trans )
838 g = ee/( 2.0d0*( 1.0d0+pp ))
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) )
878 (etype, nn, ecoord, gausses, section, edisp, &
879 ndstrain, ndstress, tt, t0, ntemp)
887 integer(kind = kint),
intent(in) :: etype
888 integer(kind = kint),
intent(in) :: nn
889 real(kind = kreal),
intent(in) :: ecoord(3, nn)
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
900 real(kind=kreal) :: stiffx(12, 12)
901 real(kind=kreal) :: tdisp(12)
902 real(kind=kreal) :: rnqm(12)
906 integer(kind = kint) :: i, j, k, jj
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
928 alp = 0.0d0; alp0 = 0.0d0
929 tempc = 0.0d0; temp0 = 0.0d0
933 pi = 4.0d0*datan( 1.0d0 )
937 if(
present( tt ) .AND.
present( t0 ) )
then
939 tempc = 0.5d0*( tt(1)+tt(2) )
940 temp0 = 0.5d0*( t0(1)+t0(2) )
946 if( ntemp .EQ. 1 )
then
950 call fetch_tabledata( mc_isoelastic, gausses(1)%pMaterial%dict, outa1, ierr, ina1 )
962 ee = gausses(1)%pMaterial%variables(m_youngs)
963 pp = gausses(1)%pMaterial%variables(m_poisson)
974 if( ntemp .EQ. 1 )
then
978 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa2(:), ierr, ina2 )
980 if( ierr ) stop
"Fails in fetching expansion coefficient!"
988 if( ntemp .EQ. 1 )
then
992 call fetch_tabledata( mc_themoexp, gausses(1)%pMaterial%dict, outa2(:), ierr, ina2 )
994 if( ierr ) stop
"Fails in fetching expansion coefficient!"
1002 refv(1) = section(1)
1003 refv(2) = section(2)
1004 refv(3) = section(3)
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)
1013 call framtr(refv, ec, le, trans)
1015 transt= transpose( trans )
1024 radius = gausses(1)%pMaterial%variables(m_beam_radius)
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)
1039 angle(k) = angle(k)/180.0d0*pi
1041 x2_hat = radius*dcos( angle(k) )
1042 x3_hat = radius*dsin( angle(k) )
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)
1056 tdisp(jj) = edisp(i,j)
1067 e_hat(1, 1) = ( edisp_hat(1, 2)-edisp_hat(1, 1) )/le
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) )
1080 if( ntemp .EQ. 1 )
then
1088 e_hat_tmp(1:3,:) = matmul( trans, e_hat(1:3,:) )
1089 t_hat_tmp(1:3,:) = matmul( trans, t_hat(1:3,:) )
1091 e(:, 1:3) = matmul( e_hat_tmp(:,1:3), transt )
1092 t(:, 1:3) = matmul( t_hat_tmp(:,1:3), transt )
1094 gausses(1)%strain(k) = e_hat(1, 1)
1095 gausses(1)%stress(k) = t_hat(1, 1)
1098 gausses(1)%strain_out(k) = gausses(1)%strain(k)
1099 gausses(1)%stress_out(k) = gausses(1)%stress(k)
1103 ndstrain(1, k) = 0.0d0
1104 ndstrain(2, k) = 0.0d0
1105 ndstrain(3, k) = 0.0d0
1106 ndstrain(4, k) = 0.0d0
1108 ndstress(1, k) = 0.0d0
1109 ndstress(2, k) = 0.0d0
1110 ndstress(3, k) = 0.0d0
1111 ndstress(4, k) = 0.0d0
1118 e_hat(1, 1) = ( edisp_hat(1, 2)-edisp_hat(1, 1) )/le
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) )
1131 if( ntemp .EQ. 1 )
then
1139 e_hat_tmp(1:3, :) = matmul( trans, e_hat(1:3, :) )
1140 t_hat_tmp(1:3, :) = matmul( trans, t_hat(1:3, :) )
1142 e(:, 1:3) = matmul( e_hat_tmp(:, 1:3), transt )
1143 t(:, 1:3) = matmul( t_hat_tmp(:, 1:3), transt )
1145 ndstrain(1, k) = e_hat(1, 1)
1146 ndstress(1, k) = t_hat(1, 1)
1153 e_hat(1, 1) = ( edisp_hat(1, 2)-edisp_hat(1, 1) )/le
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) )
1166 if( ntemp .EQ. 1 )
then
1174 e_hat_tmp(1:3, :) = matmul( trans, e_hat(1:3, :) )
1175 t_hat_tmp(1:3, :) = matmul( trans, t_hat(1:3, :) )
1177 e(:, 1:3) = matmul( e_hat_tmp(:, 1:3), transt )
1178 t(:, 1:3) = matmul( t_hat_tmp(:, 1:3), transt )
1180 ndstrain(2, k) = e_hat(1, 1)
1181 ndstress(2, k) = t_hat(1, 1)
1191 (etype, nn, ecoord, gausses, section, stiffx, tt, t0, tdisp, rnqm )
1193 gausses(1)%nqm(1:12) = rnqm(1:12)
1214 ( gausses, estrain, estress, enqm )
1223 real(kind = kreal),
intent(out) :: estrain(6)
1224 real(kind = kreal),
intent(out) :: estress(6)
1225 real(kind = kreal),
intent(out) :: enqm(12)
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)
1237 (etype, nn, ecoord, u, du, gausses, section, qf, tt, t0)
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)
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)
1255 real(kind = kreal) :: stiff(nn*3, nn*3), totaldisp(nn*3)
1256 integer(kind = kint) :: i
1258 call stf_beam_641(etype, nn, ecoord, gausses, section, stiff, tt, t0)
1262 totaldisp(3*i-2:3*i) = u(1:3,i) + du(1:3,i)
1265 qf = matmul(stiff,totaldisp)
This module encapsulate the basic functions of all elements provide by this software.
integer(kind=kind(2)) function getnumberofnodes(etype)
Obtain number of nodes of the element.
subroutine getsubface(intype, innumber, outtype, nodes)
Find the definition of surface of the element.
This module provides data for gauss quadrature.
This module defined coomon data and basic structures for analysis.
real(kind=kreal), pointer ref_temp
REFTEMP.
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.
This modules defines a structure to record history dependent parameter in static analysis.
All data should be recorded in every quadrature points.