11 real(kind=kreal) :: coef(3)
12 integer(kind=kint) :: i, in, imat, itab, ntab
13 real(kind=kreal) :: tpoi, temp(ntab), funca(ntab+1), funcb(ntab+1)
16 if(tpoi < temp(1))
then
18 elseif(temp(ntab) <= tpoi)
then
22 if(temp(in) <= tpoi .and. tpoi < temp(in+1))
then
29 coef(1) = funca(itab)*tpoi + funcb(itab)
35 ntab, temp, funcA, funcB)
38 integer(kind=kint),
intent(in) :: etype
39 integer(kind=kint),
intent(in) :: nn
40 real(kind=kreal),
intent(in) :: temperature(nn)
41 real(kind=kreal),
intent(out) :: stiff(:,:)
42 integer(kind=kint),
parameter :: ndof = 1
43 integer(kind=kint) :: i, j, IMAT
44 integer(kind=kint) :: ntab
45 real(kind=kreal) :: ecoord(3, nn)
46 real(kind=kreal) :: dx, dy, dz, surf, val, temp_i, length
47 real(kind=kreal) :: cc(3)
48 real(kind=kreal) :: temp(ntab), funca(ntab+1), funcb(ntab+1)
51 ecoord(1,2) = ecoord(1,3); ecoord(2,2) = ecoord(2,3); ecoord(3,2) = ecoord(3,3)
54 dx = ecoord(1,2) - ecoord(1,1)
55 dy = ecoord(2,2) - ecoord(2,1)
56 dz = ecoord(3,2) - ecoord(3,1)
57 length = dsqrt(dx*dx + dy*dy + dz*dz)
60 temp_i = (temperature(1) + temperature(2))*0.5d0
63 stiff(1,1) = val*cc(1)
64 stiff(2,1) = -val*cc(1)
65 stiff(1,2) = -val*cc(1)
66 stiff(2,2) = val*cc(1)
70 ntab, temp, funcA, funcB)
76 integer(kind=kint),
intent(in) :: etype
77 integer(kind=kint),
intent(in) :: nn
78 real(kind=kreal),
intent(in) :: ecoord(2,nn)
79 real(kind=kreal),
intent(out) :: stiff(:,:)
80 real(kind=kreal),
intent(in) :: temperature(nn)
81 type(tmaterial),
pointer :: matl
82 integer(kind=kint) :: i, j, lx, imat, ntab
83 real(kind=kreal) :: naturalcoord(2)
84 real(kind=kreal) :: func(nn), thick, temp_i
85 real(kind=kreal) :: det, wg, rho, diag_stiff, total_stiff
86 real(kind=kreal) :: d(2,2), n(2,nn), dn(2,nn)
87 real(kind=kreal) :: gderiv(nn,2)
88 real(kind=kreal) :: cc(3)
89 real(kind=kreal) :: temp(ntab), funca(ntab+1), funcb(ntab+1)
99 temp_i = dot_product(func, temperature)
107 dn = matmul(d, transpose(gderiv))
108 forall(i = 1:nn, j = 1:nn)
109 stiff(i,j) = stiff(i,j) + dot_product(gderiv(i,:), dn(:,j))*wg
115 ntab, temp, funcA, funcB)
121 integer(kind=kint),
intent(in) :: etype
122 integer(kind=kint),
intent(in) :: nn
123 real(kind=kreal),
intent(in) :: ecoord(3,nn)
124 real(kind=kreal),
intent(out) :: stiff(:,:)
125 real(kind=kreal),
intent(in) :: temperature(nn)
126 type(tmaterial),
pointer :: matl
127 integer(kind=kint) :: i, j, lx, imat, ntab
128 real(kind=kreal) :: naturalcoord(3)
129 real(kind=kreal) :: func(nn), temp_i
130 real(kind=kreal) :: det, wg, rho, diag_stiff, total_stiff
131 real(kind=kreal) :: d(3, 3), n(3, nn), dn(3, nn)
132 real(kind=kreal) :: gderiv(nn, 3)
133 real(kind=kreal) :: spe(3)
134 real(kind=kreal) :: temp(ntab), funca(ntab+1), funcb(ntab+1)
144 temp_i = dot_product(func, temperature)
153 dn = matmul(d, transpose(gderiv))
154 forall(i = 1:nn, j = 1:nn)
155 stiff(i,j) = stiff(i,j) + dot_product(gderiv(i,:), dn(:,j))*wg
162 ntab, temp, funcA, funcB)
169 integer(kind=kint),
intent(in) :: etype
170 integer(kind=kint),
intent(in) :: nn
171 real(kind=kreal),
intent(in) :: ecoord(3,nn)
172 real(kind=kreal),
intent(out) :: stiff(:,:)
173 real(kind=kreal),
intent(inout) :: ss(:)
174 real(kind=kreal),
intent(inout) :: tt(nn)
175 type(tmaterial),
pointer :: matl
176 integer(kind=kint) :: i, j, lx, imat, ntab
177 integer(kind=kint) :: ig1, ig2, ig3, inod
178 real(kind=kreal) :: surf, thick, temp_i
179 real(kind=kreal) :: ri, rm, rp, si, sm, sp, ti, valx, valy, valz, var
180 real(kind=kreal) :: xj11, xj12, xj13, xj21, xj22, xj23, xj31, xj32, xj33, xsum
181 real(kind=kreal) :: det, wg, rho, diag_stiff, total_stiff
182 real(kind=kreal) :: cc(3), ctemp, dum
183 real(kind=kreal) :: temp(ntab), funca(ntab+1), funcb(ntab+1)
184 real(kind=kreal) :: xg(2), wgt(2), h(4), hr(4), hs(4)
185 real(kind=kreal) :: cod(3, 4)
186 real(kind=kreal) :: g1(3), g2(3), g3(3), e1(3), e2(3), e3(3), ref(3)
187 real(kind=kreal) :: en(3, 4), the(3, 3), amat(3, 3), bv(3), wk(3)
188 real(kind=kreal) :: dtdx(4), dtdy(4)
189 data xg/-0.5773502691896d0,0.5773502691896d0/
190 data wgt/1.0d0,1.0d0/
193 cod(1,i) = ecoord(1,i)
194 cod(2,i) = ecoord(2,i)
195 cod(3,i) = ecoord(3,i)
199 cod(1,nn) = cod(1,nn-1)
200 cod(2,nn) = cod(2,nn-1)
201 cod(3,nn) = cod(3,nn-1)
204 ref(1) = 0.25*( cod(1,2) + cod(1,3) - cod(1,1) - cod(1,4) )
205 ref(2) = 0.25*( cod(2,2) + cod(2,3) - cod(2,1) - cod(2,4) )
206 ref(3) = 0.25*( cod(3,2) + cod(3,3) - cod(3,1) - cod(3,4) )
208 g1(1) = cod(1,1) - cod(1,2)
209 g1(2) = cod(2,1) - cod(2,2)
210 g1(3) = cod(3,1) - cod(3,2)
211 g2(1) = cod(1,2) - cod(1,3)
212 g2(2) = cod(2,2) - cod(2,3)
213 g2(3) = cod(3,2) - cod(3,3)
216 g3(1) = g1(2)*g2(3) - g1(3)*g2(2)
217 g3(2) = g1(3)*g2(1) - g1(1)*g2(3)
218 g3(3) = g1(1)*g2(2) - g1(2)*g2(1)
221 xsum = dsqrt( g3(1)**2 + g3(2)**2 + g3(3)**2 )
224 en(1,i) = g3(1) / xsum
225 en(2,i) = g3(2) / xsum
226 en(3,i) = g3(3) / xsum
264 var = cod(i,inod) + thick*0.5*ti*en(i,inod)
265 g1(i) = g1(i) + hr(inod)*var
266 g2(i) = g2(i) + hs(inod)*var
267 g3(i) = g3(i) + thick*0.5*h(inod)*en(i,inod)
283 det = xj11*xj22*xj33 &
292 amat(1,1) = dum*( xj22*xj33 - xj23*xj32 )
293 amat(2,1) = dum*( -xj21*xj33 + xj23*xj31 )
294 amat(3,1) = dum*( xj21*xj32 - xj22*xj31 )
295 amat(1,2) = dum*( -xj12*xj33 + xj13*xj32 )
296 amat(2,2) = dum*( xj11*xj33 - xj13*xj31 )
297 amat(3,2) = dum*( -xj11*xj32 + xj12*xj31 )
298 amat(1,3) = dum*( xj12*xj23 - xj13*xj22 )
299 amat(2,3) = dum*( -xj11*xj23 + xj13*xj21 )
300 amat(3,3) = dum*( xj11*xj22 - xj12*xj21 )
303 xsum = dsqrt( g3(1)**2 + g3(2)**2 + g3(3)**2 )
307 e2(1) = -ref(2)*e3(3) + ref(3)*e3(2)
308 e2(2) = -ref(3)*e3(1) + ref(1)*e3(3)
309 e2(3) = -ref(1)*e3(2) + ref(2)*e3(1)
310 e1(1) = -e3(2)*e2(3) + e3(3)*e2(2)
311 e1(2) = -e3(3)*e2(1) + e3(1)*e2(3)
312 e1(3) = -e3(1)*e2(2) + e3(2)*e2(1)
313 xsum = dsqrt(e2(1)**2 + e2(2)**2 + e2(3)**2)
315 if ( xsum .GT. 1.e-15 )
then
319 e1(1) = -e3(2)*e2(3) + e3(3)*e2(2)
320 e1(2) = -e3(3)*e2(1) + e3(1)*e2(3)
321 e1(3) = -e3(1)*e2(2) + e3(2)*e2(1)
322 xsum = dsqrt( e1(1)**2 + e1(2)**2 + e1(3)**2 )
349 wk(1) = amat(1,1)*bv(1) &
352 wk(2) = amat(2,1)*bv(1) &
355 wk(3) = amat(3,1)*bv(1) &
358 bv(1) = the(1,1)*wk(1) &
361 bv(2) = the(2,1)*wk(1) &
364 bv(3) = the(3,1)*wk(1) &
374 ctemp = ctemp + h(i)*tt(i)
379 valx = cc(1)*wgt(ig1)*wgt(ig2)*wgt(ig3)*det
380 valy = cc(2)*wgt(ig1)*wgt(ig2)*wgt(ig3)*det
381 ss( 1) = ss( 1) + dtdx(1)*dtdx(1)*valx &
382 + dtdy(1)*dtdy(1)*valy
383 ss( 2) = ss( 2) + dtdx(1)*dtdx(2)*valx &
384 + dtdy(1)*dtdy(2)*valy
385 ss( 3) = ss( 3) + dtdx(1)*dtdx(3)*valx &
386 + dtdy(1)*dtdy(3)*valy
387 ss( 4) = ss( 4) + dtdx(1)*dtdx(4)*valx &
388 + dtdy(1)*dtdy(4)*valy
389 ss( 6) = ss( 6) + dtdx(2)*dtdx(2)*valx &
390 + dtdy(2)*dtdy(2)*valy
391 ss( 7) = ss( 7) + dtdx(2)*dtdx(3)*valx &
392 + dtdy(2)*dtdy(3)*valy
393 ss( 8) = ss( 8) + dtdx(2)*dtdx(4)*valx &
394 + dtdy(2)*dtdy(4)*valy
395 ss(11) = ss(11) + dtdx(3)*dtdx(3)*valx &
396 + dtdy(3)*dtdy(3)*valy
397 ss(12) = ss(12) + dtdx(3)*dtdx(4)*valx &
398 + dtdy(3)*dtdy(4)*valy
399 ss(16) = ss(16) + dtdx(4)*dtdx(4)*valx &
400 + dtdy(4)*dtdy(4)*valy
405 ss( 3) = ss( 3) + ss( 4)
406 ss( 7) = ss( 7) + ss( 8)
407 ss(11) = ss(11) + ss(16) + 2.0*ss(12)
422 ntab, temp, funcA, funcB)
428 integer(kind=kint),
intent(in) :: etype
429 integer(kind=kint),
intent(in) :: nn
430 real(kind=kreal),
intent(in) :: ecoord(3,nn)
431 real(kind=kreal),
intent(out) :: stiff(:,:)
432 real(kind=kreal),
intent(inout) :: ss(:)
433 real(kind=kreal),
intent(in) :: tt(nn)
434 type(tmaterial),
pointer :: matl
435 integer(kind=kint) :: i, j, lx, ly, imat, ntab
436 integer(kind=kint) :: ig1, ig2, ig3, inod
437 real(kind=kreal) :: ti, valx, valy, valz, var
438 real(kind=kreal) :: xj11, xj12, xj13, xj21, xj22, xj23, xj31, xj32, xj33
439 real(kind=kreal) :: naturalcoord(2), xsum
440 real(kind=kreal) :: func(nn), thick, temp_i
441 real(kind=kreal) :: det, wg, rho, diag_stiff, total_stiff
442 real(kind=kreal) :: d(1,1), n(1, nn), dn(1, nn)
443 real(kind=kreal) :: gderiv(nn,2)
444 real(kind=kreal) :: cc(3)
445 real(kind=kreal) :: temp(ntab), funca(ntab+1), funcb(ntab+1)
446 real(kind=kreal) :: xg(2), ri, si, rp, sp, rm, sm, hr(4), hs(4)
447 real(kind=kreal) :: xr, xs, yr, ys, zr, zs
448 real(kind=kreal) :: h(4), x(4), y(4), z(4)
449 real(kind=kreal) :: wgt(2), ctemp, dum
450 real(kind=kreal) :: cod(3, 4)
451 real(kind=kreal) :: g1(3), g2(3), g3(3), e1(3), e2(3), e3(3), ref(3)
452 real(kind=kreal) :: en(3, 4), the(3, 3), amat(3, 3), bv(3), wk(3)
453 real(kind=kreal) :: dtdx(4), dtdy(4)
454 data xg/-0.5773502691896d0,0.5773502691896d0/
455 data wgt/1.0d0,1.0d0/
458 cod(1,i) = ecoord(1,i)
459 cod(2,i) = ecoord(2,i)
460 cod(3,i) = ecoord(3,i)
464 ref(1) = 0.25*( cod(1,2) + cod(1,3) - cod(1,1) - cod(1,4) )
465 ref(2) = 0.25*( cod(2,2) + cod(2,3) - cod(2,1) - cod(2,4) )
466 ref(3) = 0.25*( cod(3,2) + cod(3,3) - cod(3,1) - cod(3,4) )
470 g1(1) = -0.5*cod(1,1) + 0.5*cod(1,2)
471 g1(2) = -0.5*cod(2,1) + 0.5*cod(2,2)
472 g1(3) = -0.5*cod(3,1) + 0.5*cod(3,2)
473 g2(1) = -0.5*cod(1,1) + 0.5*cod(1,4)
474 g2(2) = -0.5*cod(2,1) + 0.5*cod(2,4)
475 g2(3) = -0.5*cod(3,1) + 0.5*cod(3,4)
476 else if( i .EQ. 2 )
then
477 g1(1) = -0.5*cod(1,1) + 0.5*cod(1,2)
478 g1(2) = -0.5*cod(2,1) + 0.5*cod(2,2)
479 g1(3) = -0.5*cod(3,1) + 0.5*cod(3,2)
480 g2(1) = -0.5*cod(1,2) + 0.5*cod(1,3)
481 g2(2) = -0.5*cod(2,2) + 0.5*cod(2,3)
482 g2(3) = -0.5*cod(3,2) + 0.5*cod(3,3)
483 else if( i .EQ. 3 )
then
484 g1(1) = 0.5*cod(1,3) - 0.5*cod(1,4)
485 g1(2) = 0.5*cod(2,3) - 0.5*cod(2,4)
486 g1(3) = 0.5*cod(3,3) - 0.5*cod(3,4)
487 g2(1) = -0.5*cod(1,2) + 0.5*cod(1,3)
488 g2(2) = -0.5*cod(2,2) + 0.5*cod(2,3)
489 g2(3) = -0.5*cod(3,2) + 0.5*cod(3,3)
490 else if( i .EQ. 4 )
then
491 g1(1) = 0.5*cod(1,3) - 0.5*cod(1,4)
492 g1(2) = 0.5*cod(2,3) - 0.5*cod(2,4)
493 g1(3) = 0.5*cod(3,3) - 0.5*cod(3,4)
494 g2(1) = -0.5*cod(1,1) + 0.5*cod(1,4)
495 g2(2) = -0.5*cod(2,1) + 0.5*cod(2,4)
496 g2(3) = -0.5*cod(3,1) + 0.5*cod(3,4)
500 g3(1) = g1(2)*g2(3) - g1(3)*g2(2)
501 g3(2) = g1(3)*g2(1) - g1(1)*g2(3)
502 g3(3) = g1(1)*g2(2) - g1(2)*g2(1)
505 xsum = dsqrt( g3(1)**2 + g3(2)**2 + g3(3)**2 )
547 var = cod(i,inod) + thick*0.5*ti*en(i,inod)
548 g1(i) = g1(i) + hr(inod)*var
549 g2(i) = g2(i) + hs(inod)*var
550 g3(i) = g3(i) + thick*0.5*h(inod)*en(i,inod)
566 det = xj11*xj22*xj33 &
575 amat(1,1) = dum*( xj22*xj33 - xj23*xj32 )
576 amat(2,1) = dum*( -xj21*xj33 + xj23*xj31 )
577 amat(3,1) = dum*( xj21*xj32 - xj22*xj31 )
578 amat(1,2) = dum*( -xj12*xj33 + xj13*xj32 )
579 amat(2,2) = dum*( xj11*xj33 - xj13*xj31 )
580 amat(3,2) = dum*( -xj11*xj32 + xj12*xj31 )
581 amat(1,3) = dum*( xj12*xj23 - xj13*xj22 )
582 amat(2,3) = dum*( -xj11*xj23 + xj13*xj21 )
583 amat(3,3) = dum*( xj11*xj22 - xj12*xj21 )
586 xsum = dsqrt( g3(1)**2 + g3(2)**2 + g3(3)**2 )
590 e2(1) = -ref(2)*e3(3) + ref(3)*e3(2)
591 e2(2) = -ref(3)*e3(1) + ref(1)*e3(3)
592 e2(3) = -ref(1)*e3(2) + ref(2)*e3(1)
593 e1(1) = -e3(2)*e2(3) + e3(3)*e2(2)
594 e1(2) = -e3(3)*e2(1) + e3(1)*e2(3)
595 e1(3) = -e3(1)*e2(2) + e3(2)*e2(1)
596 xsum = dsqrt( e2(1)**2 + e2(2)**2 + e2(3)**2 )
598 if ( xsum .GT. 1.e-15 )
then
602 e1(1) = -e3(2)*e2(3) + e3(3)*e2(2)
603 e1(2) = -e3(3)*e2(1) + e3(1)*e2(3)
604 e1(3) = -e3(1)*e2(2) + e3(2)*e2(1)
605 xsum = dsqrt( e1(1)**2 + e1(2)**2 + e1(3)**2 )
631 wk(1) = amat(1,1)*bv(1) &
634 wk(2) = amat(2,1)*bv(1) &
637 wk(3) = amat(3,1)*bv(1) &
640 bv(1) = the(1,1)*wk(1) &
643 bv(2) = the(2,1)*wk(1) &
646 bv(3) = the(3,1)*wk(1) &
656 ctemp = ctemp + h(i)*tt(i)
661 valx = cc(1)*wgt(ig1)*wgt(ig2)*wgt(ig3)*det
662 valy = cc(2)*wgt(ig1)*wgt(ig2)*wgt(ig3)*det
663 ss( 1) = ss( 1) + dtdx(1)*dtdx(1)*valx &
664 + dtdy(1)*dtdy(1)*valy
665 ss( 2) = ss( 2) + dtdx(1)*dtdx(2)*valx &
666 + dtdy(1)*dtdy(2)*valy
667 ss( 3) = ss( 3) + dtdx(1)*dtdx(3)*valx &
668 + dtdy(1)*dtdy(3)*valy
669 ss( 4) = ss( 4) + dtdx(1)*dtdx(4)*valx &
670 + dtdy(1)*dtdy(4)*valy
671 ss( 6) = ss( 6) + dtdx(2)*dtdx(2)*valx &
672 + dtdy(2)*dtdy(2)*valy
673 ss( 7) = ss( 7) + dtdx(2)*dtdx(3)*valx &
674 + dtdy(2)*dtdy(3)*valy
675 ss( 8) = ss( 8) + dtdx(2)*dtdx(4)*valx &
676 + dtdy(2)*dtdy(4)*valy
677 ss(11) = ss(11) + dtdx(3)*dtdx(3)*valx &
678 + dtdy(3)*dtdy(3)*valy
679 ss(12) = ss(12) + dtdx(3)*dtdx(4)*valx &
680 + dtdy(3)*dtdy(4)*valy
681 ss(16) = ss(16) + dtdx(4)*dtdx(4)*valx &
682 + dtdy(4)*dtdy(4)*valy
696 implicit real(kind=kreal) (a - h, o - z)
697 real(kind=kreal),
intent(in) :: ecoord(3,nn)
698 real(kind=kreal),
intent(out) :: stiff(:,:)
699 dimension xxx(nn), yyy(nn), zzz(nn), temp(nn), ss(nn*nn)
700 dimension xx(4), yy(4), zz(4)
742 ha1= (( rrr1 * t1z )**2 + ( rrr2 * t5z )**2 ) &
743 & * ( rrr1 * t1z + rrr2 * t5z ) * rrr1
744 ha2= (( rrr1 * t2z )**2 + ( rrr2 * t6z )**2 ) &
745 & * ( rrr1 * t2z + rrr2 * t6z ) * rrr1
746 ha3= (( rrr1 * t3z )**2 + ( rrr2 * t7z )**2 ) &
747 & * ( rrr1 * t3z + rrr2 * t7z ) * rrr1
748 ha4= (( rrr1 * t4z )**2 + ( rrr2 * t8z )**2 ) &
749 & * ( rrr1 * t4z + rrr2 * t8z ) * rrr1
751 hb1= (( rrr1 * t1z )**2 + ( rrr2 * t5z )**2 ) &
752 & * ( rrr1 * t1z + rrr2 * t5z ) * rrr2
753 hb2= (( rrr1 * t2z )**2 + ( rrr2 * t6z )**2 ) &
754 & * ( rrr1 * t2z + rrr2 * t6z ) * rrr2
755 hb3= (( rrr1 * t3z )**2 + ( rrr2 * t7z )**2 ) &
756 & * ( rrr1 * t3z + rrr2 * t7z ) * rrr2
757 hb4= (( rrr1 * t4z )**2 + ( rrr2 * t8z )**2 ) &
758 & * ( rrr1 * t4z + rrr2 * t8z ) * rrr2
762 ss( 1) = ( hhh + ha1 ) * sa * 0.25
763 ss(10) = ( hhh + ha2 ) * sa * 0.25
764 ss(19) = ( hhh + ha3 ) * sa * 0.25
765 ss(28) = ( hhh + ha4 ) * sa * 0.25
767 ss(37) = ( hhh + hb1 ) * sb * 0.25
768 ss(46) = ( hhh + hb2 ) * sb * 0.25
769 ss(55) = ( hhh + hb3 ) * sb * 0.25
770 ss(64) = ( hhh + hb4 ) * sb * 0.25
772 sm = ( sa + sb ) * 0.5
773 hh1 = ( ha1 + hb1 ) * 0.5
774 hh2 = ( ha2 + hb2 ) * 0.5
775 hh3 = ( ha3 + hb3 ) * 0.5
776 hh4 = ( ha4 + hb4 ) * 0.5
778 ss(33) = -( hhh + hh1 ) * sm * 0.25
779 ss(42) = -( hhh + hh2 ) * sm * 0.25
780 ss(51) = -( hhh + hh3 ) * sm * 0.25
781 ss(60) = -( hhh + hh4 ) * sm * 0.25
800 implicit real(kind=kreal)(a-h,o-z)
801 dimension xx(4),yy(4),zz(4)
802 dimension xg(2),h(4),hr(4),hs(4)
806 xg(1) =-0.5773502691896258d0
835 xr=hr(1)*xx(1)+hr(2)*xx(2)+hr(3)*xx(3)+hr(4)*xx(4)
836 xs=hs(1)*xx(1)+hs(2)*xx(2)+hs(3)*xx(3)+hs(4)*xx(4)
837 yr=hr(1)*yy(1)+hr(2)*yy(2)+hr(3)*yy(3)+hr(4)*yy(4)
838 ys=hs(1)*yy(1)+hs(2)*yy(2)+hs(3)*yy(3)+hs(4)*yy(4)
839 zr=hr(1)*zz(1)+hr(2)*zz(2)+hr(3)*zz(3)+hr(4)*zz(4)
840 zs=hs(1)*zz(1)+hs(2)*zz(2)+hs(3)*zz(3)+hs(4)*zz(4)
842 det=(yr*zs-zr*ys)**2+(zr*xs-xr*zs)**2+(xr*ys-yr*xs)**2
This module encapsulate the basic functions of all elements provide by this software.
subroutine getshapefunc(fetype, localcoord, func)
Calculate the shape function in natural coodinate system.
subroutine getquadpoint(fetype, np, pos)
Fetch the coordinate of gauss point.
subroutine getglobalderiv(fetype, nn, localcoord, elecoord, det, gderiv)
Calculate shape derivative in global coordinate system.
real(kind=kreal) function getweight(fetype, np)
Fetch the weight value in given gauss point.
integer function numofquadpoints(fetype)
Obtains the number of quadrature points of the element.
This module contains subroutines used in 3d eigen analysis for.
subroutine heat_conductivity_shell_741(etype, nn, ecoord, tt, imat, thick, ss, stiff, ntab, temp, funca, funcb)
subroutine heat_conductivity_c3(etype, nn, ecoord, temperature, imat, stiff, ntab, temp, funca, funcb)
subroutine heat_get_coefficient(tpoi, imat, coef, ntab, temp, funca, funcb)
subroutine heat_conductivity_c1(etype, nn, ecoord, temperature, imat, surf, stiff, ntab, temp, funca, funcb)
subroutine heat_conductivity_c2(etype, nn, ecoord, temperature, imat, thick, stiff, ntab, temp, funca, funcb)
subroutine heat_conductivity_541(nn, ecoord, temp, tzero, thick, hh, rr1, rr2, ss, stiff)
subroutine heat_conductivity_shell_731(etype, nn, ecoord, tt, imat, thick, ss, stiff, ntab, temp, funca, funcb)
CALCULATION 4 NODE SHELL ELEMENT.
subroutine heat_get_area(xx, yy, zz, aa)
This module manages calculation relates with materials.
This modules defines a structure to record history dependent parameter in static analysis.