FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_EIG_setMASS.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
11 use m_fstr
12
13contains
14
15 subroutine setmass(fstrSOLID, hecMESH, hecMAT, fstrEIG)
16 use hecmw_util
17 use m_fstr
18 use mmaterial
20 implicit none
21 type(hecmwst_matrix) :: hecMAT
22 type(hecmwst_local_mesh) :: hecMESH
23 type(fstr_solid) :: fstrSOLID
24 type(fstr_eigen) :: fstrEIG
25 integer(kind=kint) :: i, iS, iE, ii, nn, j, jn, jS, k
26 integer(kind=kint) :: N, NP, NDOF
27 integer(kind=kint) :: icel, ic_type, itype, isect, ihead, sec_opt, cid
28 integer(kind=kint) :: nodLOCAL(20)
29 real(kind=kreal) :: val, surf, chkmass
30 real(kind=kreal) :: rho, thick, length
31 real(kind=kreal) :: lumped(120), mass(20*6, 20*6), ecoord(3,20)
32
33 if(myrank == 0)then
34 write(imsg,*)"* mass matrix generation"
35 endif
36
37 n = hecmat%N
38 np = hecmat%NP
39 ndof = hecmesh%n_dof
40
41 allocate(fstreig%mass(np*ndof))
42 fstreig%mass = 0.0d0
43 val = 0.0d0
44
45 do itype = 1, hecmesh%n_elem_type
46 is = hecmesh%elem_type_index(itype-1) + 1
47 ie = hecmesh%elem_type_index(itype )
48 ic_type = hecmesh%elem_type_item(itype)
49 nn = hecmw_get_max_node(ic_type)
50
51 if(hecmw_is_etype_link(ic_type)) cycle
52 if(hecmw_is_etype_patch(ic_type)) cycle
53 if(ic_type == 3414) cycle
54
55 do icel = is, ie
56 js = hecmesh%elem_node_index(icel-1)
57 do j = 1, nn
58 nodlocal(j) = hecmesh%elem_node_item(js+j)
59 do i = 1, 3
60 ecoord(i,j) = hecmesh%node(3*(nodlocal(j)-1)+i)
61 enddo
62 enddo
63
64 isect = hecmesh%section_ID(icel)
65 ihead = hecmesh%section%sect_R_index(isect-1)
66 cid = hecmesh%section%sect_mat_ID_item(isect)
67 sec_opt = hecmesh%section%sect_opt(isect)
68 rho = fstrsolid%materials(cid)%variables(m_density)
69 thick = fstrsolid%materials(cid)%variables(m_thick)
70
71 lumped = 0.0d0
72 if(ic_type == 231 .or. ic_type == 232 .or. ic_type == 241 .or. ic_type == 242)then
73 call mass_c2(ic_type, nn, ecoord(1:2,1:nn), fstrsolid%elements(icel)%gausses, sec_opt, thick, mass, lumped)
74
75 elseif(ic_type == 341 .or. ic_type == 342 .or. ic_type == 351 .or. ic_type == 352 .or. &
76 & ic_type == 361 .or. ic_type == 362 )then
77 call mass_c3(ic_type, nn, ecoord(1:3,1:nn), fstrsolid%elements(icel)%gausses, mass, lumped)
78
79 elseif(ic_type==731 .or. ic_type==741 .or. ic_type==743) then
80 rho = fstrsolid%materials(cid)%variables(m_density)
81 thick = fstrsolid%materials(cid)%variables(m_thick)
82 call mass_shell(ic_type, nn, ecoord(1:3,1:nn), rho, thick, fstrsolid%elements(icel)%gausses, mass, lumped)
83
84 elseif(ic_type == 761)then
85 surf = get_face3(ecoord(1:3,1:nn))
86 rho = fstrsolid%materials(cid)%variables(m_density)
87 thick = fstrsolid%materials(cid)%variables(m_thick)
88 val = surf*thick*rho/3.0d0
89
90 elseif(ic_type == 781)then
91 surf = get_face4(ecoord(1:3,1:nn))
92 rho = fstrsolid%materials(cid)%variables(m_density)
93 thick = fstrsolid%materials(cid)%variables(m_thick)
94 val = surf*thick*rho/4.0d0
95
96 elseif(ic_type == 611 .or. ic_type == 641)then
97 surf = hecmesh%section%sect_R_item(ihead+4)
98 length = get_length(ecoord(1:3,1:nn))
99 rho = fstrsolid%materials(cid)%variables(m_density)
100 val = 0.5d0*surf*length*rho
101
102 elseif(ic_type == 301)then
103 surf = hecmesh%section%sect_R_item(ihead+1)
104 length = get_length(ecoord(1:3,1:nn))
105 rho = fstrsolid%materials(cid)%variables(m_density)
106 val = 0.5d0*surf*length*rho
107
108 else
109 write(*,*)"** error setMASS"
110 endif
111
112 do j = 1,nn
113 jn = nodlocal(j)
114 js = ndof*(jn-1)
115
116 if(ic_type == 611)then
117 fstreig%mass(js+1) = fstreig%mass(js+1) + val
118 fstreig%mass(js+2) = fstreig%mass(js+2) + val
119 fstreig%mass(js+3) = fstreig%mass(js+3) + val
120 fstreig%mass(js+4) = fstreig%mass(js+4) + 0.0d0
121 fstreig%mass(js+5) = fstreig%mass(js+5) + 0.0d0
122 fstreig%mass(js+6) = fstreig%mass(js+6) + 0.0d0
123
124 elseif(ic_type == 641)then
125 if(j == 1 .or. j == 2)then
126 fstreig%mass(js+1) = fstreig%mass(js+1) + val
127 fstreig%mass(js+2) = fstreig%mass(js+2) + val
128 fstreig%mass(js+3) = fstreig%mass(js+3) + val
129 elseif(j == 3 .or. j == 4)then
130 fstreig%mass(js+1) = fstreig%mass(js+1) + 0.0d0
131 fstreig%mass(js+2) = fstreig%mass(js+2) + 0.0d0
132 fstreig%mass(js+3) = fstreig%mass(js+3) + 0.0d0
133 end if
134
135 elseif(ic_type == 761)then
136 if(j == 1 .or. j == 2 .or. j == 3)then
137 fstreig%mass(js+1) = fstreig%mass(js+1) + val
138 fstreig%mass(js+2) = fstreig%mass(js+2) + val
139 fstreig%mass(js+3) = fstreig%mass(js+3) + val
140 elseif(j == 4 .or. j == 5 .or. j == 6)then
141 fstreig%mass(js+1) = fstreig%mass(js+1) + 0.0d0
142 fstreig%mass(js+2) = fstreig%mass(js+2) + 0.0d0
143 fstreig%mass(js+3) = fstreig%mass(js+3) + 0.0d0
144 endif
145
146 else if( ic_type == 781 ) then
147 if(j == 1 .or. j == 2 .or. j == 3 .or. j == 4)then
148 fstreig%mass(js+1) = fstreig%mass(js+1) + val
149 fstreig%mass(js+2) = fstreig%mass(js+2) + val
150 fstreig%mass(js+3) = fstreig%mass(js+3) + val
151 elseif(j == 5 .or. j == 6 .or. j == 7 .or. j == 8)then
152 fstreig%mass(js+1) = fstreig%mass(js+1) + 0.0d0
153 fstreig%mass(js+2) = fstreig%mass(js+2) + 0.0d0
154 fstreig%mass(js+3) = fstreig%mass(js+3) + 0.0d0
155 endif
156
157 elseif(ic_type == 301) then
158 fstreig%mass(js+1) = fstreig%mass(js+1) + val
159 fstreig%mass(js+2) = fstreig%mass(js+2) + val
160 fstreig%mass(js+3) = fstreig%mass(js+3) + val
161
162 else
163 do k = 1, ndof
164 fstreig%mass(js+k) = fstreig%mass(js+k) + lumped(ndof*(j-1)+k)
165 enddo
166 endif
167 enddo
168 enddo
169 enddo
170
171 call hecmw_update_m_r(hecmesh, fstreig%mass, np, ndof)
172
173 chkmass = 0.0d0
174 do i = 1, n
175 ii = ndof*(i-1) + 1
176 chkmass = chkmass + fstreig%mass(ii)
177 end do
178 call hecmw_allreduce_r1(hecmesh, chkmass, hecmw_sum)
179 fstreig%totalmass = chkmass
180
181 if(myrank == 0)then
182 write(imsg,"(a,1pe12.5)")"** Total mass: ", chkmass
183 endif
184 end subroutine setmass
185
186end module m_fstr_eig_setmass
187
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=kint), parameter hecmw_sum
integer(kind=4), parameter kreal
This module contains subroutines used in 3d eigen analysis for.
Definition: dynamic_mass.f90:6
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)
This modules just summarizes all modules used in eigen analysis.
Definition: eigen_LIB.f90:6
Set up lumped mass matrix.
subroutine setmass(fstrsolid, hecmesh, hecmat, fstreig)
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint) myrank
PARALLEL EXECUTION.
Definition: m_fstr.f90:80
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
This module provide a function to fetch material properties from hecmw.
subroutine mass_shell(etype, nn, elem, rho, thick, gausses, mass, lumped)
This module summarizes all infomation of material properties.
Definition: material.f90:6
integer(kind=kint), parameter m_density
Definition: material.f90:86
integer(kind=kint), parameter m_thick
Definition: material.f90:87
Package of data used by Lanczos eigenvalue solver.
Definition: m_fstr.f90:562