FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
fstr_solve_dynamic.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 
8 
9  use m_fstr
13  use fstr_frequency_analysis !Frequency analysis module
14 
15 contains
16 
17  !C================================================================C
19  !C================================================================C
20  subroutine fstr_solve_dynamic(hecMESH,hecMAT,fstrSOLID,fstrEIG &
21  ,fstrDYNAMIC,fstrRESULT,fstrPARAM &
22  ,fstrCPL,fstrFREQ, fstrMAT &
23  ,conMAT )
24  use m_fstr_setup
25  implicit none
26  type(hecmwst_local_mesh) :: hecMESH
27  type(hecmwst_matrix) :: hecMAT
28  type(fstr_eigen) :: fstrEIG
29  type(fstr_solid) :: fstrSOLID
30  type(hecmwst_result_data) :: fstrRESULT
31  type(fstr_param) :: fstrPARAM
32  type(fstr_dynamic) :: fstrDYNAMIC
33  type(fstr_couple) :: fstrCPL !for COUPLE
34  type(fstr_freqanalysis) :: fstrFREQ
35  type(fstrst_matrix_contact_lagrange) :: fstrMAT
36  type(fstr_info_contactchange) :: infoCTChange
37  type(hecmwst_matrix), optional :: conMAT
38  integer(kind=kint) :: i, j, num_monit, ig, is, iE, ik, in, ing, iunitS, iunit, ierror, flag, limit
39  character(len=HECMW_FILENAME_LEN) :: fname, header
40  integer(kind=kint) :: restrt_step_num, ndof
41  integer(kind=kint) :: restrt_step(1)
42 
43  num_monit = 0
44 
45  if(dabs(fstrdynamic%t_delta) < 1.0e-20) then
46  if( hecmesh%my_rank == 0 ) then
47  write(imsg,*) 'stop due to fstrDYNAMIC%t_delta = 0'
48  end if
49  call hecmw_abort( hecmw_comm_get_comm())
50  end if
51  call fstr_dynamic_alloc( hecmesh, fstrdynamic )
52 
53  !C-- file open for local use
54  !C
55  if(fstrdynamic%idx_resp == 1) then ! time history analysis
56  call hecmw_ctrl_is_subdir( flag, limit )
57  if( flag == 0)then
58  header = ""
59  else
60  header = adjustl("MONITOR/")
61  call hecmw_ctrl_make_subdir( adjustl("MONITOR/test.txt"), limit )
62  endif
63  ig = fstrdynamic%ngrp_monit
64  is = hecmesh%node_group%grp_index(ig-1)+1
65  ie = hecmesh%node_group%grp_index(ig)
66  do ik=is,ie
67  in = hecmesh%node_group%grp_item(ik)
68  if (in > hecmesh%nn_internal) cycle
69  num_monit = num_monit+1
70  ing = hecmesh%global_node_id(in)
71  iunits = 10*(num_monit-1)
72 
73  iunit = iunits + fstrdynamic%dynamic_IW4
74  write(fname,'(a,i0,a)') trim(header)//'dyna_disp_',ing,'.txt'
75  if(fstrdynamic%restart_nout < 0 ) then
76  open(iunit,file=fname, position = 'append', iostat=ierror)
77  else
78  open(iunit,file=fname, status = 'replace', iostat=ierror)
79  endif
80  if( ierror /= 0 ) then
81  write(*,*) 'stop due to file opening error:',trim(fname)
82  call hecmw_abort( hecmw_comm_get_comm())
83  end if
84 
85  iunit = iunits + fstrdynamic%dynamic_IW5
86  write(fname,'(a,i0,a)') trim(header)//'dyna_velo_',ing,'.txt'
87  if(fstrdynamic%restart_nout < 0 ) then
88  open(iunit,file=fname, position = 'append', iostat=ierror)
89  else
90  open(iunit,file=fname, status = 'replace', iostat=ierror)
91  endif
92  if( ierror /= 0 ) then
93  write(*,*) 'stop due to file opening error',trim(fname)
94  call hecmw_abort( hecmw_comm_get_comm())
95  end if
96 
97  iunit = iunits + fstrdynamic%dynamic_IW6
98  write(fname,'(a,i0,a)') trim(header)//'dyna_acce_',ing,'.txt'
99  if(fstrdynamic%restart_nout < 0 ) then
100  open(iunit,file=fname, position = 'append', iostat=ierror)
101  else
102  open(iunit,file=fname, status = 'replace', iostat=ierror)
103  endif
104  if( ierror /= 0 ) then
105  write(*,*) 'stop due to file opening error',trim(fname)
106  call hecmw_abort( hecmw_comm_get_comm())
107  end if
108 
109  iunit = iunits + fstrdynamic%dynamic_IW7
110  write(fname,'(a,i0,a)') trim(header)//'dyna_force_',ing,'.txt'
111  if(fstrdynamic%restart_nout < 0 ) then
112  open(iunit,file=fname, position = 'append', iostat=ierror)
113  else
114  open(iunit,file=fname, status = 'replace', iostat=ierror)
115  endif
116  if( ierror /= 0 ) then
117  write(*,*) 'stop due to file opening error',trim(fname)
118  call hecmw_abort( hecmw_comm_get_comm())
119  end if
120  iunit = iunits + fstrdynamic%dynamic_IW8
121  write(fname,'(a,i0,a)') trim(header)//'dyna_strain_',ing,'.txt'
122  if(fstrdynamic%restart_nout < 0 ) then
123  open(iunit,file=fname, position = 'append', iostat=ierror)
124  else
125  open(iunit,file=fname, status = 'replace', iostat=ierror)
126  endif
127  if( ierror /= 0 ) then
128  write(*,*) 'stop due to file opening error',trim(fname)
129  call hecmw_abort( hecmw_comm_get_comm())
130  end if
131 
132  iunit = iunits + fstrdynamic%dynamic_IW9
133  write(fname,'(a,i0,a)') trim(header)//'dyna_stress_',ing,'.txt'
134  if(fstrdynamic%restart_nout < 0 ) then
135  open(iunit,file=fname, position = 'append', iostat=ierror)
136  else
137  open(iunit,file=fname, status = 'replace', iostat=ierror)
138  endif
139  if( ierror /= 0 ) then
140  write(*,*) 'stop due to file opening error',trim(fname)
141  call hecmw_abort( hecmw_comm_get_comm())
142  end if
143  enddo
144  endif
145 
146  !C
147  !C-- initial condition
148  !C
149  fstrdynamic%DISP = 0.d0
150  fstrdynamic%VEL = 0.d0
151  fstrdynamic%ACC = 0.d0
152  fstrsolid%unode(:) =0.d0
153  fstrsolid%QFORCE(:) =0.d0
154  fstrdynamic%kineticEnergy=0.d0
155  fstrdynamic%strainEnergy=0.d0
156  fstrdynamic%totalEnergy=0.d0
157  call fstr_updatestate( hecmesh, fstrsolid, 0.d0 )
158 
159  ! ---- restart
160 
161  restrt_step_num = 1
162  fstrdynamic%i_step = 0
163  infoctchange%contactNode_previous = 0
164 
165  if(associated(g_initialcnd))then
166  ndof = hecmat%NDOF
167  do j = 1, size(g_initialcnd)
168  if(g_initialcnd(j)%cond_name == "velocity")then
169  do i= 1, hecmesh%n_node
170  ing = g_initialcnd(j)%intval(i)
171  if(ing <= 0) cycle
172  fstrdynamic%VEL(ndof*i-(ndof-ing),1) = g_initialcnd(j)%realval(i)
173  end do
174  elseif(g_initialcnd(j)%cond_name == "acceleration")then
175  do i = 1, hecmesh%n_node
176  ing = g_initialcnd(j)%intval(i)
177  if(ing <= 0) cycle
178  fstrdynamic%ACC(ndof*i-(ndof-ing),1) = g_initialcnd(j)%realval(i)
179  enddo
180  endif
181  enddo
182  endif
183 
184  if(fstrdynamic%restart_nout >= 0 ) then
185  call dynamic_bc_init (hecmesh, hecmat, fstrsolid, fstrdynamic)
186  call dynamic_bc_init_vl(hecmesh, hecmat, fstrsolid, fstrdynamic)
187  call dynamic_bc_init_ac(hecmesh, hecmat, fstrsolid, fstrdynamic)
188  endif
189 
190  !restart
191  if(fstrdynamic%restart_nout < 0 ) then
192  if( .not. associated( fstrsolid%contacts ) ) then
193  call fstr_read_restart_dyna_nl(restrt_step_num,hecmesh,fstrsolid,fstrdynamic,fstrparam)
194  else
195  call fstr_read_restart_dyna_nl(restrt_step_num,hecmesh,fstrsolid,fstrdynamic,fstrparam,&
196  infoctchange%contactNode_previous)
197  endif
198  restrt_step_num = restrt_step_num + 1
199  fstrdynamic%restart_nout = - fstrdynamic%restart_nout
200  hecmat%Iarray(98) = 1
201  end if
202 
203  if(fstrdynamic%idx_resp == 1) then ! time history analysis
204 
205  if(fstrdynamic%idx_eqa == 1) then ! implicit dynamic analysis
206  if(.not. associated( fstrsolid%contacts ) ) then
207  call fstr_solve_dynamic_nlimplicit(1, hecmesh,hecmat,fstrsolid,fstreig &
208  ,fstrdynamic,fstrresult,fstrparam &
209  ,fstrcpl, restrt_step_num )
210  elseif( fstrparam%contact_algo == kcaslagrange ) then
211  ! ---- For Parallel Contact with Multi-Partition Domains
212  if(paracontactflag.and.present(conmat)) then
213  call fstr_solve_dynamic_nlimplicit_contactslag(1, hecmesh,hecmat,fstrsolid,fstreig &
214  ,fstrdynamic,fstrresult,fstrparam &
215  ,fstrcpl,fstrmat,restrt_step_num,infoctchange &
216  ,conmat )
217  else
218  call fstr_solve_dynamic_nlimplicit_contactslag(1, hecmesh,hecmat,fstrsolid,fstreig &
219  ,fstrdynamic,fstrresult,fstrparam &
220  ,fstrcpl,fstrmat,restrt_step_num,infoctchange )
221  endif
222  endif
223 
224  else if(fstrdynamic%idx_eqa == 11) then ! explicit dynamic analysis
225  call fstr_solve_dynamic_nlexplicit(hecmesh,hecmat,fstrsolid,fstreig &
226  ,fstrdynamic,fstrresult,fstrparam,infoctchange &
227  ,fstrcpl, restrt_step_num )
228  endif
229 
230  else if(fstrdynamic%idx_resp == 2) then ! frequency response analysis
231 
232  if( fstrparam%nlgeom ) then
233  if( hecmesh%my_rank .eq. 0 ) then
234  write(imsg,*) 'stop: steady-state harmonic response analysis is not yet available !'
235  end if
236  call hecmw_abort( hecmw_comm_get_comm())
237  end if
238 
239  if( hecmesh%my_rank .eq. 0 ) then
240  call fstr_solve_frequency_analysis(hecmesh, hecmat, fstrsolid, fstreig, fstrdynamic, &
241  fstrresult, fstrparam, fstrcpl, fstrfreq, fstrmat, &
242  restrt_step_num)
243  end if
244  end if
245 
246  !C-- file close for local use
247  if(fstrdynamic%idx_resp == 1) then ! time history analysis
248  do i=1,num_monit
249  iunits = 10*(i-1)
250  close(iunits + fstrdynamic%dynamic_IW4)
251  close(iunits + fstrdynamic%dynamic_IW5)
252  close(iunits + fstrdynamic%dynamic_IW6)
253  close(iunits + fstrdynamic%dynamic_IW7)
254  close(iunits + fstrdynamic%dynamic_IW8)
255  close(iunits + fstrdynamic%dynamic_IW9)
256  enddo
257  endif
258  !C-- end of finalization
259 
260  call fstr_dynamic_finalize( fstrdynamic )
261  call hecmat_finalize( hecmat )
262 
263  deallocate( fstreig%mass )
264 
265  end subroutine fstr_solve_dynamic
266 
267 end module fstr_solver_dynamic
void hecmw_ctrl_is_subdir(int *flag, int *limit)
void hecmw_ctrl_make_subdir(char *filename, int *err, int len)
This module contains subroutines for nonlinear explicit dynamic analysis.
subroutine fstr_solve_dynamic_nlexplicit(hecMESH, hecMAT, fstrSOLID, fstrEIG, fstrDYN, fstrRESULT, fstrPARAM, infoCTChange, fstrCPL, restrt_step_num)
This module contains subroutines for nonlinear implicit dynamic analysis.
subroutine fstr_solve_dynamic_nlimplicit_contactslag(cstep, hecMESH, hecMAT, fstrSOLID, fstrEIG, fstrDYNAMIC, fstrRESULT, fstrPARAM, fstrCPL, fstrMAT, restrt_step_num, infoCTChange, conMAT)
This subroutine provides function of nonlinear implicit dynamic analysis using the Newmark method....
subroutine fstr_solve_dynamic_nlimplicit(cstep, hecMESH, hecMAT, fstrSOLID, fstrEIG, fstrDYNAMIC, fstrRESULT, fstrPARAM, fstrCPL, restrt_step_num)
This subroutine provides function of nonlinear implicit dynamic analysis using the Newmark method.
subroutine fstr_solve_frequency_analysis(hecMESH, hecMAT, fstrSOLID, fstrEIG, fstrDYNAMIC, fstrRESULT, fstrPARAM, fstrCPL, fstrFREQ, fstrMAT, restart_step_num)
This module provides functions of reconstructing.
This module contains subroutines controlling dynamic calculation.
subroutine fstr_solve_dynamic(hecMESH, hecMAT, fstrSOLID, fstrEIG, fstrDYNAMIC, fstrRESULT, fstrPARAM, fstrCPL, fstrFREQ, fstrMAT, conMAT)
Master subroutine for dynamic analysis.
This module provides functions to read in data from control file and do neccessary preparation for fo...
Definition: fstr_setup.f90:7
subroutine fstr_dynamic_alloc(hecMESH, fstrDYNAMIC)
Initial setting of dynamic calculation.
subroutine fstr_dynamic_finalize(fstrDYNAMIC)
Finalizer of fstr_solid.
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
subroutine hecmat_finalize(hecMAT)
Definition: m_fstr.f90:891
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
integer(kind=kint), parameter kcaslagrange
contact analysis algorithm
Definition: m_fstr.f90:53
type(tinitialcondition), dimension(:), pointer, save g_initialcnd
Definition: m_fstr.f90:135
logical paracontactflag
PARALLEL CONTACT FLAG.
Definition: m_fstr.f90:84
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