1!$Id:$ 2 subroutine tint3d(ll,lint,s) 3 4! * * F E A P * * A Finite Element Analysis Program 5 6!.... Copyright (c) 1984-2017: Regents of the University of California 7! All rights reserved 8 9!-----[--.----+----.----+----.-----------------------------------------] 10! Purpose: Gauss quadrature for 3-d tetrahedral element 11 12! Inputs: 13! ll - Type of quadrature 14 15! Outputs: 16! lint - Number of quadrature points 17! s(5,*) - Values of volume coordinates and weights 18!-----[--.----+----.----+----.-----------------------------------------] 19 implicit none 20 21 integer i, j, ll, lint 22 real*8 s(5,*) 23 24 save 25 26! 1 pt. quadrature O(h^2) 27 28 if(ll.eq.1) then 29 lint = 1 30 do i = 1,4 31 s(i,1) = 0.25d0 32 end do ! i 33 s(5,1) = 1.0d0/6.d0 34 35! 4 pt. quadrature O(h^3) 36 37 elseif(ll.eq.2) then 38 lint = 4 39 s(5,4) = 0.25d0/6.d0 40 do i = 1,4 41 do j = 1,4 42 s(i,j) = 0.1381966011250105d+00 43 end do ! j 44 s(i,i) = 0.5854101966249658d+00 45 s(5,i) = s(5,4) 46 end do ! i 47 48! 11 pt. quadrature O(h^4) 49 50 elseif(ll.eq.3) then 51 lint = 11 52 do i = 1,3 53 do j = 1,10 54 s(i,j) = 0.0d0 55 end do ! j 56 s(i, i ) = 1.00d0 57 s(i, i+4) = 0.50d0 58 s(i, i+7) = 0.50d0 59 s(i, 11 ) = 0.25d0 60 end do ! i 61 s(2, 5) = 0.50d0 62 s(3, 6) = 0.50d0 63 s(1,10) = 0.50d0 64 do j = 1,4 65 s(5,j) = 1.d0/360.d0 66 end do ! j 67 do j = 5,10 68 s(5,j) = 1.d0/90.d0 69 end do ! j 70 s(5,11) = 4.d0/45.d0 71 else 72 73! 16 pt. quadrature O(h^5) 74 75 lint = 16 76 s(5,4) = 0.8395632516687135d-02 77 do i = 1,3 78 do j = 1,4 79 s(i,j) = 0.7611903264425430d-01 80 end do ! j 81 s(i,i) = 0.7716429020672371d+00 82 s(5,i) = s(5,4) 83 end do ! i 84 do i = 5,16 85 s(5,i) = 0.1109034477221540d-01 86 end do ! i 87 88 s(1, 5) = 0.1197005277978019d+00 89 s(2, 5) = 0.7183164526766925d-01 90 s(3, 5) = 0.4042339134672644d+00 91 s(1, 6) = 0.4042339134672644d+00 92 s(2, 6) = 0.1197005277978019d+00 93 s(3, 6) = 0.7183164526766925d-01 94 s(1, 7) = 0.4042339134672644d+00 95 s(2, 7) = 0.4042339134672644d+00 96 s(3, 7) = 0.1197005277978019d+00 97 s(1, 8) = 0.7183164526766925d-01 98 s(2, 8) = 0.4042339134672644d+00 99 s(3, 8) = 0.4042339134672644d+00 100 101 s(1, 9) = 0.1197005277978019d+00 102 s(2, 9) = 0.4042339134672644d+00 103 s(3, 9) = 0.7183164526766925d-01 104 s(1,10) = 0.4042339134672644d+00 105 s(2,10) = 0.1197005277978019d+00 106 s(3,10) = 0.4042339134672644d+00 107 s(1,11) = 0.7183164526766925d-01 108 s(2,11) = 0.4042339134672644d+00 109 s(3,11) = 0.1197005277978019d+00 110 s(1,12) = 0.4042339134672644d+00 111 s(2,12) = 0.7183164526766925d-01 112 s(3,12) = 0.4042339134672644d+00 113 114 s(1,13) = 0.1197005277978019d+00 115 s(2,13) = 0.4042339134672644d+00 116 s(3,13) = 0.4042339134672644d+00 117 s(1,14) = 0.7183164526766925d-01 118 s(2,14) = 0.1197005277978019d+00 119 s(3,14) = 0.4042339134672644d+00 120 s(1,15) = 0.4042339134672644d+00 121 s(2,15) = 0.7183164526766925d-01 122 s(3,15) = 0.1197005277978019d+00 123 s(1,16) = 0.4042339134672644d+00 124 s(2,16) = 0.4042339134672644d+00 125 s(3,16) = 0.7183164526766925d-01 126 127 endif 128 129! Compute fourth points 130 131 do j = 1,lint 132 s(4,j) = 1.d0 - (s(1,j) + s(2,j) + s(3,j)) 133 end do ! j 134 135 end 136