15 logical,
private :: is_cutback_active = .false.
25 type(hecmwst_local_mesh) :: hecMESH
29 integer(kind=kint) :: istep, i, j
30 integer(kind=kint) :: ng
31 integer(kind=kint) :: ncont, nstate
33 do istep=1,fstrsolid%nstep_tot
34 if( fstrsolid%step_ctrl(istep)%inc_type == stepautoinc ) is_cutback_active = .true.
37 if( .not. is_cutback_active )
return
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)))
48 allocate(fstrsolid%elements_bkup(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
53 allocate( fstrsolid%elements_bkup(i)%gausses(ng) )
55 fstrsolid%elements_bkup(i)%gausses(j)%pMaterial => fstrsolid%elements(i)%gausses(j)%pMaterial
58 if(
associated( fstrsolid%elements(i)%aux ) )
then
59 allocate( fstrsolid%elements_bkup(i)%aux(3,3) )
64 ncont =
size(fstrsolid%contacts)
65 if(
associated(fstrsolid%contacts) )
then
66 allocate(fstrsolid%contacts_bkup(ncont))
68 nstate =
size(fstrsolid%contacts(i)%states)
69 allocate(fstrsolid%contacts_bkup(i)%states(nstate))
79 integer(kind=kint) :: i, j
80 integer(kind=kint) :: ng
81 integer(kind=kint) :: ncont
83 if( .not. is_cutback_active )
return
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)
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
101 deallocate( fstrsolid%elements_bkup(i)%gausses )
102 if(
associated( fstrsolid%elements_bkup(i)%aux ) )
then
103 deallocate( fstrsolid%elements_bkup(i)%aux )
106 deallocate(fstrsolid%elements_bkup)
109 ncont =
size(fstrsolid%contacts)
110 if(
associated(fstrsolid%contacts) )
then
112 deallocate(fstrsolid%contacts_bkup(i)%states)
114 deallocate(fstrsolid%contacts_bkup)
124 integer(kind=kint) :: i, j
125 integer(kind=kint) :: ng
126 integer(kind=kint) :: ncont, nstate
128 if( .not. is_cutback_active )
return
131 do i=1,
size(fstrsolid%unode)
132 fstrsolid%unode_bkup(i) = fstrsolid%unode(i)
134 do i=1,
size(fstrsolid%QFORCE)
135 fstrsolid%QFORCE_bkup(i) = fstrsolid%QFORCE(i)
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)
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
149 call fstr_copy_gauss( fstrsolid%elements(i)%gausses(j), fstrsolid%elements_bkup(i)%gausses(j) )
151 if(
associated( fstrsolid%elements(i)%aux ) )
then
152 fstrsolid%elements_bkup(i)%aux(:,:) = fstrsolid%elements(i)%aux(:,:)
157 ncont =
size(fstrsolid%contacts)
158 if(
associated(fstrsolid%contacts) )
then
160 nstate =
size(fstrsolid%contacts(i)%states)
162 call contact_state_copy(fstrsolid%contacts(i)%states(j), fstrsolid%contacts_bkup(i)%states(j))
165 infoctchange_bak = infoctchange
175 integer(kind=kint) :: i, j
176 integer(kind=kint) :: ng
177 integer(kind=kint) :: ncont, nstate
179 if( .not. is_cutback_active )
return
182 do i=1,
size(fstrsolid%unode)
183 fstrsolid%unode(i) = fstrsolid%unode_bkup(i)
185 do i=1,
size(fstrsolid%QFORCE)
186 fstrsolid%QFORCE(i) = fstrsolid%QFORCE_bkup(i)
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)
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
200 call fstr_copy_gauss( fstrsolid%elements_bkup(i)%gausses(j), fstrsolid%elements(i)%gausses(j) )
202 if(
associated( fstrsolid%elements(i)%aux ) )
then
203 fstrsolid%elements(i)%aux(:,:) = fstrsolid%elements_bkup(i)%aux(:,:)
208 ncont =
size(fstrsolid%contacts)
209 if(
associated(fstrsolid%contacts) )
then
211 nstate =
size(fstrsolid%contacts(i)%states)
213 call contact_state_copy(fstrsolid%contacts_bkup(i)%states(j), fstrsolid%contacts(i)%states(j))
216 infoctchange = infoctchange_bak
This module encapsulate the basic functions of all elements provide by this software.
integer function numofquadpoints(fetype)
Obtains the number of quadrature points of the element.
This module provides functions to deal with cutback.
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.
This modules defines a structure to record history dependent parameter in static analysis.
subroutine fstr_init_gauss(gauss)
Initializer.
subroutine fstr_finalize_gauss(gauss)
Finializer.
subroutine fstr_copy_gauss(gauss1, gauss2)
Copy.
FSTR INNER CONTROL PARAMETERS (fstrPARAM)