FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
fstr_Cutback.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  use m_fstr
9  use elementinfo
10  use mmechgauss
11  use mcontactdef
12 
13  implicit none
14 
15  logical, private :: is_cutback_active = .false.
16 
17 contains
18 
19  logical function fstr_cutback_active()
20  fstr_cutback_active = is_cutback_active
21  end function
22 
24  subroutine fstr_cutback_init( hecMESH, fstrSOLID, fstrPARAM )
25  type(hecmwst_local_mesh) :: hecMESH
26  type(fstr_param) :: fstrPARAM
27  type(fstr_solid) :: fstrSOLID
28 
29  integer(kind=kint) :: istep, i, j
30  integer(kind=kint) :: ng
31  integer(kind=kint) :: ncont, nstate
32 
33  do istep=1,fstrsolid%nstep_tot
34  if( fstrsolid%step_ctrl(istep)%inc_type == stepautoinc ) is_cutback_active = .true.
35  end do
36 
37  if( .not. is_cutback_active ) return
38 
39  !allocate variables to store analysis state
40  !(1)nodal values
41  allocate(fstrsolid%unode_bkup(size(fstrsolid%unode)))
42  allocate(fstrsolid%QFORCE_bkup(size(fstrsolid%QFORCE)))
43  if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres > 0 ) then
44  allocate(fstrsolid%last_temp_bkup(size(fstrsolid%last_temp)))
45  endif
46 
47  !(2)elemental values
48  allocate(fstrsolid%elements_bkup(hecmesh%n_elem))
49  do i=1,hecmesh%n_elem
50  if (hecmw_is_etype_link( fstrsolid%elements(i)%etype )) cycle
51  if (hecmw_is_etype_patch( fstrsolid%elements(i)%etype )) cycle
52  ng = numofquadpoints( fstrsolid%elements(i)%etype )
53  allocate( fstrsolid%elements_bkup(i)%gausses(ng) )
54  do j=1,ng
55  fstrsolid%elements_bkup(i)%gausses(j)%pMaterial => fstrsolid%elements(i)%gausses(j)%pMaterial
56  call fstr_init_gauss( fstrsolid%elements_bkup(i)%gausses(j) )
57  end do
58  if( associated( fstrsolid%elements(i)%aux ) ) then
59  allocate( fstrsolid%elements_bkup(i)%aux(3,3) )
60  endif
61  end do
62 
63  !(3)contact values
64  ncont = size(fstrsolid%contacts)
65  if( associated(fstrsolid%contacts) ) then
66  allocate(fstrsolid%contacts_bkup(ncont))
67  do i=1,ncont
68  nstate = size(fstrsolid%contacts(i)%states)
69  allocate(fstrsolid%contacts_bkup(i)%states(nstate))
70  end do
71  end if
72 
73  end subroutine
74 
76  subroutine fstr_cutback_finalize( fstrSOLID )
77  type(fstr_solid) :: fstrSOLID
78 
79  integer(kind=kint) :: i, j
80  integer(kind=kint) :: ng
81  integer(kind=kint) :: ncont
82 
83  if( .not. is_cutback_active ) return
84 
85  !deallocate variables to store analysis state
86  !(1)nodal values
87  if( associated(fstrsolid%unode_bkup) ) deallocate(fstrsolid%unode_bkup)
88  if( associated(fstrsolid%QFORCE_bkup) ) deallocate(fstrsolid%QFORCE_bkup)
89  if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres > 0 ) then
90  if( associated(fstrsolid%last_temp_bkup) ) deallocate(fstrsolid%last_temp_bkup)
91  endif
92 
93  !(2)elemental values
94  do i=1,size(fstrsolid%elements)
95  if (hecmw_is_etype_link( fstrsolid%elements(i)%etype )) cycle
96  if (hecmw_is_etype_patch( fstrsolid%elements(i)%etype )) cycle
97  ng = numofquadpoints( fstrsolid%elements(i)%etype )
98  do j=1,ng
99  call fstr_finalize_gauss( fstrsolid%elements_bkup(i)%gausses(j) )
100  end do
101  deallocate( fstrsolid%elements_bkup(i)%gausses )
102  if( associated( fstrsolid%elements_bkup(i)%aux ) ) then
103  deallocate( fstrsolid%elements_bkup(i)%aux )
104  endif
105  end do
106  deallocate(fstrsolid%elements_bkup)
107 
108  !(3)contact values
109  ncont = size(fstrsolid%contacts)
110  if( associated(fstrsolid%contacts) ) then
111  do i=1,ncont
112  deallocate(fstrsolid%contacts_bkup(i)%states)
113  end do
114  deallocate(fstrsolid%contacts_bkup)
115  end if
116  end subroutine
117 
119  subroutine fstr_cutback_save( fstrSOLID, infoCTChange, infoCTChange_bak )
120  type(fstr_solid), intent(inout) :: fstrSOLID
121  type(fstr_info_contactchange), intent(inout) :: infoCTChange
122  type(fstr_info_contactchange), intent(inout) :: infoCTChange_bak
123 
124  integer(kind=kint) :: i, j
125  integer(kind=kint) :: ng
126  integer(kind=kint) :: ncont, nstate
127 
128  if( .not. is_cutback_active ) return
129 
130  !(1)nodal values
131  do i=1,size(fstrsolid%unode)
132  fstrsolid%unode_bkup(i) = fstrsolid%unode(i)
133  end do
134  do i=1,size(fstrsolid%QFORCE)
135  fstrsolid%QFORCE_bkup(i) = fstrsolid%QFORCE(i)
136  end do
137  if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres > 0 ) then
138  do i=1,size(fstrsolid%last_temp)
139  fstrsolid%last_temp_bkup(i) = fstrsolid%last_temp(i)
140  end do
141  endif
142 
143  !(2)elemental values
144  do i=1,size(fstrsolid%elements)
145  if (hecmw_is_etype_link( fstrsolid%elements(i)%etype )) cycle
146  if (hecmw_is_etype_patch( fstrsolid%elements(i)%etype )) cycle
147  ng = numofquadpoints( fstrsolid%elements(i)%etype )
148  do j=1,ng
149  call fstr_copy_gauss( fstrsolid%elements(i)%gausses(j), fstrsolid%elements_bkup(i)%gausses(j) )
150  end do
151  if( associated( fstrsolid%elements(i)%aux ) ) then
152  fstrsolid%elements_bkup(i)%aux(:,:) = fstrsolid%elements(i)%aux(:,:)
153  endif
154  end do
155 
156  !(3)contact values
157  ncont = size(fstrsolid%contacts)
158  if( associated(fstrsolid%contacts) ) then
159  do i=1,ncont
160  nstate = size(fstrsolid%contacts(i)%states)
161  do j=1,nstate
162  call contact_state_copy(fstrsolid%contacts(i)%states(j), fstrsolid%contacts_bkup(i)%states(j))
163  enddo
164  end do
165  infoctchange_bak = infoctchange
166  end if
167  end subroutine
168 
170  subroutine fstr_cutback_load( fstrSOLID, infoCTChange, infoCTChange_bak )
171  type(fstr_solid), intent(inout) :: fstrSOLID
172  type(fstr_info_contactchange), intent(inout) :: infoCTChange
173  type(fstr_info_contactchange), intent(inout) :: infoCTChange_bak
174 
175  integer(kind=kint) :: i, j
176  integer(kind=kint) :: ng
177  integer(kind=kint) :: ncont, nstate
178 
179  if( .not. is_cutback_active ) return
180 
181  !(1)nodal values
182  do i=1,size(fstrsolid%unode)
183  fstrsolid%unode(i) = fstrsolid%unode_bkup(i)
184  end do
185  do i=1,size(fstrsolid%QFORCE)
186  fstrsolid%QFORCE(i) = fstrsolid%QFORCE_bkup(i)
187  end do
188  if( fstrsolid%TEMP_ngrp_tot > 0 .or. fstrsolid%TEMP_irres > 0 ) then
189  do i=1,size(fstrsolid%last_temp)
190  fstrsolid%last_temp(i) = fstrsolid%last_temp_bkup(i)
191  end do
192  endif
193 
194  !(2)elemental values
195  do i=1,size(fstrsolid%elements)
196  if (hecmw_is_etype_link( fstrsolid%elements(i)%etype )) cycle
197  if (hecmw_is_etype_patch( fstrsolid%elements(i)%etype )) cycle
198  ng = numofquadpoints( fstrsolid%elements(i)%etype )
199  do j=1,ng
200  call fstr_copy_gauss( fstrsolid%elements_bkup(i)%gausses(j), fstrsolid%elements(i)%gausses(j) )
201  end do
202  if( associated( fstrsolid%elements(i)%aux ) ) then
203  fstrsolid%elements(i)%aux(:,:) = fstrsolid%elements_bkup(i)%aux(:,:)
204  endif
205  end do
206 
207  !(3)contact values
208  ncont = size(fstrsolid%contacts)
209  if( associated(fstrsolid%contacts) ) then
210  do i=1,ncont
211  nstate = size(fstrsolid%contacts(i)%states)
212  do j=1,nstate
213  call contact_state_copy(fstrsolid%contacts_bkup(i)%states(j), fstrsolid%contacts(i)%states(j))
214  enddo
215  end do
216  infoctchange = infoctchange_bak
217  end if
218 
219  end subroutine
220 
221 end module m_fstr_cutback
This module encapsulate the basic functions of all elements provide by this software.
Definition: element.f90:43
integer function numofquadpoints(fetype)
Obtains the number of quadrature points of the element.
Definition: element.f90:445
This module provides functions to deal with cutback.
Definition: fstr_Cutback.f90:7
subroutine fstr_cutback_save(fstrSOLID, infoCTChange, infoCTChange_bak)
Save analysis status.
subroutine fstr_cutback_load(fstrSOLID, infoCTChange, infoCTChange_bak)
Load analysis status.
subroutine fstr_cutback_init(hecMESH, fstrSOLID, fstrPARAM)
Initializer of cutback variables.
subroutine fstr_cutback_finalize(fstrSOLID)
Finalizer of cutback variables.
logical function fstr_cutback_active()
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
This module manage the data structure for contact calculation.
This modules defines a structure to record history dependent parameter in static analysis.
Definition: mechgauss.f90:6
subroutine fstr_init_gauss(gauss)
Initializer.
Definition: mechgauss.f90:41
subroutine fstr_finalize_gauss(gauss)
Finializer.
Definition: mechgauss.f90:79
subroutine fstr_copy_gauss(gauss1, gauss2)
Copy.
Definition: mechgauss.f90:86
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138