FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
dynamic_mat_ass_bc_ac.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 contains
9 
10  !C***
12  !C***
13 
14  subroutine dynamic_mat_ass_bc_ac(hecMESH, hecMAT, fstrSOLID ,fstrDYNAMIC, fstrPARAM, fstrMAT, iter, conMAT)
15  use m_fstr
16  use m_table_dyn
19  use mcontact
20 
21  implicit none
22  type(hecmwst_matrix) :: hecMAT
23  type(hecmwst_local_mesh) :: hecMESH
24  type(fstr_solid) :: fstrSOLID
25  type(fstr_dynamic) :: fstrDYNAMIC
26  type(fstr_param) :: fstrPARAM
27  type(fstrst_matrix_contact_lagrange) :: fstrMAT
28  integer, optional :: iter
29  type(hecmwst_matrix), optional :: conMAT
30 
31  integer(kind=kint) :: ig0, ig, ityp, NDOF, iS0, iE0, ik, in, idofS, idofE, idof
32  integer(kind=kint) :: dyn_step, flag_u
33  real(kind=kreal) :: b2, b3, b4, c1
34  real(kind=kreal) :: rhs, rhs0, f_t
35 
36  if( fstrsolid%ACCELERATION_type == kbcinitial )return
37 
38  dyn_step = fstrdynamic%i_step
39  flag_u = 3
40 
41  b2 = fstrdynamic%t_delta
42  b3 = fstrdynamic%t_delta**2*(0.5d0-fstrdynamic%beta)
43  b4 = fstrdynamic%t_delta**2*fstrdynamic%beta
44  c1 = fstrdynamic%t_delta**2
45 
46  ndof = hecmat%NDOF
47 
48  !C=============================C
49  !C-- implicit dynamic analysis
50  !C=============================C
51  if( fstrdynamic%idx_eqa == 1 ) then
52 
53  do ig0 = 1, fstrsolid%ACCELERATION_ngrp_tot
54  ig = fstrsolid%ACCELERATION_ngrp_ID(ig0)
55  rhs = fstrsolid%ACCELERATION_ngrp_val(ig0)
56 
57  call table_dyn(hecmesh, fstrsolid, fstrdynamic, ig0, f_t, flag_u)
58  rhs = rhs * f_t
59  rhs0 = rhs
60 
61  ityp = fstrsolid%ACCELERATION_ngrp_type(ig0)
62 
63  idofs = ityp/10
64  idofe = ityp - idofs*10
65 
66  is0 = hecmesh%node_group%grp_index(ig-1) + 1
67  ie0 = hecmesh%node_group%grp_index(ig )
68 
69  do ik = is0, ie0
70  in = hecmesh%node_group%grp_item(ik)
71  do idof = idofs, idofe
72 
73  if( present(iter) ) then ! increment
74  if( iter>1 ) then
75  rhs = 0.d0
76  else
77  rhs = &
78  + b2*fstrdynamic%VEL (ndof*in-(ndof-idof),1) &
79  + b3*fstrdynamic%ACC (ndof*in-(ndof-idof),1) &
80  + b4*rhs0
81  endif
82  else
83  rhs = fstrdynamic%DISP(ndof*in-(ndof-idof),1) &
84  + b2*fstrdynamic%VEL (ndof*in-(ndof-idof),1) &
85  + b3*fstrdynamic%ACC (ndof*in-(ndof-idof),1) &
86  + b4*rhs0
87  endif
88  if(present(conmat)) then
89  call hecmw_mat_ass_bc(hecmat, in, idof, rhs, conmat)
90  else
91  call hecmw_mat_ass_bc(hecmat, in, idof, rhs)
92  endif
93  if( fstr_is_contact_active() .and. fstrparam%contact_algo == kcaslagrange &
94  .and. fstrparam%nlgeom .and. fstrdynamic%idx_resp == 1 ) then
95  if(present(conmat)) then
96  call fstr_mat_ass_bc_contact(conmat,fstrmat,in,idof,rhs)
97  else
98  call fstr_mat_ass_bc_contact(hecmat,fstrmat,in,idof,rhs)
99  endif
100  endif
101 
102  !for output reaction force
103  fstrsolid%REACTION(ndof*(in-1)+idof) = fstrsolid%QFORCE(ndof*(in-1)+idof)
104  enddo
105  enddo
106 
107  enddo
108  !C
109  !C-- end of implicit dynamic analysis
110  !C
111 
112  !C=============================C
113  !C-- explicit dynamic analysis
114  !C=============================C
115  else if( fstrdynamic%idx_eqa == 11 ) then
116  !C
117  do ig0 = 1, fstrsolid%ACCELERATION_ngrp_tot
118  ig = fstrsolid%ACCELERATION_ngrp_ID(ig0)
119  rhs = fstrsolid%ACCELERATION_ngrp_val(ig0)
120 
121  call table_dyn(hecmesh, fstrsolid, fstrdynamic, ig0, f_t, flag_u)
122  rhs = rhs * f_t
123  rhs0 = rhs
124 
125  ityp = fstrsolid%ACCELERATION_ngrp_type(ig0)
126 
127  is0 = hecmesh%node_group%grp_index(ig-1) + 1
128  ie0 = hecmesh%node_group%grp_index(ig )
129  idofs = ityp/10
130  idofe = ityp - idofs*10
131 
132  do ik = is0, ie0
133  in = hecmesh%node_group%grp_item(ik)
134 
135  do idof = idofs, idofe
136  rhs = 2.0*fstrdynamic%DISP(ndof*in-(ndof-idof),1) &
137  - fstrdynamic%DISP(ndof*in-(ndof-idof),3) &
138  + c1*rhs0
139  hecmat%B (ndof*in-(ndof-idof)) = rhs
140  fstrdynamic%VEC1(ndof*in-(ndof-idof)) = 1.0d0
141 
142  !for output reaction force
143  fstrsolid%REACTION(ndof*(in-1)+idof) = fstrsolid%QFORCE(ndof*(in-1)+idof)
144  end do
145  enddo
146  enddo
147  !C
148  !C-- end of explicit dynamic analysis
149  !C
150  end if
151  !
152  return
153  end subroutine dynamic_mat_ass_bc_ac
154 
155 
156  !C***
158  !C***
159  subroutine dynamic_bc_init_ac(hecMESH, hecMAT, fstrSOLID ,fstrDYNAMIC)
160  use m_fstr
161  use m_table_dyn
162 
163  implicit none
164  type(hecmwst_matrix) :: hecmat
165  type(hecmwst_local_mesh) :: hecMESH
166  type(fstr_solid) :: fstrSOLID
167  type(fstr_dynamic) :: fstrDYNAMIC
168 
169  integer(kind=kint) :: NDOF, ig0, ig, ityp, iS0, iE0, ik, in, idofS, idofE, idof
170  !!!
171  integer(kind=kint) :: flag_u
172  real(kind=kreal) :: rhs, f_t
173 
174  if( fstrsolid%ACCELERATION_type == kbctransit )return
175 
176  !!!
177  flag_u = 3
178  !!!
179 
180  ndof = hecmat%NDOF
181  do ig0 = 1, fstrsolid%ACCELERATION_ngrp_tot
182  ig = fstrsolid%ACCELERATION_ngrp_ID(ig0)
183  rhs = fstrsolid%ACCELERATION_ngrp_val(ig0)
184 
185  !!!!!! time history
186  call table_dyn(hecmesh, fstrsolid, fstrdynamic, ig0, f_t, flag_u)
187  rhs = rhs * f_t
188  !!!!!!
189  ityp = fstrsolid%ACCELERATION_ngrp_type(ig0)
190 
191  is0 = hecmesh%node_group%grp_index(ig-1) + 1
192  ie0 = hecmesh%node_group%grp_index(ig )
193  idofs = ityp/10
194  idofe = ityp - idofs*10
195 
196  do ik = is0, ie0
197  in = hecmesh%node_group%grp_item(ik)
198  !C
199  do idof = idofs, idofe
200  fstrdynamic%ACC (ndof*in-(ndof-idof),1) = rhs
201  end do
202  enddo
203  !*End ig0 (main) loop
204  enddo
205 
206  return
207  end subroutine dynamic_bc_init_ac
208 
209  subroutine dynamic_explicit_ass_ac(hecMESH, hecMAT, fstrSOLID ,fstrDYNAMIC, iter)
210  use m_fstr
211  use m_table_dyn
214  use mcontact
215  type(hecmwst_matrix) :: hecmat
216  type(hecmwst_local_mesh) :: hecMESH
217  type(fstr_solid) :: fstrSOLID
218  type(fstr_dynamic) :: fstrDYNAMIC
219  integer, optional :: iter
220 
221  integer(kind=kint) :: ig0, ig, ityp, NDOF, iS0, iE0, ik, in, idofS, idofE, idof
222  integer(kind=kint) :: dyn_step, flag_u
223  real(kind=kreal) :: b2, b3, b4, c1
224  real(kind=kreal) :: rhs, rhs0, f_t
225 
226  if( fstrsolid%ACCELERATION_type == kbcinitial )return
227 
228  dyn_step = fstrdynamic%i_step
229  flag_u = 3
230 
231  b2 = fstrdynamic%t_delta
232  b3 = fstrdynamic%t_delta**2*(0.5d0-fstrdynamic%beta)
233  b4 = fstrdynamic%t_delta**2*fstrdynamic%beta
234  c1 = fstrdynamic%t_delta**2
235 
236  ndof = hecmat%NDOF
237 
238  do ig0 = 1, fstrsolid%ACCELERATION_ngrp_tot
239  ig = fstrsolid%ACCELERATION_ngrp_ID(ig0)
240  rhs = fstrsolid%ACCELERATION_ngrp_val(ig0)
241 
242  call table_dyn(hecmesh, fstrsolid, fstrdynamic, ig0, f_t, flag_u)
243  rhs = rhs * f_t
244  rhs0 = rhs
245 
246  ityp = fstrsolid%ACCELERATION_ngrp_type(ig0)
247 
248  is0 = hecmesh%node_group%grp_index(ig-1) + 1
249  ie0 = hecmesh%node_group%grp_index(ig )
250  idofs = ityp/10
251  idofe = ityp - idofs*10
252 
253  do ik = is0, ie0
254  in = hecmesh%node_group%grp_item(ik)
255 
256  do idof = idofs, idofe
257  rhs = 2.0*fstrdynamic%DISP(ndof*in-(ndof-idof),1) &
258  - fstrdynamic%DISP(ndof*in-(ndof-idof),3) &
259  + c1*rhs0
260  hecmat%B(ndof*in-(ndof-idof)) = rhs* fstrdynamic%VEC1(ndof*in-(ndof-idof))
261  ! fstrDYNAMIC%VEC1(NDOF*in-(NDOF-idof)) = 1.0d0
262 
263  !for output reaction force
264  fstrsolid%REACTION(ndof*(in-1)+idof) = fstrsolid%QFORCE(ndof*(in-1)+idof)
265  end do
266  enddo
267  enddo
268  end subroutine dynamic_explicit_ass_ac
269 
270 end module m_dynamic_mat_ass_bc_ac
This module provides functions of reconstructing.
This module provides functions: 1) obtain contact stiffness matrix of each contact pair and assemble ...
subroutine, public fstr_mat_ass_bc_contact(hecMAT, fstrMAT, inode, idof, RHS)
Modify Lagrange multiplier-related part of stiffness matrix and right-hand side vector for dealing wi...
This module contains functions to set acceleration boundary condition in dynamic analysis.
subroutine dynamic_mat_ass_bc_ac(hecMESH, hecMAT, fstrSOLID, fstrDYNAMIC, fstrPARAM, fstrMAT, iter, conMAT)
This subrouitne set acceleration boundary condition in dynamic analysis.
subroutine dynamic_explicit_ass_ac(hecMESH, hecMAT, fstrSOLID, fstrDYNAMIC, iter)
subroutine dynamic_bc_init_ac(hecMESH, hecMAT, fstrSOLID, fstrDYNAMIC)
This function sets initial condition of acceleration.
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint), parameter kcaslagrange
contact analysis algorithm
Definition: m_fstr.f90:53
integer(kind=kint), parameter kbcinitial
Definition: m_fstr.f90:60
integer(kind=kint), parameter kbctransit
Definition: m_fstr.f90:61
Table of lading step in dynamic analysis.
Definition: table_dyn.f90:6
subroutine table_dyn(hecMESH, fstrSOLID, fstrDYNAMIC, ig0, f_t, flag_u)
Definition: table_dyn.f90:19
This module provides functions to calcualte contact stiff matrix.
Definition: fstr_contact.f90:6
logical function fstr_is_contact_active()
Structure for Lagrange multiplier-related part of stiffness matrix (Lagrange multiplier-related matri...
Data for DYNAMIC ANSLYSIS (fstrDYNAMIC)
Definition: m_fstr.f90:473
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138