FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_EIG_output.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!-------------------------------------------------------------------------------
6
7contains
8
9 subroutine fstr_eigen_output(hecMESH, hecMAT, fstrEIG)
10 use m_fstr
13 implicit none
14 type(hecmwst_local_mesh) :: hecMESH
15 type(hecmwst_matrix) :: hecMAT
16 type(fstr_eigen) :: fstrEIG
17
18 integer(kind=kint) :: N, NP, NDOF, NNDOF, NPNDOF
19 integer(kind=kint) :: i, j, k, in, jn, nget
20 real(kind=kreal) :: chk, gm
21 real(kind=kreal), allocatable :: s(:), t(:), u(:), r(:)
22 real(kind=kreal), pointer :: mass(:), eigval(:), eigvec(:,:)
23
24 n = hecmat%N
25 np = hecmat%NP
26 ndof = hecmesh%n_dof
27 nndof = n *ndof
28 npndof = np*ndof
29
30 nget = fstreig%nget
31 mass => fstreig%mass
32 eigval => fstreig%eigval
33 eigvec => fstreig%eigvec
34
35 allocate(fstreig%effmass(ndof*nget))
36 allocate(fstreig%partfactor(ndof*nget))
37 allocate(s(npndof))
38 allocate(t(npndof))
39 allocate(u(npndof))
40 allocate(r(ndof))
41 fstreig%effmass = 0.0d0
42 fstreig%partfactor = 0.0d0
43
44 do i = 1, nget
45 r = 0.0d0
46 gm = 0.0d0
47 do j = 1, n
48 do k = 1, ndof
49 in = ndof*(j-1) + k
50 r(k) = r(k) + mass(in)*eigvec(in,i)
51 gm = gm + mass(in)*eigvec(in,i)*eigvec(in,i)
52 enddo
53 enddo
54
55 call hecmw_allreduce_r(hecmesh, r, ndof, hecmw_sum)
56 call hecmw_allreduce_r1(hecmesh, gm, hecmw_sum)
57
58 gm = 1.0d0/gm
59 do j = 1, ndof
60 in = ndof*(i-1) + j
61 fstreig%partfactor(in) = gm*r(j)
62 fstreig%effmass(in) = gm*r(j)*r(j)
63 enddo
64 enddo
65
66 call eglist(hecmesh, hecmat, fstreig)
67
68 if(myrank == 0)then
69 write(imsg,*) ''
70 write(imsg,*) '*----------------------------------------------*'
71 write(imsg,*) '*--E I G E N V A L U E C O N V E R G E N C E--*'
72 write(imsg,*) 'Absolute residual = |(||Kx - lambda*Mx||)|'
73 write(imsg,*) ' Iter.# Eigenvalue Abs. Residual '
74 write(ilog,*) ' Iter.# Eigenvalue Abs. Residual '
75 write(imsg,*) ' *-----* *---------* *--------------*'
76 endif
77
78 do j = 1, fstreig%nget
79 do i = 1, nndof
80 u(i) = eigvec(i,j)
81 enddo
82 call hecmw_matvec(hecmesh, hecmat, u, t)
83
84 s = 0.0d0
85 do i = 1, nndof
86 s(i) = mass(i)*eigvec(i,j)
87 enddo
88
89 chk = 0.0d0
90 do i = 1,nndof
91 chk = chk + (t(i) - (eigval(j)-fstreig%sigma)*s(i))**2
92 enddo
93 call hecmw_allreduce_r1(hecmesh, chk, hecmw_sum)
94 chk = dsqrt(chk)
95
96 if(myrank == 0)then
97 write(imsg,'(2x,i5,2x,1p5e15.6)') j, eigval(j), chk
98 write(ilog,'(i5,1p5e12.4)') j, eigval(j), chk
99 endif
100 enddo
101
102 if(myrank == 0)then
103 write(imsg,*)'* ---END Eigenvalue listing--- *'
104 endif
105
106 deallocate(s)
107 deallocate(t)
108 deallocate(u)
109 deallocate(r)
110 end subroutine fstr_eigen_output
111
112 subroutine fstr_eigen_make_result(hecMESH, hecMAT, fstrEIG, fstrRESULT)
113 use m_fstr
115 use hecmw_util
116 implicit none
117 type(hecmwst_local_mesh) :: hecmesh
118 type(hecmwst_matrix) :: hecMAT
119 type(fstr_eigen) :: fstrEIG
120 type(hecmwst_result_data) :: fstrRESULT
121
122 integer(kind=kint) :: i, istep, nget, NP, NDOF, NPNDOF, totalmpc, MPC_METHOD
123 real(kind=kreal) :: t1
124 real(kind=kreal), pointer :: eigvec(:,:)
125 real(kind=kreal), allocatable :: x(:), egval(:)
126 character(len=HECMW_HEADER_LEN) :: header
127 character(len=HECMW_MSG_LEN) :: comment
128 character(len=HECMW_NAME_LEN) :: label
129 character(len=HECMW_NAME_LEN) :: nameID
130
131 nget = fstreig%nget
132 ndof = hecmat%NDOF
133 npndof = hecmat%NP*hecmat%NDOF
134 !totalmpc = hecMESH%mpc%n_mpc
135 !call hecmw_allreduce_I1 (hecMESH, totalmpc, hecmw_sum)
136
137 eigvec => fstreig%eigvec
138
139 allocate(x(npndof))
140 x = 0.0d0
141 allocate(egval(1))
142
143 do istep = 1, nget
144 egval(1) = fstreig%eigval(istep)
145 do i=1,npndof
146 x(i) = eigvec(i,istep)
147 enddo
148
149 !if (totalmpc > 0) then
150 ! MPC_METHOD = hecmw_mat_get_mpc_method(hecMAT)
151 ! if (MPC_METHOD < 1 .or. 3 < MPC_METHOD) MPC_METHOD = 3
152 ! if (MPC_METHOD == 3) then ! elimination
153 ! call hecmw_tback_x_33(hecMESH, X, t1)
154 ! else
155 ! if (hecMESH%my_rank.eq.0) write(0,*) "### ERROR: MPC_METHOD must set to 3"
156 ! stop
157 ! endif
158 !endif
159
160 call hecmw_update_m_r(hecmesh, x, hecmat%NP, ndof)
161
162 if( iresult.eq.1 ) then
163 header = "*fstrresult"
164 comment = "eigen_result"
165 call hecmw_result_init(hecmesh,istep,header,comment)
166 label = "EIGENVALUE"
167 call hecmw_result_add(3,1,label,egval)
168 label = "DISPLACEMENT"
169 call hecmw_result_add(1,ndof,label,x)
170 nameid = "fstrRES"
171 call hecmw_result_write_by_name(nameid)
172 call hecmw_result_finalize
173 endif
174
175 if( ivisual.eq.1 ) then
176 call hecmw_nullify_result_data(fstrresult)
177 fstrresult%ng_component = 1
178 fstrresult%nn_component = 1
179 fstrresult%ne_component = 0
180 allocate(fstrresult%ng_dof(1))
181 allocate(fstrresult%global_label(1))
182 allocate(fstrresult%global_val_item(1))
183 fstrresult%ng_dof(1) = 1
184 fstrresult%global_label(1) = 'EIGENVALUE'
185 fstrresult%global_val_item(1) = egval(1)
186 allocate(fstrresult%nn_dof(1))
187 allocate(fstrresult%node_label(1))
188 allocate(fstrresult%node_val_item(ndof*hecmat%NP))
189 fstrresult%nn_dof(1) = ndof
190 fstrresult%node_label(1) = 'DISPLACEMENT'
191 fstrresult%node_val_item = x
192 call fstr2hecmw_mesh_conv(hecmesh)
193 call hecmw_visualize_init
194 call hecmw_visualize( hecmesh, fstrresult, istep )
195 call hecmw_visualize_finalize
196 call hecmw2fstr_mesh_conv(hecmesh)
197 call hecmw_result_free(fstrresult)
198 endif
199 enddo
200
201 deallocate(x)
202
203 end subroutine fstr_eigen_make_result
204
206 subroutine eglist(hecMESH, hecMAT, fstrEIG)
207 use m_fstr
208 use hecmw_util
210 implicit none
211 type(hecmwst_local_mesh) :: hecMESH
212 type(hecmwst_matrix) :: hecMAT
213 type(fstr_eigen) :: fstrEIG
214
215 integer(kind=kint) :: NDOF
216 integer(kind=kint) :: i, j, in, iter ,nget
217 real(kind=kreal) :: pi, angle, freq, pf(3), em(3)
218 real(kind=kreal), pointer :: eigval(:)
219
220 ndof = hecmat%NDOF
221 nget = fstreig%nget
222 iter = fstreig%iter
223 pi = 4.0d0 * datan(1.0d0)
224 eigval => fstreig%eigval
225
226 if(myrank == 0)then
227 write(ilog,*)""
228 write(ilog,"(a)")"********************************"
229 write(ilog,"(a)")"*RESULT OF EIGEN VALUE ANALYSIS*"
230 write(ilog,"(a)")"********************************"
231 write(ilog,"(a)")""
232 write(ilog,"(a,i8)")"NUMBER OF ITERATIONS = ",iter
233 write(ilog,"(a,1pe12.4)")"TOTAL MASS = ",fstreig%totalmass
234 write(ilog,"(a)")""
235 write(ilog,"(3a)")" ANGLE FREQUENCY ",&
236 "PARTICIPATION FACTOR EFFECTIVE MASS"
237 write(ilog,"(3a)")" NO. EIGENVALUE FREQUENCY (HZ) ",&
238 "X Y Z X Y Z"
239 write(ilog,"(3a)")" --- ---------- ---------- ---------- ",&
240 "---------- ---------- ---------- ---------- ---------- ----------"
241 write(*,*)""
242 write(*,"(a)")"#----------------------------------#"
243 write(*,"(a)")"# RESULT OF EIGEN VALUE ANALYSIS #"
244 write(*,"(a)")"#----------------------------------#"
245 write(*,"(a)")""
246 write(*,"(a,i8)")"### NUMBER OF ITERATIONS = ",iter
247 write(*,"(a,1pe12.4)")"### TOTAL MASS = ",fstreig%totalmass
248 write(*,"(a)")""
249 write(*,"(3a)")" PERIOD FREQUENCY ",&
250 "PARTICIPATION FACTOR EFFECTIVE MASS"
251 write(*,"(3a)")" NO. [Sec] [HZ] ",&
252 "X Y Z X Y Z"
253 write(*,"(3a)")" --- --------- --------- ",&
254 "--------- --------- --------- --------- --------- ---------"
255
256 in = 0
257 do i = 1, nget
258 in = in + 1
259 if(eigval(i) < 0.0d0) eigval(i) = 0.0d0
260 angle = dsqrt(eigval(i))
261 freq = angle*0.5d0/pi
262
263 pf = 0.0d0
264 em = 0.0d0
265 do j = 1, min(ndof, 3)
266 pf(j) = fstreig%partfactor(ndof*(i-1) + j)
267 em(j) = fstreig%effmass(ndof*(i-1) + j)
268 enddo
269
270 write(ilog,'(I5,1P9E12.4)') in, eigval(i), angle, freq, pf(1), pf(2), pf(3), em(1), em(2), em(3)
271 write(* ,'(I5,1P8E11.3)') in, 1.0d0/freq, freq, pf(1), pf(2), pf(3), em(1), em(2), em(3)
272 enddo
273 write(ilog,*)""
274 write(*,*)""
275 endif
276
277 end subroutine eglist
278end module m_fstr_eig_output
279
280
subroutine, public hecmw_matvec(hecmesh, hecmat, x, y, commtime)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
subroutine fstr_eigen_output(hecmesh, hecmat, fstreig)
subroutine fstr_eigen_make_result(hecmesh, hecmat, fstreig, fstrresult)
subroutine eglist(hecmesh, hecmat, fstreig)
Output eigenvalues and vectors.
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint), pointer iresult
Definition: m_fstr.f90:106
integer(kind=kint) myrank
PARALLEL EXECUTION.
Definition: m_fstr.f90:80
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
integer(kind=kint), parameter ilog
FILE HANDLER.
Definition: m_fstr.f90:91
integer(kind=kint), pointer ivisual
Definition: m_fstr.f90:107
HECMW to FSTR Mesh Data Converter. Convering Conectivity of Element Type 232, 342 and 352.
subroutine hecmw2fstr_mesh_conv(hecmesh)
subroutine fstr2hecmw_mesh_conv(hecmesh)
Package of data used by Lanczos eigenvalue solver.
Definition: m_fstr.f90:562