10 subroutine mass_c2(etype, nn, ecoord, gausses, sec_opt, thick, mass, lumped, temperature)
16 integer(kind=kint),
intent(in) :: etype
17 integer(kind=kint),
intent(in) :: nn
18 real(kind=kreal),
intent(in) :: ecoord(2,nn)
19 real(kind=kreal),
intent(out) :: mass(:,:)
20 real(kind=kreal),
intent(out) :: lumped(:)
21 real(kind=kreal),
intent(in),
optional :: temperature(nn)
22 type(tmaterial),
pointer :: matl
23 integer(kind=kint),
parameter :: ndof = 2
24 integer(kind=kint) :: i, j, lx, sec_opt
25 real(kind=kreal) :: naturalcoord(2)
26 real(kind=kreal) :: func(nn), thick
27 real(kind=kreal) :: det, wg, rho
28 real(kind=kreal) :: d(2,2), n(2, nn*ndof), dn(2, nn*ndof)
29 real(kind=kreal) :: gderiv(nn,2)
34 matl => gausses(1)%pMaterial
36 if(sec_opt == 2) thick = 1.0d0
43 if(
present(temperature))
then
47 rho = matl%variables(m_density)
54 rho = matl%variables(m_density)
72 dn(1:2, 1:nn*ndof) = matmul(d, n(1:2, 1:nn*ndof))
73 forall(i = 1:nn*ndof, j = 1:nn*ndof)
74 mass(i,j) = mass(i,j) + dot_product(n(:,i), dn(:,j))*wg
82 subroutine mass_c3(etype, nn, ecoord, gausses, mass, lumped, temperature)
88 integer(kind=kint),
intent(in) :: etype
89 integer(kind=kint),
intent(in) :: nn
90 real(kind=kreal),
intent(in) :: ecoord(3,nn)
91 real(kind=kreal),
intent(out) :: mass(:,:)
92 real(kind=kreal),
intent(out) :: lumped(:)
93 real(kind=kreal),
intent(in),
optional :: temperature(nn)
94 type(tmaterial),
pointer :: matl
95 integer(kind=kint),
parameter :: ndof = 3
96 integer(kind=kint) :: i, j, lx
97 real(kind=kreal) :: naturalcoord(3)
98 real(kind=kreal) :: func(nn)
99 real(kind=kreal) :: det, wg, rho
100 real(kind=kreal) :: d(3, 3), n(3, nn*ndof), dn(3, nn*ndof)
101 real(kind=kreal) :: gderiv(nn, 3)
106 matl => gausses(1)%pMaterial
113 if(
present(temperature))
then
117 rho = matl%variables(m_density)
124 rho = matl%variables(m_density)
144 dn(1:3, 1:nn*ndof) = matmul(d, n(1:3, 1:nn*ndof))
145 forall(i = 1:nn*ndof, j = 1:nn*ndof)
146 mass(i,j) = mass(i,j) + dot_product(n(:,i), dn(:,j))*wg
157 integer(kind=kint) :: i, j, nn, ndof
158 real(kind=kreal) :: lumped(:), mass(:,:)
159 real(kind=kreal) :: diag_mass, total_mass
162 do i = 1, nn*ndof, ndof
163 do j = 1, nn*ndof, ndof
164 total_mass = total_mass + mass(j,i)
169 do i = 1, nn*ndof, ndof
170 diag_mass = diag_mass + mass(i,i)
173 diag_mass = 1.0d0/diag_mass
175 lumped(i) = lumped(i) + mass(i,i)*total_mass*diag_mass
180 mass(i,i) = lumped(i)
190 (ecoord(1,2) - ecoord(1,1))**2 + &
191 (ecoord(2,2) - ecoord(2,1))**2 + &
192 (ecoord(3,2) - ecoord(3,1))**2 )
198 real(kind=kreal) ::
get_face3, ecoord(3,20)
199 real(kind=kreal) :: a1, a2, a3
200 real(kind=kreal) :: x(3), y(3), z(3)
202 x(1) = ecoord(1,1); y(1) = ecoord(2,1); z(1) = ecoord(3,1)
203 x(2) = ecoord(1,2); y(2) = ecoord(2,2); z(2) = ecoord(3,2)
204 x(3) = ecoord(1,3); y(3) = ecoord(2,3); z(3) = ecoord(3,3)
206 a1 = (x(2) - x(1))**2 + (y(2) - y(1))**2 + (z(2) - z(1))**2
207 a2 = (x(1) - x(3))*(x(2) - x(1)) &
208 & + (y(1) - y(3))*(y(2) - y(1)) &
209 & + (z(1) - z(3))*(z(2) - z(1))
210 a3 = (x(3) - x(1))**2 + (y(3) - y(1))**2 + (z(3) - z(1))**2
218 integer(kind=kint) :: lx, ly
219 real(kind=kreal) ::
get_face4, ecoord(3,20)
220 real(kind=kreal) :: xg(2), ri, si, rp, sp, rm, sm, hr(4), hs(4)
221 real(kind=kreal) :: xr, xs, yr, ys, zr, zs
222 real(kind=kreal) :: x(4), y(4), z(4), det
224 x(1) = ecoord(1,1); y(1) = ecoord(2,1); z(1) = ecoord(3,1)
225 x(2) = ecoord(1,2); y(2) = ecoord(2,2); z(2) = ecoord(3,2)
226 x(3) = ecoord(1,3); y(3) = ecoord(2,3); z(3) = ecoord(3,3)
227 x(4) = ecoord(1,4); y(4) = ecoord(2,4); z(4) = ecoord(3,4)
229 xg(1) = -0.5773502691896258d0
255 xr = hr(1)*x(1) + hr(2)*x(2) + hr(3)*x(3) + hr(4)*x(4)
256 xs = hs(1)*x(1) + hs(2)*x(2) + hs(3)*x(3) + hs(4)*x(4)
257 yr = hr(1)*y(1) + hr(2)*y(2) + hr(3)*y(3) + hr(4)*y(4)
258 ys = hs(1)*y(1) + hs(2)*y(2) + hs(3)*y(3) + hs(4)*y(4)
259 zr = hr(1)*z(1) + hr(2)*z(2) + hr(3)*z(3) + hr(4)*z(4)
260 zs = hs(1)*z(1) + hs(2)*z(2) + hs(3)*z(3) + hs(4)*z(4)
262 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.
real(kind=kreal) function get_length(ecoord)
subroutine mass_c2(etype, nn, ecoord, gausses, sec_opt, thick, mass, lumped, temperature)
real(kind=kreal) function get_face4(ecoord)
real(kind=kreal) function get_face3(ecoord)
subroutine mass_c3(etype, nn, ecoord, gausses, mass, lumped, temperature)
subroutine get_lumped_mass(nn, ndof, mass, lumped)
This module manages calculation relates with materials.
This modules defines a structure to record history dependent parameter in static analysis.
All data should be recorded in every quadrature points.