FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
fstr_frequency_analysis.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  use m_fstr
8 
9  implicit none
10 contains
11  subroutine fstr_freq_result_init(hecMESH, numcomp, fstrRESULT)
12  !---- args
13  type(hecmwst_local_mesh), intent(in) :: hecMESH
14  integer(kind=kint), intent(in) :: numcomp
15  type(hecmwst_result_data), intent(inout) :: fstrRESULT
16  !---- vals
17  !---- body
18 
19  call hecmw_nullify_result_data(fstrresult)
20  fstrresult%ng_component = 0
21  fstrresult%nn_component = numcomp
22  fstrresult%ne_component = 0
23  allocate( fstrresult%nn_dof(numcomp) )
24  allocate( fstrresult%node_label(numcomp) )
25  allocate( fstrresult%node_val_item(numcomp*hecmesh%n_dof*hecmesh%n_node)) ! Should we use nn_internal?
26  end subroutine
27 
28  subroutine fstr_freq_result_add(fstrRESULT, hecMESH, comp_index, ndof, label, vect)
29  !---- args
30  type(hecmwst_result_data), intent(inout) :: fstrRESULT
31  type(hecmwst_local_mesh), intent(in) :: hecMESH
32  integer(kind=kint), intent(in) :: comp_index
33  integer(kind=kint), intent(in) :: ndof
34  character(len=HECMW_NAME_LEN), intent(in) :: label
35  real(kind=kreal), intent(in) :: vect(:)
36  !---- vals
37  integer(kind=kint) :: i, k, alldof, offset
38  !---- body
39 
40  fstrresult%nn_dof(comp_index) = ndof
41  fstrresult%node_label(comp_index) = label
42  alldof = fstrresult%nn_component*ndof
43  offset = ndof*(comp_index-1)
44  do i=1, hecmesh%n_node
45  do k=1, ndof
46  fstrresult%node_val_item(alldof*(i-1) + k + offset) = vect(ndof*(i-1) + k)
47  end do
48  end do
49 
50  end subroutine
51 
52 end module
53 
54 
56 
57  use m_fstr
59  use m_fstr_addbc
64 
65  implicit none
66 
67 contains
68 
69  subroutine fstr_solve_frequency_analysis(hecMESH, hecMAT, fstrSOLID, fstrEIG, &
70  fstrDYNAMIC, fstrRESULT, fstrPARAM, &
71  fstrCPL, fstrFREQ, fstrMAT, restart_step_num)
72  !C
73  !C-- global variable
74  !C
75  type(hecmwst_local_mesh) :: hecMESH
76  type(hecmwst_matrix) :: hecMAT
77  type(fstr_eigen) :: fstrEIG
78  type(fstr_solid) :: fstrSOLID
79  type(hecmwst_result_data) :: fstrRESULT
80  type(fstr_param) :: fstrPARAM
81  type(fstr_dynamic) :: fstrDYNAMIC
82  type(fstr_couple) :: fstrCPL
83  type(fstr_freqanalysis) :: fstrFREQ
84  type(fstrst_matrix_contact_lagrange) :: fstrMAT
85  integer(kind=kint) :: restart_step_num
86 
87  !C
88  !C-- local variable
89  !C
90  integer(kind=kint) :: numnode, numelm, startmode, endmode, nummode, ndof, im, in, ntotal, vistype
91  integer(kind=kint) :: numfreq, idnode, numdisp
92  integer(kind=kint) :: freqiout(3)
93  integer(kind=kint), parameter :: ilogin = 9056
94  real(kind=kreal), allocatable :: eigenvalue(:), loadvecre(:), loadvecim(:)
95  real(kind=kreal), allocatable :: bjre(:), bjim(:), dvare(:), dvaim(:), disp(:), vel(:), acc(:)
96  real(kind=kreal), allocatable :: dispre(:), dispim(:), velre(:), velim(:), accre(:), accim(:)
97  real(kind=kreal) :: freq, omega, val, dx, dy, dz, f_start, f_end
98  real(kind=kreal) :: t_start, t_end, time, dxi, dyi, dzi
99  type(fstr_freqanalysis_data) :: freqdata
100 
101  numnode = hecmesh%nn_internal
102  numelm = hecmesh%n_elem
103  ndof = hecmesh%n_dof
104  startmode = fstrfreq%start_mode
105  endmode = fstrfreq%end_mode
106  nummode = endmode - startmode +1
107 
108  freqdata%numMode = nummode
109  freqdata%numNodeDOF = numnode*ndof
110  allocate(freqdata%eigOmega(nummode))
111  allocate(freqdata%eigVector(numnode*ndof, nummode))
112 
113  call setupfreqparam(fstrdynamic, f_start, f_end, numfreq, freqdata%rayAlpha, freqdata%rayBeta, idnode, vistype, freqiout)
114  write(*,*) "Rayleigh alpha:", freqdata%rayAlpha
115  write(*,*) "Rayleigh beta:", freqdata%rayBeta
116  write(ilog,*) "Rayleigh alpha:", freqdata%rayAlpha
117  write(ilog,*) "Rayleigh beta:", freqdata%rayBeta
118 
119 
120  allocate(eigenvalue(nummode))
121  allocate(loadvecre(numnode*ndof))
122  allocate(loadvecim(numnode*ndof))
123  allocate(bjre(nummode))
124  allocate(bjim(nummode))
125  allocate(dvare(numnode*ndof))
126  allocate(dvaim(numnode*ndof))
127  allocate(disp(numnode*ndof))
128  allocate(vel(numnode*ndof))
129  allocate(acc(numnode*ndof))
130  allocate(dispre(numnode*ndof))
131  allocate(dispim(numnode*ndof))
132  allocate(velre(numnode*ndof))
133  allocate(velim(numnode*ndof))
134  allocate(accre(numnode*ndof))
135  allocate(accim(numnode*ndof))
136 
137  loadvecre(:) = 0.0d0
138  loadvecim(:) = 0.0d0
139 
140  write(*,*) "--frequency analysis--"
141  write(*, *) "read from=", trim(fstrfreq%eigenlog_filename)
142  write(ilog,*) "read from=", trim(fstrfreq%eigenlog_filename)
143  write(*, *) "start mode=", startmode
144  write(ilog,*) "start mode=", startmode
145  write(*, *) "end mode=", endmode
146  write(ilog,*) "end mode=", endmode
147  open(unit=ilogin, file=trim(fstrfreq%eigenlog_filename), status="OLD", action="READ")
148  call read_eigen_values(ilogin, startmode, endmode, eigenvalue, freqdata%eigOmega)
149  !call read_eigen_vector(ilogin, startmode, endmode, ndof, numnode, freqData%eigVector)
150  close(ilogin)
151 
152  call read_eigen_vector_res(hecmesh, startmode, endmode, ndof, numnode, freqdata%eigVector)
153 
154  call extract_surf2node(hecmesh, fstrfreq, ndof, loadvecre, loadvecim)
155  call assemble_nodeload(hecmesh, fstrfreq, ndof, loadvecre, loadvecim)
156 
157  write(*,*) "calc mass matrix"
158  call calcmassmatrix(fstrparam, hecmesh, hecmat, fstrsolid, fstreig, fstrmat)
159  write(*,*) "scale eigenvector"
160  call scaleeigenvector(fstreig, ndof*numnode, nummode, freqdata%eigVector)
161 
162  write(*, *) "start frequency:", f_start
163  write(ilog,*) "start frequency:", f_start
164  write(*, *) "end frequency:", f_end
165  write(ilog,*) "end frequency:", f_end
166  write(* ,*) "number of the sampling points", numfreq
167  write(ilog,*) "number of the sampling points", numfreq
168  write(* ,*) "monitor nodeid=", idnode
169  write(ilog,*) "monitor nodeid=", idnode
170 
171  do im=1, numfreq
172  freq = (f_end-f_start)/dble(numfreq)*dble(im) + f_start
173  omega = 2.0d0 * 3.14159265358979d0 * freq
174 
175  call calcfreqcoeff(freqdata, loadvecre, loadvecim, omega, bjre, bjim)
176  call calcdispvector(freqdata, bjre, bjim, dvare, dvaim)
177 
178  dx = sqrt(dvare(3*(idnode-1)+1)**2 + dvaim(3*(idnode-1)+1)**2)
179  dy = sqrt(dvare(3*(idnode-1)+2)**2 + dvaim(3*(idnode-1)+2)**2)
180  dz = sqrt(dvare(3*(idnode-1)+3)**2 + dvaim(3*(idnode-1)+3)**2)
181  val = sqrt(dx**2 + dy**2 + dz**2)
182  write(*, *) freq, "[Hz] : ", val
183  write(ilog, *) freq, "[Hz] : ", val
184  disp(:) = abs(cmplx(dvare(:), dvaim(:)))
185 
186  call calcvelvector(freqdata, omega, bjre, bjim, dvare, dvaim)
187  vel(:) = abs(cmplx(dvare(:), dvaim(:)))
188 
189  call calcaccvector(freqdata, omega, bjre, bjim, dvare, dvaim)
190  acc(:) = abs(cmplx(dvare(:), dvaim(:)))
191 
192  if(iresult==1) then
193  write(*, *) freq, "[Hz] : ", im, ".res"
194  write(ilog,*) freq, "[Hz] : ", im, ".res"
195  call output_resfile(hecmesh, freq, im, disp, vel, acc, freqiout)
196  end if
197  if(ivisual==1 .and. vistype==1) then
198  write(*, *) freq, "[Hz] : ", im, ".vis"
199  write(ilog,*) freq, "[Hz] : ", im, ".vis"
200  call output_visfile(hecmesh, im, disp, vel, acc, freqiout)
201  end if
202  end do
203 
204  call setupdynaparam(fstrdynamic, t_start, t_end, freq, numdisp)
205  write(*, *) "start time:", t_start
206  write(ilog,*) "start time:", t_start
207  write(*, *) "end time:", t_end
208  write(ilog,*) "end time:", t_end
209  write(*, *) "freqency:", freq
210  write(ilog,*) "freqency:", freq
211  write(*, *) "node id:", idnode
212  write(ilog,*) "node id:", idnode
213  write(*, *) "num disp:", numdisp
214  write(ilog,*) "num disp:", numdisp
215 
216  omega = 2.0d0 * 3.14159265358979d0 * freq
217  call calcfreqcoeff(freqdata, loadvecre, loadvecim, omega, bjre, bjim)
218  call calcdispvector(freqdata, bjre, bjim, dvare, dvaim)
219 
220  do im=1, numdisp
221  time = (t_end-t_start)/dble(numdisp)*dble(im-1) + t_start
222  call calcdispvectortime(freqdata, time, omega, bjre, bjim, dvare, dvaim)
223  dx = dvare(3*(idnode-1)+1)
224  dy = dvare(3*(idnode-1)+2)
225  dz = dvare(3*(idnode-1)+3)
226  dxi = dvaim(3*(idnode-1)+1)
227  dyi = dvaim(3*(idnode-1)+2)
228  dzi = dvaim(3*(idnode-1)+3)
229 
230  call calcvelvectortime(freqdata, time, omega, bjre, bjim, velre, velim)
231  call calcaccvectortime(freqdata, time, omega, bjre, bjim, accre, accim)
232  if(iresult==1) then
233  write(*, *) "time=", time, " : ", im, ".res"
234  write(ilog,*) "time=", time, " : ", im, ".res"
235  call outputdyna_resfile(hecmesh, time, im, dvare, dvaim, velre, velim, accre, accim, freqiout)
236  end if
237  if(ivisual==1 .and. vistype==2) then
238  write(*, *) "time=", time, " : ", im, ".vis"
239  write(ilog,*) "time=", time, " : ", im, ".vis"
240  call outputdyna_visfile(hecmesh, im, dvare, dvaim, velre, velim, accre, accim, freqiout)
241  end if
242  end do
243 
244  deallocate(freqdata%eigOmega)
245  deallocate(freqdata%eigVector)
246  deallocate(eigenvalue)
247  deallocate(loadvecre)
248  deallocate(loadvecim)
249  deallocate(bjre)
250  deallocate(bjim)
251  deallocate(dvare)
252  deallocate(dvaim)
253  deallocate(disp)
254  deallocate(vel)
255  deallocate(acc)
256  deallocate(dispre)
257  deallocate(dispim)
258  deallocate(velre)
259  deallocate(velim)
260  deallocate(accre)
261  deallocate(accim)
262 
263  end subroutine
264 
265  subroutine read_eigen_values(logfile, startmode, endmode, eigenvalue, anglfreq)
266  !---- args
267  integer(kind=kint), intent(in) :: logfile
268  integer(kind=kint), intent(in) :: startmode
269  integer(kind=kint), intent(in) :: endmode
270  real(kind=kreal), intent(inout) :: eigenvalue(:) !intend(endmode-startmode+1)
271  real(kind=kreal), intent(inout) :: anglfreq(:)
272  !---- vals
273  integer(kind=kint) :: im, endflag, id
274  character(len=HECMW_MSG_LEN) :: line
275  real(kind=kreal) :: freq
276  !---- body
277 
278  rewind(logfile)
279  endflag = 0
280  !Find eigenvalue header.
281  do
282  read(logfile, '(A80)', err=119) line
283  if(trim(adjustl(line)) == "NO. EIGENVALUE FREQUENCY (HZ) X Y Z X") then
284  endflag = 1
285  exit
286  end if
287  end do
288  read(logfile, '(A80)') line
289  !read eigenvalue
290  do im=1, startmode-1
291  read(logfile, '(A80)') line
292  end do
293  do im=1, (endmode-startmode+1)
294  read(logfile, '(i5,3e12.4,a)', err=119) id, eigenvalue(im), anglfreq(im), freq, line
295  end do
296  return
297 
298  !error handling
299  119 write(*,*) "Error to find eigenvalue information from logfile"
300  write(ilog,*) "Error to find eigenvalue information from logfile"
301  stop
302  end subroutine
303 
304  subroutine read_eigen_vector(logfile, startmode, endmode, numdof, numnode, eigenvector)
305  !---- args
306  integer(kind=kint), intent(in) :: logfile
307  integer(kind=kint), intent(in) :: startmode
308  integer(kind=kint), intent(in) :: endmode
309  integer(kind=kint), intent(in) :: numdof
310  integer(kind=kint), intent(in) :: numnode
311  real(kind=kreal), intent(inout) :: eigenvector(:, :) !intend (numdof*NN,nmode)
312  !---- vals
313  integer(kind=kint) :: im, in, gblid, j, idx
314  real(kind=kreal) :: vec(6)
315  character(len=HECMW_MSG_LEN) :: line
316  !---- body
317 
318  rewind(logfile)
319  !Find first eigenvector header
320  do im=1, startmode-1
321  do
322  read(logfile, '(a80)', err=119, end=119) line
323  if(line(1:9) == " Mode No.") then
324  exit ! goto next mode
325  end if
326  end do
327  end do
328 
329  !read eigenvector
330  do im=1, (endmode-startmode+1)
331  !find header
332  do
333  read(logfile, '(a80)', err=119, end=119) line
334  if(line(1:9) == " Mode No.") then
335  exit !find eigenmode
336  end if
337  end do
338 
339  read(logfile, '(a80)', err=119) line
340  read(logfile, '(a80)', err=119) line
341  read(logfile, '(a80)', err=119) line
342  read(logfile, '(a80)', err=119) line
343  !read eigenvector component
344  do in=1, numnode
345  select case(numdof)
346  case(2)
347  read(logfile, '(i10,2e12.4)', err=119) gblid, (vec(j), j=1,2)
348 
349  case(3)
350  !read(logfile, '(i10,3e12.4)', ERR=119) gblid, (vec(j), j=1,3)
351  read(logfile, '(i10,3e16.8)', err=119) gblid, (vec(j), j=1,3)
352 
353  case(6)
354  read(logfile, '(i10,6e12.4)', err=119) gblid, (vec(j), j=1,6)
355  case default
356  !error
357  goto 119
358  end select
359 
360  do j=1, numdof
361  idx = (in-1)*numdof + j
362  eigenvector(idx,im) = vec(j)
363  end do
364  end do
365  end do
366  return
367 
368  !error handling
369  119 write(*,*) "Error to find eigenvector from logfile"
370  write(ilog,*) "Error to find eigenvector from logfile"
371  stop
372 
373  end subroutine
374 
375  subroutine read_eigen_vector_res(hecMESH, startmode, endmode, numdof, numnode, eigenvector)
376  !---- args
377  type(hecmwst_local_mesh), intent(in) :: hecMESH
378  integer(kind=kint), intent(in) :: startmode
379  integer(kind=kint), intent(in) :: endmode
380  integer(kind=kint), intent(in) :: numdof
381  integer(kind=kint), intent(in) :: numnode
382  real(kind=kreal), intent(inout) :: eigenvector(:, :) !intend (numdof*NN,nmode)
383  !---- vals
384  integer(kind=kint), parameter :: compidx = 1 !Component index of displacement
385  integer(kind=kint) :: imode, idx, ind, a, b, nallcomp, j
386  type(hecmwst_result_data) :: eigenres
387  character(len=HECMW_NAME_LEN) :: name
388  !---- body
389 
390  name = 'result-in'
391  do imode=startmode, endmode
392  call nullify_result_data(eigenres)
393  call hecmw_result_read_by_name(hecmesh, name, imode, eigenres)
394 
395  nallcomp = 0
396  do ind=1,eigenres%nn_component
397  nallcomp = nallcomp + eigenres%nn_dof(ind)
398  end do
399 
400  idx = imode - startmode + 1
401  do ind=1, numnode
402  do j=1, numdof
403  a = (ind-1)*nallcomp + j !src vector index
404  b = (ind-1)*numdof + j
405  eigenvector(b,imode) = eigenres%node_val_item(a)
406  end do
407  end do
408  call free_result_data(eigenres)
409  end do
410 
411  contains
412 
413  subroutine free_result_data(res)
414  !---- args
415  type(hecmwst_result_data) :: res
416  !---- vals
417  !---- body
418  if(associated(res%nn_dof)) deallocate(res%nn_dof)
419  if(associated(res%ne_dof)) deallocate(res%ne_dof)
420  if(associated(res%node_label)) deallocate(res%node_label)
421  if(associated(res%elem_label)) deallocate(res%elem_label)
422  if(associated(res%node_val_item)) deallocate(res%node_val_item)
423  if(associated(res%elem_val_item)) deallocate(res%elem_val_item)
424  end subroutine
425 
426  subroutine nullify_result_data(res)
427  !---- args
428  type(hecmwst_result_data) :: res
429  !---- vals
430  !---- body
431  nullify(res%nn_dof)
432  nullify(res%ne_dof)
433  nullify(res%node_label)
434  nullify(res%elem_label)
435  nullify(res%node_val_item)
436  nullify(res%elem_val_item)
437  end subroutine
438 
439 
440  end subroutine
441 
442  subroutine output_resfile(hecMESH, freq, ifreq, disp, vel, acc, iout)
443  !---- args
444  type(hecmwst_local_mesh), intent(in) :: hecMESH
445  real(kind=kreal), intent(in) :: freq
446  integer(kind=kint), intent(in) :: ifreq
447  real(kind=kreal), intent(in) :: disp(:) !intend (numnodeDOF)
448  real(kind=kreal), intent(in) :: vel(:) !intend (numnodeDOF)
449  real(kind=kreal), intent(in) :: acc(:) !intend (numnodeDOF)
450  integer(kind=kint), intent(in) :: iout(3)
451  !---- vals
452  integer(kind=kint) :: im
453  character(len=HECMW_HEADER_LEN) :: header
454  character(len=HECMW_MSG_LEN) :: comment
455  character(len=HECMW_NAME_LEN) :: label, nameid
456  real(kind=kreal) :: freqval(1)
457  !---- body
458 
459  nameid='fstrRES'
460  header='*fstrresult'
461  comment='frequency_result'
462  call hecmw_result_init(hecmesh, ifreq, header, comment)
463 
464  label = "frequency"
465  freqval(1) = freq
466  call hecmw_result_add(3, 1, label, freqval)
467 
468  if(iout(1) == 1) then
469  label='displacement'
470  call hecmw_result_add(1, 3, label, disp) !mode=node, ndof=3
471  end if
472  if(iout(2) == 1) then
473  label='velocity'
474  call hecmw_result_add(1, 3, label, vel) !mode=node, ndof=3
475  end if
476  if(iout(3) == 1) then
477  label='acceleration'
478  call hecmw_result_add(1, 3, label, acc) !mode=node, ndof=3
479  end if
480  call hecmw_result_write_by_name(nameid)
481  call hecmw_result_finalize()
482  return
483  end subroutine
484 
485  subroutine output_visfile(hecMESH, ifreq, disp, vel, acc, iout)
486  !---- args
487  type(hecmwst_local_mesh), intent(in) :: hecmesh
488  integer(kind=kint), intent(in) :: ifreq
489  real(kind=kreal), intent(in) :: disp(:) !intend (numnodeDOF)
490  real(kind=kreal), intent(in) :: vel(:) !intend (numnodeDOF)
491  real(kind=kreal), intent(in) :: acc(:) !intend (numnodeDOF)
492  integer(kind=kint), intent(in) :: iout(3)
493  !---- vals
494  type(hecmwst_result_data) :: fstrRESULT
495  character(len=HECMW_NAME_LEN) :: label
496  integer(kind=kint) :: ncomp, i
497  !---- body
498  ncomp = 0
499  do i=1, 3
500  if(iout(i) == 1) then
501  ncomp = ncomp + 1
502  end if
503  end do
504 
505  call fstr_freq_result_init(hecmesh, ncomp, fstrresult)
506  ncomp=1
507  if(iout(1) == 1) then
508  label = 'displace_abs'
509  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, disp)
510  ncomp = ncomp + 1
511  end if
512  if(iout(2) == 1) then
513  label = 'velocity_abs'
514  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, vel)
515  ncomp = ncomp + 1
516  end if
517  if(iout(3) == 1) then
518  label = 'acceleration_abs'
519  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, acc)
520  ncomp = ncomp + 1
521  end if
522 
523  call fstr2hecmw_mesh_conv(hecmesh)
524  call hecmw_visualize_init
525  call hecmw_visualize( hecmesh, fstrresult, ifreq )
526  call hecmw2fstr_mesh_conv(hecmesh)
527  call hecmw_result_free(fstrresult)
528  end subroutine
529 
530  subroutine extract_surf2node(hecMESH, freqData, numdof, loadvecRe, loadvecIm)
531  !---- args
532  type(hecmwst_local_mesh), intent(in) :: hecMESH
533  type(fstr_freqanalysis), intent(in) :: freqData
534  integer(kind=kint), intent(in) :: numdof
535  real(kind=kreal), intent(inout) :: loadvecre(:) !intend(numnode*ndof)
536  real(kind=kreal), intent(inout) :: loadvecim(:) !intend(numnode*ndof)
537  !---- vals
538  integer(kind=kint), parameter :: MAXNODE = 100
539  integer(kind=kint) :: sgrpID, is, ie, ic, nsurf, ic_type, outtype, node_index(MAXNODE)
540  integer(kind=kint) :: nn, iss, nodeid, dof_index, ndof
541  integer(kind=kint) :: i, j, k, l, m, isn, nsize
542  integer(kind=kint) :: iwk(60), nodLOCAL(20)
543  real(kind=kreal) :: vect(60), xx(20), yy(20), zz(20), forcere(3), forceim(3)
544  !---- body
545 
546  ndof = 3
547  do i=1,freqdata%FLOAD_ngrp_tot
548  if(freqdata%FLOAD_ngrp_TYPE(i) == kfloadtype_surf) then !FLOAD type=surface
549  sgrpid = freqdata%FLOAD_ngrp_ID(i)
550  dof_index = freqdata%FLOAD_ngrp_DOF(i)
551  forcere(:) = 0.0d0
552  forceim(:) = 0.0d0
553  forcere(dof_index) = freqdata%FLOAD_ngrp_valre(i)
554  forceim(dof_index) = freqdata%FLOAD_ngrp_valim(i)
555 
556  is = hecmesh%surf_group%grp_index(sgrpid-1) + 1
557  ie = hecmesh%surf_group%grp_index(sgrpid)
558  do j=is, ie
559  ic = hecmesh%surf_group%grp_item(2*j-1)
560  nsurf = hecmesh%surf_group%grp_item(2*j)
561  ic_type = hecmesh%elem_type(ic)
562  nn = hecmw_get_max_node(ic_type)
563  isn = hecmesh%elem_node_index(ic-1)
564  do k=1, nn
565  nodlocal(k) = hecmesh%elem_node_item(isn+k)
566  xx(k) = hecmesh%node(3*nodlocal(k)-2)
567  yy(k) = hecmesh%node(3*nodlocal(k)-1)
568  zz(k) = hecmesh%node(3*nodlocal(k) )
569  do l=1, ndof
570  iwk(ndof*(k-1)+l) = ndof*(nodlocal(k)-1)+l
571  end do
572  end do
573 
574  call dl_c3_freq(ic_type, nn, xx, yy, zz, nsurf, forcere, vect, nsize)
575  do k=1,nsize
576  loadvecre(iwk(k)) = loadvecre(iwk(k)) + vect(k)
577  end do
578 
579  call dl_c3_freq(ic_type, nn, xx, yy, zz, nsurf, forceim, vect, nsize)
580  do k=1,nsize
581  loadvecim(iwk(k)) = loadvecim(iwk(k)) + vect(k)
582  end do
583  end do
584  end if
585  end do
586 
587  return
588  end subroutine
589 
590  subroutine dl_c3_freq(ETYPE, NN, XX, YY, ZZ, LTYPE, force, VECT, nsize)
591  !---- args
592  integer(kind=kint), intent(in) :: ETYPE !--solid element type
593  integer(kind=kint), intent(in) :: NN !--node num
594  integer(kind=kint), intent(in) :: LTYPE !--solid element face
595  real(kind=kreal), intent(in) :: xx(:) !--node x pos
596  real(kind=kreal), intent(in) :: yy(:) !--node y pos
597  real(kind=kreal), intent(in) :: zz(:) !--node z pos
598  real(kind=kreal), intent(in) :: force(3) !--node surfforce
599  real(kind=kreal), intent(inout) :: vect(:)
600  integer(kind=kint), intent(inout) :: nsize
601  !---- vals
602  integer(kind=kint), parameter :: NDOF = 3
603  real(kind=kreal) :: wg
604  integer(kind=kint) :: NOD(NN)
605  real(kind=kreal) :: elecoord(3, nn), localcoord(3)
606  real(kind=kreal) :: h(nn)
607  integer(kind=kint) :: I, IG2, NSUR, SURTYPE
608  !---- body
609 
610  call getsubface( etype, ltype, surtype, nod )
611  nsur = getnumberofnodes( surtype )
612 
613  do i=1,nsur
614  elecoord(1,i)=xx(nod(i))
615  elecoord(2,i)=yy(nod(i))
616  elecoord(3,i)=zz(nod(i))
617  end do
618  nsize = nn*ndof
619  vect(1:nsize) = 0.0d0
620  do ig2=1,numofquadpoints( surtype )
621  call getquadpoint( surtype, ig2, localcoord(1:2) )
622  call getshapefunc( surtype, localcoord(1:2), h(1:nsur) )
623 
624  wg=getweight( surtype, ig2 )
625  do i=1,nsur
626  vect(3*nod(i)-2)=vect(3*nod(i)-2)+wg*h(i)*force(1)
627  vect(3*nod(i)-1)=vect(3*nod(i)-1)+wg*h(i)*force(2)
628  vect(3*nod(i) )=vect(3*nod(i) )+wg*h(i)*force(3)
629  end do
630  end do
631  end subroutine
632 
633  subroutine assemble_nodeload(hecMESH, freqData, numdof, loadvecRe, loadvecIm)
634  !---- args
635  type(hecmwst_local_mesh), intent(in) :: hecMESH
636  type(fstr_freqanalysis), intent(in) :: freqData
637  integer(kind=kint), intent(in) :: numdof
638  real(kind=kreal), intent(inout) :: loadvecre(:)
639  real(kind=kreal), intent(inout) :: loadvecim(:)
640  !---- vals
641  integer(kind=kint) :: i, vecsize, ig, is, ie, in, nodeid, dof_index
642 
643  !---- body
644 
645  do i=1, freqdata%FLOAD_ngrp_tot
646  if(freqdata%FLOAD_ngrp_TYPE(i) == kfloadtype_node) then
647  ig = freqdata%FLOAD_ngrp_ID(i)
648  is = hecmesh%node_group%grp_index(ig-1) + 1
649  ie = hecmesh%node_group%grp_index(ig)
650  do in=is, ie
651  nodeid = hecmesh%node_group%grp_item(in)
652  dof_index = freqdata%FLOAD_ngrp_DOF(i)
653  loadvecre((nodeid-1)*numdof + dof_index) = loadvecre((nodeid-1)*numdof + dof_index) + freqdata%FLOAD_ngrp_valre(i)
654  loadvecim((nodeid-1)*numdof + dof_index) = loadvecim((nodeid-1)*numdof + dof_index) + freqdata%FLOAD_ngrp_valim(i)
655  end do
656  end if
657  end do
658  return
659  end subroutine
660 
661  subroutine calcmassmatrix(fstrPARAM, hecMESH, hecMAT, fstrSOLID, fstrEIG, fstrMAT)
662  !---- args
663  type(fstr_param), intent(in) :: fstrPARAM
664  type(hecmwst_local_mesh), intent(in) :: hecMESH
665  type(hecmwst_matrix), intent(inout) :: hecMAT
666  type(fstr_solid), intent(inout) :: fstrSOLID
667  type(fstr_eigen), intent(inout) :: fstrEIG
668  type(fstrst_matrix_contact_lagrange), intent(inout) :: fstrMAT
669  !---- vals
670  integer(kind=kint) :: ntotal
671  !---- body
672 
673 
674  fstrsolid%dunode = 0.d0
675  call fstr_stiffmatrix( hecmesh, hecmat, fstrsolid, 0.d0, 0.d0 )
676  call fstr_addbc(1, hecmesh, hecmat, fstrsolid, fstrparam, fstrmat, 2)
677 
678  call setmass(fstrsolid, hecmesh, hecmat, fstreig)
679 
680  end subroutine
681 
682  subroutine scaleeigenvector(fstrEIG, ntotaldof, nmode, eigenvector)
683  !---- args
684  type(fstr_eigen), intent(in) :: fstrEIG
685  integer(kind=kint), intent(in) :: ntotaldof
686  integer(kind=kint), intent(in) :: nmode
687  real(kind=kreal), intent(inout) :: eigenvector(:, :)
688  !---- vals
689  integer(kind=kint) :: imode, idof
690  real(kind=kreal) :: mas
691  !---- body
692 
693  do imode=1,nmode
694  mas = 0.0d0
695  do idof=1,ntotaldof
696  mas = mas + fstreig%mass(idof)*eigenvector(idof,imode)**2
697  end do
698  do idof=1,ntotaldof
699  eigenvector(idof,imode) = eigenvector(idof,imode) / sqrt(mas)
700  end do
701  end do
702  end subroutine
703 
704  subroutine checkorthvector(fstrEIG, eigenvector, imode, jmode, prod)
705  !---- args
706  type(fstr_eigen), intent(in) :: fstreig
707  real(kind=kreal), intent(in) :: eigenvector(:, :)
708  integer(kind=kint), intent(in) :: imode
709  integer(kind=kint), intent(in) :: jmode
710  real(kind=kreal), intent(inout) :: prod
711  !---- vals
712  integer(kind=kint) :: idof, s
713  !---- body
714  s = size(eigenvector(:,1))
715  prod = 0.0d0
716 
717  do idof=1, s
718  prod = prod + eigenvector(idof,imode)*fstreig%mass(idof)*eigenvector(idof,jmode)
719  end do
720  return
721  end subroutine
722 
723  subroutine writeoutvector(im, vector)
724  !---- args
725  integer(kind=kint), intent(in) :: im
726  real(kind=kreal), intent(in) :: vector(:)
727  !---- vals
728  integer(kind=kint) :: i, s
729  !---- body
730  s = size(vector)
731  do i=1, s
732  if(i == 1) then
733  write(*,'("eigenvec",i2.2,":[",e12.5", ")') im, vector(i)
734  else if(i /= s) then
735  write(*,'(e12.5,", ")') vector(i)
736  else
737  write(*,'(e12.5,"];")') vector(i)
738  end if
739  end do
740  write(*,*)
741  return
742  end subroutine
743 
744  subroutine calcdotproduct(a, b, c)
745  !---- args
746  real(kind=kreal), intent(in) :: a(:)
747  real(kind=kreal), intent(in) :: b(:)
748  real(kind=kreal), intent(inout) :: c
749  !---- vals
750  !---- body
751  c = dot_product(a, b)
752  !Next, we need allreduce operation to implement distribute mesh.
753  return
754  end subroutine
755 
756  subroutine calcfreqcoeff(freqData, loadRe, loadIm, inpOmega, bjRe, bjIm)
757  !---- args
758  type(fstr_freqanalysis_data), intent(in) :: freqdata
759  real(kind=kreal), intent(in) :: loadre(:) !intend (numNodeDOF)
760  real(kind=kreal), intent(in) :: loadim(:) !intend (numNodeDOF)
761  real(kind=kreal), intent(in) :: inpomega
762  real(kind=kreal), intent(inout) :: bjre(:) !intend (numMode)
763  real(kind=kreal), intent(inout) :: bjim(:) !intend (numMode)
764  !---- vals
765  integer(kind=kint) :: imode
766  real(kind=kreal) :: ujfr, ujfi, a, b, alp, beta
767  !---- body
768 
769  alp = freqdata%rayAlpha
770  beta = freqdata%rayBeta
771 
772  do imode=1, freqdata%numMode
773  call calcdotproduct(freqdata%eigVector(:,imode), loadre, ujfr)
774  call calcdotproduct(freqdata%eigVector(:,imode), loadim, ujfi)
775 
776  a = ujfr*(freqdata%eigOmega(imode)**2 - inpomega**2) + ujfi*(alp + beta*freqdata%eigOmega(imode)**2)*inpomega
777  b = (freqdata%eigOmega(imode)**2 - inpomega**2)**2 + ((alp + beta*freqdata%eigOmega(imode)**2)*inpomega)**2
778  bjre(imode) = a / b
779 
780  a = ujfi*(freqdata%eigOmega(imode)**2 -inpomega**2) - ujfr*(alp + beta*freqdata%eigOmega(imode)**2)*inpomega
781  b = (freqdata%eigOmega(imode)**2 - inpomega**2)**2 + ((alp + beta*freqdata%eigOmega(imode)**2)*inpomega)**2
782  bjim(imode) = a / b
783  end do
784  return
785  end subroutine
786 
787  subroutine calcdispvector(freqData, bjRe, bjIm, dispRe, dispIm)
788  !---- args
789  type(fstr_freqanalysis_data), intent(in) :: freqdata
790  real(kind=kreal), intent(in) :: bjre(:) !intend (numMode)
791  real(kind=kreal), intent(in) :: bjim(:) !intend (numMode)
792  real(kind=kreal), intent(inout) :: dispre(:) !intend (numNodeDOF)
793  real(kind=kreal), intent(inout) :: dispim(:) !intend (numNodeDOF)
794  !---- vals
795  integer(kind=kint) :: imode
796  !---- body
797 
798  dispre(:) = 0.0d0
799  dispim(:) = 0.0d0
800 
801  do imode=1, freqdata%numMode
802  dispre(:) = dispre(:) + bjre(imode)*freqdata%eigVector(:,imode)
803  dispim(:) = dispim(:) + bjim(imode)*freqdata%eigVector(:,imode)
804  end do
805  return
806  end subroutine
807 
808  subroutine calcvelvector(freqData, omega, bjRe, bjIm, velRe, velIm)
809  !---- args
810  type(fstr_freqanalysis_data), intent(in) :: freqData
811  real(kind=kreal), intent(in) :: omega
812  real(kind=kreal), intent(in) :: bjre(:) !intend (numMode)
813  real(kind=kreal), intent(in) :: bjim(:) !intend (numMode)
814  real(kind=kreal), intent(inout) :: velre(:) !intend (numNodeDOF)
815  real(kind=kreal), intent(inout) :: velim(:) !intend (numNodeDOF)
816  !---- vals
817  integer(kind=kint) :: imode
818  !---- body
819 
820  velre(:) = 0.0d0
821  velim(:) = 0.0d0
822 
823  do imode=1, freqdata%numMode
824  velre(:) = velre(:) - omega * bjim(imode) * freqdata%eigVector(:,imode)
825  velim(:) = velim(:) + omega * bjre(imode) * freqdata%eigVector(:,imode)
826  end do
827  end subroutine
828 
829  subroutine calcaccvector(freqData, omega, bjRe, bjIm, accRe, accIm)
830  !---- args
831  type(fstr_freqanalysis_data), intent(in) :: freqData
832  real(kind=kreal), intent(in) :: omega
833  real(kind=kreal), intent(in) :: bjre(:) !intend (numMode)
834  real(kind=kreal), intent(in) :: bjim(:) !intend (numMode)
835  real(kind=kreal), intent(inout) :: accre(:) !intend (numNodeDOF)
836  real(kind=kreal), intent(inout) :: accim(:) !intend (numNodeDOF)
837  !---- vals
838  integer(kind=kint) :: imode
839  !---- body
840 
841  accre(:) = 0.0d0
842  accim(:) = 0.0d0
843 
844  do imode=1, freqdata%numMode
845  accre(:) = accre(:) - omega**2 * bjre(imode) * freqdata%eigVector(:,imode)
846  accim(:) = accim(:) - omega**2 * bjim(imode) * freqdata%eigVector(:,imode)
847  end do
848 
849  end subroutine
850 
851  subroutine setupfreqparam(fstrDYNAMIC, f_start, f_end, numfreq, raym, rayk, idnode, vistype, ioutl)
852  !---- args
853  type(fstr_dynamic), intent(in) :: fstrDYNAMIC
854  real(kind=kreal), intent(inout) :: f_start
855  real(kind=kreal), intent(inout) :: f_end
856  integer(kind=kint), intent(inout) :: numfreq
857  real(kind=kreal), intent(inout) :: raym
858  real(kind=kreal), intent(inout) :: rayk
859  integer(kind=kint), intent(inout) :: idnode
860  integer(kind=kint), intent(inout) :: vistype
861  integer(kind=kint), intent(inout) :: ioutl(3)
862  !---- vals
863 
864  !---- body
865  f_start = fstrdynamic%t_start
866  f_end = fstrdynamic%t_end
867  numfreq = fstrdynamic%n_step
868  raym = fstrdynamic%ray_m
869  rayk = fstrdynamic%ray_k
870  idnode = fstrdynamic%nout_monit
871  vistype = fstrdynamic%ngrp_monit
872  ioutl(1:3) = fstrdynamic%iout_list(1:3)
873  return
874  end subroutine
875 
876  subroutine calcdispvectortime(freqData, time, omega, bjRe, bjIm, dispRe, dispIm)
877  !---- args
878  type(fstr_freqanalysis_data), intent(in) :: freqData
879  real(kind=kreal), intent(in) :: time
880  real(kind=kreal), intent(in) :: omega
881  real(kind=kreal), intent(in) :: bjre(:) !intend (numMode)
882  real(kind=kreal), intent(in) :: bjim(:) !intend (numMode)
883  real(kind=kreal), intent(inout) :: dispre(:) !intend (numNodeDOF)
884  real(kind=kreal), intent(inout) :: dispim(:) !intend (numNodeDOF)
885  !---- vals
886  integer(kind=kint) :: imode, idf, s
887  complex(kind=kreal) :: a, b, c
888  !---- body
889 
890  dispre(:) = 0.0d0
891  dispim(:) = 0.0d0
892  a = exp(cmplx(0.0d0, omega*time))
893 
894  do imode=1, freqdata%numMode
895  s = size(freqdata%eigvector(:,imode))
896  b = cmplx(bjre(imode), bjim(imode)) * a
897  do idf=1, s
898  c = b*cmplx(freqdata%eigVector(idf,imode), 0.0d0)
899  dispre(idf) = dispre(idf) + dble(c)
900  dispim(idf) = dispim(idf) + imag(c)
901  end do
902  end do
903  return
904  end subroutine
905 
906  subroutine calcvelvectortime(freqData, time, omega, bjRe, bjIm, velRe, velIm)
907  !---- args
908  type(fstr_freqanalysis_data), intent(in) :: freqData
909  real(kind=kreal), intent(in) :: time
910  real(kind=kreal), intent(in) :: omega
911  real(kind=kreal), intent(in) :: bjre(:) !intend (numMode)
912  real(kind=kreal), intent(in) :: bjim(:) !intend (numMode)
913  real(kind=kreal), intent(inout) :: velre(:) !intend (numNodeDOF)
914  real(kind=kreal), intent(inout) :: velim(:) !intend (numNodeDOF)
915  !---- vals
916  integer(kind=kint) :: imode, idf, s
917  complex(kind=kreal) :: a, b, c
918  !---- body
919 
920  velre(:) = 0.0d0
921  velim(:) = 0.0d0
922  a = cmplx(0.0d0, 1.0d0)*cmplx(omega, 0.0d0)*exp(cmplx(0.0d0, omega*time))
923 
924  do imode=1, freqdata%numMode
925  s = size(freqdata%eigvector(:,imode))
926  b = cmplx(bjre(imode), bjim(imode)) * a
927  do idf=1, s
928  c = b*cmplx(freqdata%eigVector(idf,imode), 0.0d0)
929  velre(idf) = velre(idf) + dble(c)
930  velim(idf) = velim(idf) + imag(c)
931  end do
932  end do
933  return
934  end subroutine
935 
936  subroutine calcaccvectortime(freqData, time, omega, bjRe, bjIm, accRe, accIm)
937  !---- args
938  type(fstr_freqanalysis_data), intent(in) :: freqData
939  real(kind=kreal), intent(in) :: time
940  real(kind=kreal), intent(in) :: omega
941  real(kind=kreal), intent(in) :: bjre(:) !intend (numMode)
942  real(kind=kreal), intent(in) :: bjim(:) !intend (numMode)
943  real(kind=kreal), intent(inout) :: accre(:) !intend (numNodeDOF)
944  real(kind=kreal), intent(inout) :: accim(:) !intend (numNodeDOF)
945  !---- vals
946  integer(kind=kint) :: imode, idf, s
947  complex(kind=kreal) :: a, b, c
948  !---- body
949 
950  accre(:) = 0.0d0
951  accim(:) = 0.0d0
952  a = cmplx(-1.0d0, 0.0d0)*cmplx(omega**2, 0.0d0)*exp(cmplx(0.0d0, omega*time))
953 
954  do imode=1, freqdata%numMode
955  s = size(freqdata%eigvector(:,imode))
956  b = cmplx(bjre(imode), bjim(imode)) * a
957  do idf=1, s
958  c = b*cmplx(freqdata%eigVector(idf,imode), 0.0d0)
959  accre(idf) = accre(idf) + dble(c)
960  accim(idf) = accim(idf) + imag(c)
961  end do
962  end do
963  return
964  end subroutine
965 
966  subroutine setupdynaparam(fstrDYNAMIC, t_start, t_end, dynafreq, numdisp)
967  !---- args
968  type(fstr_dynamic), intent(in) :: fstrDYNAMIC
969  real(kind=kreal), intent(inout) :: t_start
970  real(kind=kreal), intent(inout) :: t_end
971  real(kind=kreal), intent(inout) :: dynafreq
972  integer(kind=kint), intent(inout) :: numdisp
973  !---- vals
974  !---- body
975  t_start = fstrdynamic%ganma
976  t_end = fstrdynamic%beta
977  dynafreq = fstrdynamic%t_delta
978  numdisp = fstrdynamic%nout
979  return
980  end subroutine
981 
982  subroutine outputdyna_resfile(hecMESH, time, istp, dispre, dispim, velre, velim, accre, accim, iout)
983  !---- args
984  type(hecmwst_local_mesh), intent(in) :: hecMESH
985  real(kind=kreal), intent(in) :: time
986  integer(kind=kint), intent(in) :: istp
987  real(kind=kreal), intent(in) :: dispre(:) !intend (numnodeDOF)
988  real(kind=kreal), intent(in) :: dispim(:) !intend (numnodeDOF)
989  real(kind=kreal), intent(in) :: velre(:) !intend (numnodeDOF)
990  real(kind=kreal), intent(in) :: velim(:) !intend (numnodeDOF)
991  real(kind=kreal), intent(in) :: accre(:) !intend (numnodeDOF)
992  real(kind=kreal), intent(in) :: accim(:) !intend (numnodeDOF)
993  integer(kind=kint), intent(in) :: iout(3)
994  !---- vals
995  integer(kind=kint) :: im, s
996  character(len=HECMW_HEADER_LEN) :: header
997  character(len=HECMW_MSG_LEN) :: comment
998  character(len=HECMW_NAME_LEN) :: label, nameid
999  real(kind=kreal), allocatable :: absval(:)
1000  !---- body
1001 
1002  s = size(dispre)
1003  allocate(absval(s))
1004 
1005  nameid='fstrDYNA'
1006  header='*fstrresult'
1007  comment='frequency_result'
1008 
1009  call hecmw_result_init(hecmesh, istp, header, comment)
1010 
1011  label = "time"
1012  absval(1) = time
1013  call hecmw_result_add(3, 1, label, absval)
1014 
1015  if(iout(1) == 1) then
1016  label='displacement_real'
1017  call hecmw_result_add(1, 3, label, dispre) !mode=node, ndof=3
1018  label='displacement_imag'
1019  call hecmw_result_add(1, 3, label, dispim)
1020  label='displacement_abs'
1021  absval(:) = abs(cmplx(dispre(:), dispim(:)))
1022  call hecmw_result_add(1, 3, label, absval)
1023  end if
1024 
1025  if(iout(2) == 1) then
1026  label='velocity_real'
1027  call hecmw_result_add(1, 3, label, velre) !mode=node, ndof=3
1028  label='velocity_imag'
1029  call hecmw_result_add(1, 3, label, velim)
1030  label='velocity_abs'
1031  absval(:) = abs(cmplx(velre(:), velim(:)))
1032  call hecmw_result_add(1, 3, label, absval)
1033  end if
1034 
1035  if(iout(3) == 1) then
1036  label='acceleration_real'
1037  call hecmw_result_add(1, 3, label, accre) !mode=node, ndof=3
1038  label='acceleration_imag'
1039  call hecmw_result_add(1, 3, label, accim)
1040  label='acceleration_abs'
1041  absval(:) = abs(cmplx(velre(:), velim(:)))
1042  call hecmw_result_add(1, 3, label, absval)
1043  end if
1044 
1045  call hecmw_result_write_by_name(nameid)
1046  call hecmw_result_finalize()
1047 
1048  deallocate(absval)
1049  return
1050  end subroutine
1051 
1052  subroutine outputdyna_visfile(hecMESH, istp, dispre, dispim, velre, velim, accre, accim, iout)
1053  !---- args
1054  type(hecmwst_local_mesh), intent(inout) :: hecmesh
1055  integer(kind=kint), intent(in) :: istp
1056  real(kind=kreal), intent(in) :: dispre(:)
1057  real(kind=kreal), intent(in) :: dispim(:)
1058  real(kind=kreal), intent(in) :: velre(:)
1059  real(kind=kreal), intent(in) :: velim(:)
1060  real(kind=kreal), intent(in) :: accre(:)
1061  real(kind=kreal), intent(in) :: accim(:)
1062  integer(kind=kint), intent(in) :: iout(3)
1063  !---- vals
1064  type(hecmwst_result_data) :: fstrRESULT
1065  character(len=HECMW_NAME_LEN) :: label
1066  integer(kind=kint) :: s, ncomp, i
1067  real(kind=kreal), allocatable :: absval(:)
1068  !---- body
1069 
1070  s = size(dispre)
1071  allocate(absval(s))
1072 
1073  ncomp = 0
1074  do i=1, 3
1075  if(iout(i) == 1) then
1076  ncomp = ncomp + 3 !re, im, abs
1077  end if
1078  end do
1079 
1080  call fstr_freq_result_init(hecmesh, ncomp, fstrresult) !disp, vel, acc
1081 
1082  ncomp = 1
1083 
1084  if(iout(1) == 1) then
1085  label = 'displace_real'
1086  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, dispre)
1087  ncomp = ncomp + 1
1088 
1089  label = 'displace_imag'
1090  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, dispim)
1091  ncomp = ncomp + 1
1092 
1093  label = 'displace_abs'
1094  absval(:) = abs(cmplx(dispre(:), dispim(:)))
1095  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, absval)
1096  ncomp = ncomp + 1
1097  end if
1098 
1099  if(iout(2) == 1) then
1100  label = 'velocity_real'
1101  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, velre)
1102  ncomp = ncomp + 1
1103 
1104  label = 'velocity_imag'
1105  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, velim)
1106  ncomp = ncomp + 1
1107 
1108  label = 'velocity_abs'
1109  absval(:) = abs(cmplx(velre(:), velim(:)))
1110  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, absval)
1111  ncomp = ncomp + 1
1112  end if
1113 
1114  if(iout(3) == 1) then
1115  label = 'acceleration_real'
1116  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, accre)
1117  ncomp = ncomp + 1
1118 
1119  label = 'acceleration_imag'
1120  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, accim)
1121  ncomp = ncomp + 1
1122 
1123  label = 'acceleration_abs'
1124  absval(:) = abs(cmplx(accre(:), accim(:)))
1125  call fstr_freq_result_add(fstrresult, hecmesh, ncomp, 3, label, absval)
1126  ncomp = ncomp + 1
1127  end if
1128 
1129  call fstr2hecmw_mesh_conv(hecmesh)
1130  call hecmw_visualize_init
1131  call hecmw_visualize( hecmesh, fstrresult, istp )
1132  call hecmw2fstr_mesh_conv(hecmesh)
1133  call hecmw_result_free(fstrresult)
1134 
1135  deallocate(absval)
1136  end subroutine
1137 
1138 end module
subroutine nullify_result_data(res)
subroutine free_result_data(res)
subroutine fstr_solve_frequency_analysis(hecMESH, hecMAT, fstrSOLID, fstrEIG, fstrDYNAMIC, fstrRESULT, fstrPARAM, fstrCPL, fstrFREQ, fstrMAT, restart_step_num)
subroutine calcdispvector(freqData, bjRe, bjIm, dispRe, dispIm)
subroutine read_eigen_vector(logfile, startmode, endmode, numdof, numnode, eigenvector)
subroutine output_resfile(hecMESH, freq, ifreq, disp, vel, acc, iout)
subroutine read_eigen_vector_res(hecMESH, startmode, endmode, numdof, numnode, eigenvector)
subroutine calcmassmatrix(fstrPARAM, hecMESH, hecMAT, fstrSOLID, fstrEIG, fstrMAT)
subroutine calcaccvector(freqData, omega, bjRe, bjIm, accRe, accIm)
subroutine extract_surf2node(hecMESH, freqData, numdof, loadvecRe, loadvecIm)
subroutine writeoutvector(im, vector)
subroutine calcdispvectortime(freqData, time, omega, bjRe, bjIm, dispRe, dispIm)
subroutine checkorthvector(fstrEIG, eigenvector, imode, jmode, prod)
subroutine calcfreqcoeff(freqData, loadRe, loadIm, inpOmega, bjRe, bjIm)
subroutine scaleeigenvector(fstrEIG, ntotaldof, nmode, eigenvector)
subroutine output_visfile(hecMESH, ifreq, disp, vel, acc, iout)
subroutine setupdynaparam(fstrDYNAMIC, t_start, t_end, dynafreq, numdisp)
subroutine calcvelvectortime(freqData, time, omega, bjRe, bjIm, velRe, velIm)
subroutine setupfreqparam(fstrDYNAMIC, f_start, f_end, numfreq, raym, rayk, idnode, vistype, ioutl)
subroutine outputdyna_resfile(hecMESH, time, istp, dispre, dispim, velre, velim, accre, accim, iout)
subroutine read_eigen_values(logfile, startmode, endmode, eigenvalue, anglfreq)
subroutine assemble_nodeload(hecMESH, freqData, numdof, loadvecRe, loadvecIm)
subroutine calcaccvectortime(freqData, time, omega, bjRe, bjIm, accRe, accIm)
subroutine outputdyna_visfile(hecMESH, istp, dispre, dispim, velre, velim, accre, accim, iout)
subroutine dl_c3_freq(ETYPE, NN, XX, YY, ZZ, LTYPE, force, VECT, nsize)
subroutine calcvelvector(freqData, omega, bjRe, bjIm, velRe, velIm)
This module contains steady state frequency analysis.
subroutine fstr_freq_result_init(hecMESH, numcomp, fstrRESULT)
subroutine fstr_freq_result_add(fstrRESULT, hecMESH, comp_index, ndof, label, vect)
This module provides functions of reconstructing.
This module provides a function to deal with prescribed displacement.
Definition: fstr_AddBC.f90:7
Set up lumped mass matrix.
This module provides function to calcualte tangent stiffness matrix.
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), 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.
Structure for Lagrange multiplier-related part of stiffness matrix (Lagrange multiplier-related matri...
Data for coupling analysis.
Definition: m_fstr.f90:580
Data for DYNAMIC ANSLYSIS (fstrDYNAMIC)
Definition: m_fstr.f90:473
Package of data used by Lanczos eigenvalue solver.
Definition: m_fstr.f90:562
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138