FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
fstr_Restart.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 
9  use m_utilities
10  use m_fstr
11  implicit none
12 
13 contains
14 
16  !----------------------------------------------------------------------*
17  subroutine fstr_read_restart(cstep,substep,step_count,ctime,dtime,hecMESH,fstrSOLID,fstrPARAM,contactNode)
18  !----------------------------------------------------------------------*
19  integer, intent(out) :: cstep
20  integer, intent(out) :: substep
21  integer, intent(out) :: step_count
22  real(kind=kreal), intent(out) :: ctime
23  real(kind=kreal), intent(out) :: dtime
24  integer, intent(out) :: contactNode
25  type (hecmwST_local_mesh), intent(in) :: hecMESH
26  type (fstr_solid),intent(inout) :: fstrSOLID
27  type(fstr_param), intent(in) :: fstrPARAM
28 
29  integer :: i,j,restrt_step(3),nif(2),istat(1),nload_prev(1)
30  real(kind=kreal) :: times(3)
31 
32  call hecmw_restart_open()
33 
34  call hecmw_restart_read_int(restrt_step)
35  if( fstrparam%restart_version >= 5 ) then
36  if( myrank == 0 ) write(*,*) 'Reading restart file as new format(>=ver5.0)'
37  call hecmw_restart_read_real(times)
38  call hecmw_restart_read_int(fstrsolid%NRstat_i)
39  call hecmw_restart_read_real(fstrsolid%NRstat_r)
40  call hecmw_restart_read_int(istat)
41  else
42  if( myrank == 0 ) write(*,*) 'Reading restart file as old format(<ver5.0)'
43  endif
44  call hecmw_restart_read_int(nload_prev) !load info at previous step
45  if( nload_prev(1)>0 ) then
46  allocate(fstrsolid%step_ctrl_restart%Load(nload_prev(1)))
47  call hecmw_restart_read_int(fstrsolid%step_ctrl_restart%Load)
48  endif
49 
50  call hecmw_restart_read_real(fstrsolid%unode)
51  call hecmw_restart_read_real(fstrsolid%QFORCE)
52 
53  do i= 1, hecmesh%n_elem
54  if (hecmw_is_etype_link( fstrsolid%elements(i)%etype )) cycle
55  if (hecmw_is_etype_patch( fstrsolid%elements(i)%etype )) cycle
56  do j= 1, size(fstrsolid%elements(i)%gausses)
57  call hecmw_restart_read_int(nif)
58  call hecmw_restart_read_real(fstrsolid%elements(i)%gausses(j)%strain)
59  call hecmw_restart_read_real(fstrsolid%elements(i)%gausses(j)%stress)
60  if( nif(1)>0 ) call hecmw_restart_read_int(fstrsolid%elements(i)%gausses(j)%istatus)
61  if( nif(2)>0 ) call hecmw_restart_read_real(fstrsolid%elements(i)%gausses(j)%fstatus)
62  enddo
63  enddo
64 
65  if( associated( fstrsolid%contacts ) ) then
66  call hecmw_restart_read_int(nif)
67  contactnode = nif(1)
68  do i= 1, size(fstrsolid%contacts)
69  do j= 1, size(fstrsolid%contacts(i)%slave)
70  call hecmw_restart_read_int(nif)
71  fstrsolid%contacts(i)%states(j)%surface = nif(1)
72  fstrsolid%contacts(i)%states(j)%state = nif(2)
73  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%lpos)
74  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%direction)
75  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%multiplier)
76  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%tangentForce_trial)
77  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%tangentForce_final)
78  enddo
79  enddo
80  endif
81 
82  call hecmw_restart_close()
83 
84  cstep = restrt_step(1)
85  substep = restrt_step(2) + 1
86  step_count = restrt_step(3)
87  if( fstrparam%restart_version >= 5 ) then
88  ctime = times(1)
89  dtime = times(2)
90  fstrsolid%AutoINC_stat = istat(1)
91  if( dabs(times(1)-times(3)) < 1.d-10 ) then
92  cstep = cstep + 1
93  substep = 1
94  endif
95  do i=1,size(fstrsolid%step_ctrl) !shift start time
96  fstrsolid%step_ctrl(i)%starttime = fstrsolid%step_ctrl(i)%starttime + times(3)
97  end do
98  else
99  ctime = fstrsolid%step_ctrl(cstep)%starttime
100  ctime = ctime + dble(substep-1)*fstrsolid%step_ctrl(cstep)%initdt
101  dtime = fstrsolid%step_ctrl(cstep)%initdt
102  if( dabs(ctime-fstrsolid%step_ctrl(cstep)%starttime-fstrsolid%step_ctrl(cstep)%elapsetime) < 1.d-10 ) then
103  cstep = cstep + 1
104  substep = 1
105  endif
106  endif
107 
108  end subroutine fstr_read_restart
109 
111  !----------------------------------------------------------------------*
112  subroutine fstr_write_restart(cstep,cstep_ext,substep,step_count,ctime,dtime,hecMESH, &
113  & fstrSOLID,fstrPARAM,is_StepFinished,contactNode)
114  !----------------------------------------------------------------------*
115  integer, intent(in) :: cstep
116  integer, intent(in) :: cstep_ext
117  integer, intent(in) :: substep
118  integer, intent(in) :: step_count
119  real(kind=kreal), intent(in) :: ctime
120  real(kind=kreal), intent(in) :: dtime
121  logical, intent(in) :: is_StepFinished
122  integer, intent(in) :: contactNode
123  type (hecmwST_local_mesh), intent(in) :: hecMESH
124  type (fstr_solid), intent(in) :: fstrSOLID
125  type(fstr_param), intent(in) :: fstrPARAM
126 
127  integer :: i,j,restrt_step(3),nif(2),istat(1),nload_prev(1)
128  real(kind=kreal) :: times(3)
129 
130  restrt_step(1) = cstep_ext
131  restrt_step(2) = substep
132  restrt_step(3) = step_count
133  times(1) = ctime
134  times(2) = dtime
135  if( is_stepfinished ) then
136  times(3) = ctime
137  else
138  times(3) = fstrsolid%step_ctrl(cstep)%starttime
139  end if
140  istat(1) = fstrsolid%AutoINC_stat
141  call hecmw_restart_add_int(restrt_step,size(restrt_step))
142  if( fstrparam%restart_version >= 5 ) then
143  call hecmw_restart_add_real(times,size(times))
144  call hecmw_restart_add_int(fstrsolid%NRstat_i,size(fstrsolid%NRstat_i))
145  call hecmw_restart_add_real(fstrsolid%NRstat_r,size(fstrsolid%NRstat_r))
146  call hecmw_restart_add_int(istat,1)
147  endif
148  nload_prev(1) = 0
149  if( is_stepfinished ) then
150  if( associated(fstrsolid%step_ctrl(cstep)%Load) ) nload_prev(1) = size(fstrsolid%step_ctrl(cstep)%Load)
151  call hecmw_restart_add_int(nload_prev,1)
152  if( nload_prev(1)>0 ) call hecmw_restart_add_int(fstrsolid%step_ctrl(cstep)%Load,nload_prev(1))
153  else
154  if( cstep>1 ) then
155  if( associated(fstrsolid%step_ctrl(cstep-1)%Load) ) nload_prev(1) = size(fstrsolid%step_ctrl(cstep-1)%Load)
156  call hecmw_restart_add_int(nload_prev,1)
157  if( nload_prev(1)>0 ) call hecmw_restart_add_int(fstrsolid%step_ctrl(cstep-1)%Load,nload_prev(1))
158  else
159  if( associated(fstrsolid%step_ctrl_restart%Load) ) nload_prev(1) = size(fstrsolid%step_ctrl_restart%Load)
160  call hecmw_restart_add_int(nload_prev,1)
161  if( nload_prev(1)>0 ) call hecmw_restart_add_int(fstrsolid%step_ctrl_restart%Load,nload_prev(1))
162  endif
163  end if
164 
165  call hecmw_restart_add_real(fstrsolid%unode,size(fstrsolid%unode))
166  call hecmw_restart_add_real(fstrsolid%QFORCE,size(fstrsolid%QFORCE))
167 
168  do i= 1, hecmesh%n_elem
169  if (hecmw_is_etype_link( fstrsolid%elements(i)%etype )) cycle
170  if (hecmw_is_etype_patch( fstrsolid%elements(i)%etype )) cycle
171  do j= 1, size(fstrsolid%elements(i)%gausses)
172  nif = 0
173  if( associated(fstrsolid%elements(i)%gausses(j)%istatus) ) nif(1)=size(fstrsolid%elements(i)%gausses(j)%istatus)
174  if( associated(fstrsolid%elements(i)%gausses(j)%fstatus) ) nif(2)=size(fstrsolid%elements(i)%gausses(j)%fstatus)
175  call hecmw_restart_add_int(nif,size(nif))
176  call hecmw_restart_add_real(fstrsolid%elements(i)%gausses(j)%strain,size(fstrsolid%elements(i)%gausses(j)%strain))
177  call hecmw_restart_add_real(fstrsolid%elements(i)%gausses(j)%stress,size(fstrsolid%elements(i)%gausses(j)%stress))
178  if( nif(1)>0 ) then
179  call hecmw_restart_add_int(fstrsolid%elements(i)%gausses(j)%istatus,size(fstrsolid%elements(i)%gausses(j)%istatus))
180  endif
181  if( nif(2)>0 ) then
182  call hecmw_restart_add_real(fstrsolid%elements(i)%gausses(j)%fstatus,size(fstrsolid%elements(i)%gausses(j)%fstatus))
183  endif
184  enddo
185  enddo
186 
187  if( associated( fstrsolid%contacts ) ) then
188  nif(1) = contactnode
189  call hecmw_restart_add_int(nif,size(nif))
190  do i= 1, size(fstrsolid%contacts)
191  do j= 1, size(fstrsolid%contacts(i)%slave)
192  nif(1) = fstrsolid%contacts(i)%states(j)%surface
193  nif(2) = fstrsolid%contacts(i)%states(j)%state
194  call hecmw_restart_add_int(nif,size(nif))
195  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%lpos,size(fstrsolid%contacts(i)%states(j)%lpos))
196  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%direction,size(fstrsolid%contacts(i)%states(j)%direction))
197  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%multiplier,size(fstrsolid%contacts(i)%states(j)%multiplier))
198  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%tangentForce_trial, &
199  size(fstrsolid%contacts(i)%states(j)%tangentForce_trial))
200  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%tangentForce_final, &
201  size(fstrsolid%contacts(i)%states(j)%tangentForce_final))
202  enddo
203  enddo
204  endif
205 
206  call hecmw_restart_write()
207 
208  end subroutine fstr_write_restart
209 
211  !----------------------------------------------------------------------*
212  subroutine fstr_read_restart_dyna_nl(cstep,hecMESH,fstrSOLID,fstrDYNAMIC,fstrPARAM,contactNode)
213  !----------------------------------------------------------------------*
214  integer, intent(out) :: cstep
215  integer, intent(out), optional :: contactNode
216  type (hecmwST_local_mesh), intent(in) :: hecMESH
217  type (fstr_solid),intent(inout) :: fstrSOLID
218  type ( fstr_dynamic), intent(inout) :: fstrDYNAMIC
219  type(fstr_param), intent(in) :: fstrPARAM
220 
221  integer :: i,j,restrt_step(1),nif(2)
222  real(kind=kreal) :: data(2)
223 
224  call hecmw_restart_open()
225 
226  call hecmw_restart_read_int(restrt_step)
227  cstep = restrt_step(1)
228  call hecmw_restart_read_real(fstrsolid%unode)
229  call hecmw_restart_read_real(fstrsolid%QFORCE)
230 
231  do i= 1, hecmesh%n_elem
232  if (hecmw_is_etype_link( fstrsolid%elements(i)%etype )) cycle
233  if (hecmw_is_etype_patch( fstrsolid%elements(i)%etype )) cycle
234  do j= 1, size(fstrsolid%elements(i)%gausses)
235  call hecmw_restart_read_int(nif)
236  call hecmw_restart_read_real(fstrsolid%elements(i)%gausses(j)%strain)
237  call hecmw_restart_read_real(fstrsolid%elements(i)%gausses(j)%stress)
238  if( nif(1)>0 ) call hecmw_restart_read_int(fstrsolid%elements(i)%gausses(j)%istatus)
239  if( nif(2)>0 ) call hecmw_restart_read_real(fstrsolid%elements(i)%gausses(j)%fstatus)
240  enddo
241  enddo
242 
243  if(present(contactnode)) then
244  call hecmw_restart_read_int(nif)
245  contactnode = nif(1)
246  do i= 1, size(fstrsolid%contacts)
247  do j= 1, size(fstrsolid%contacts(i)%slave)
248  call hecmw_restart_read_int(nif)
249  fstrsolid%contacts(i)%states(j)%surface = nif(1)
250  fstrsolid%contacts(i)%states(j)%state = nif(2)
251  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%lpos)
252  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%direction)
253  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%multiplier)
254  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%tangentForce_trial)
255  call hecmw_restart_read_real(fstrsolid%contacts(i)%states(j)%tangentForce_final)
256  enddo
257  enddo
258  endif
259 
260  call hecmw_restart_read_int(restrt_step)
261  fstrdynamic%idx_eqa = restrt_step(1)
262  call hecmw_restart_read_real(data)
263  fstrdynamic%t_curr = data(1)
264  fstrdynamic%strainEnergy = data(2)
265  if( fstrdynamic%idx_eqa == 1 ) then
266  call hecmw_restart_read_real(fstrdynamic%DISP(:,1))
267  call hecmw_restart_read_real(fstrdynamic%VEL(:,1))
268  call hecmw_restart_read_real(fstrdynamic%ACC(:,1))
269  else
270  call hecmw_restart_read_real(fstrdynamic%DISP(:,1))
271  call hecmw_restart_read_real(fstrdynamic%DISP(:,3))
272  endif
273  do i= 1, hecmesh%n_elem
274  call hecmw_restart_read_real(fstrsolid%elements(i)%equiForces)
275  enddo
276 
277  call hecmw_restart_close()
278 
279  end subroutine fstr_read_restart_dyna_nl
280 
282  !----------------------------------------------------------------------*
283  subroutine fstr_write_restart_dyna_nl(cstep,hecMESH,fstrSOLID,fstrDYNAMIC,fstrPARAM,contactNode)
284  !----------------------------------------------------------------------*
285  integer, intent(in) :: cstep
286  integer, intent(in), optional :: contactNode
287  type (hecmwST_local_mesh), intent(in) :: hecMESH
288  type (fstr_solid), intent(in) :: fstrSOLID
289  type ( fstr_dynamic), intent(in) :: fstrDYNAMIC
290  type(fstr_param), intent(in) :: fstrPARAM
291 
292  integer :: i,j,restrt_step(1),nif(2)
293  real(kind=kreal) :: data(2)
294 
295  restrt_step(1) = cstep
296  call hecmw_restart_add_int(restrt_step,size(restrt_step))
297  call hecmw_restart_add_real(fstrsolid%unode,size(fstrsolid%unode))
298  call hecmw_restart_add_real(fstrsolid%QFORCE,size(fstrsolid%QFORCE))
299 
300  do i= 1, hecmesh%n_elem
301  if (hecmw_is_etype_link( fstrsolid%elements(i)%etype )) cycle
302  if (hecmw_is_etype_patch( fstrsolid%elements(i)%etype )) cycle
303  do j= 1, size(fstrsolid%elements(i)%gausses)
304  nif = 0
305  if( associated(fstrsolid%elements(i)%gausses(j)%istatus) ) nif(1)=size(fstrsolid%elements(i)%gausses(j)%istatus)
306  if( associated(fstrsolid%elements(i)%gausses(j)%fstatus) ) nif(2)=size(fstrsolid%elements(i)%gausses(j)%fstatus)
307  call hecmw_restart_add_int(nif,size(nif))
308  call hecmw_restart_add_real(fstrsolid%elements(i)%gausses(j)%strain,size(fstrsolid%elements(i)%gausses(j)%strain))
309  call hecmw_restart_add_real(fstrsolid%elements(i)%gausses(j)%stress,size(fstrsolid%elements(i)%gausses(j)%stress))
310  if( nif(1)>0 ) then
311  call hecmw_restart_add_int(fstrsolid%elements(i)%gausses(j)%istatus,size(fstrsolid%elements(i)%gausses(j)%istatus))
312  endif
313  if( nif(2)>0 ) then
314  call hecmw_restart_add_real(fstrsolid%elements(i)%gausses(j)%fstatus,size(fstrsolid%elements(i)%gausses(j)%fstatus))
315  endif
316  enddo
317  enddo
318 
319  if(present(contactnode)) then
320  nif(1) = contactnode
321  call hecmw_restart_add_int(nif,size(nif))
322  do i= 1, size(fstrsolid%contacts)
323  do j= 1, size(fstrsolid%contacts(i)%slave)
324  nif(1) = fstrsolid%contacts(i)%states(j)%surface
325  nif(2) = fstrsolid%contacts(i)%states(j)%state
326  call hecmw_restart_add_int(nif,size(nif))
327  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%lpos,size(fstrsolid%contacts(i)%states(j)%lpos))
328  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%direction,size(fstrsolid%contacts(i)%states(j)%direction))
329  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%multiplier,size(fstrsolid%contacts(i)%states(j)%multiplier))
330  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%tangentForce_trial, &
331  size(fstrsolid%contacts(i)%states(j)%tangentForce_trial))
332  call hecmw_restart_add_real(fstrsolid%contacts(i)%states(j)%tangentForce_final, &
333  size(fstrsolid%contacts(i)%states(j)%tangentForce_final))
334  enddo
335  enddo
336  endif
337 
338  restrt_step(1) = fstrdynamic%idx_eqa
339  call hecmw_restart_add_int(restrt_step,size(restrt_step))
340  data(1) = fstrdynamic%t_curr
341  data(2) = fstrdynamic%strainEnergy
342  call hecmw_restart_add_real(data,size(data))
343  if( fstrdynamic%idx_eqa == 1 ) then
344  call hecmw_restart_add_real(fstrdynamic%DISP(:,1),size(fstrdynamic%DISP(:,1)))
345  call hecmw_restart_add_real(fstrdynamic%VEL(:,1),size(fstrdynamic%VEL(:,1)))
346  call hecmw_restart_add_real(fstrdynamic%ACC(:,1),size(fstrdynamic%ACC(:,1)))
347  else
348  call hecmw_restart_add_real(fstrdynamic%DISP(:,1),size(fstrdynamic%DISP(:,1)))
349  call hecmw_restart_add_real(fstrdynamic%DISP(:,3),size(fstrdynamic%DISP(:,3)))
350  endif
351  do i= 1, hecmesh%n_elem
352  call hecmw_restart_add_real(fstrsolid%elements(i)%equiForces,size(fstrsolid%elements(i)%equiForces))
353  enddo
354 
355  call hecmw_restart_write()
356 
357  end subroutine fstr_write_restart_dyna_nl
358 
359 end module m_fstr_restart
This module provides functions to read in and write out restart fiels.
Definition: fstr_Restart.f90:8
subroutine fstr_read_restart(cstep, substep, step_count, ctime, dtime, hecMESH, fstrSOLID, fstrPARAM, contactNode)
Read in restart file.
subroutine fstr_write_restart_dyna_nl(cstep, hecMESH, fstrSOLID, fstrDYNAMIC, fstrPARAM, contactNode)
write out restart file for nonlinear dynamic analysis
subroutine fstr_read_restart_dyna_nl(cstep, hecMESH, fstrSOLID, fstrDYNAMIC, fstrPARAM, contactNode)
Read in restart file for nonlinear dynamic analysis.
subroutine fstr_write_restart(cstep, cstep_ext, substep, step_count, ctime, dtime, hecMESH, fstrSOLID, fstrPARAM, is_StepFinished, contactNode)
write out restart file
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
This module provides aux functions.
Definition: utilities.f90:6
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138