FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
dynamic_mat_ass_bc_vl.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 contains
8 
9  !C***
11  !C***
12 
13  subroutine dynamic_mat_ass_bc_vl(hecMESH, hecMAT, fstrSOLID, fstrDYNAMIC, fstrPARAM, fstrMAT, iter, conMAT)
14  use m_fstr
15  use m_table_dyn
18  use mcontact
19 
20  implicit none
21  type(hecmwst_matrix) :: hecMAT
22  type(hecmwst_local_mesh) :: hecMESH
23  type(fstr_solid) :: fstrSOLID
24  type(fstr_dynamic) :: fstrDYNAMIC
25  type(fstr_param) :: fstrPARAM
26  type(fstrst_matrix_contact_lagrange) :: fstrMAT
27  type(hecmwst_matrix), optional :: conMAT
28 
29  integer, optional :: iter
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%VELOCITY_type == kbcinitial )return
37 
38  dyn_step = fstrdynamic%i_step
39  flag_u = 2
40 
41  if(dabs(fstrdynamic%ganma) .lt. 1.0e-20) then
42  if( hecmesh%my_rank == 0 ) then
43  write(imsg,*) 'stop due to fstrDYNAMIC%ganma = 0'
44  end if
45  call hecmw_abort( hecmw_comm_get_comm())
46  end if
47 
48  b2 = fstrdynamic%t_delta &
49  *(fstrdynamic%ganma-fstrdynamic%beta)/fstrdynamic%ganma
50  b3 = fstrdynamic%t_delta**2 &
51  *(fstrdynamic%ganma-2.0*fstrdynamic%beta) &
52  /(2.0*fstrdynamic%ganma)
53  b4 = fstrdynamic%t_delta*fstrdynamic%beta/fstrdynamic%ganma
54  c1 = 2.0*fstrdynamic%t_delta
55 
56  ndof = hecmat%NDOF
57 
58  !C=============================C
59  !C-- implicit dynamic analysis
60  !C=============================C
61  if( fstrdynamic%idx_eqa == 1 ) then
62 
63  do ig0 = 1, fstrsolid%VELOCITY_ngrp_tot
64  ig = fstrsolid%VELOCITY_ngrp_ID(ig0)
65  rhs = fstrsolid%VELOCITY_ngrp_val(ig0)
66 
67  call table_dyn(hecmesh, fstrsolid, fstrdynamic, ig0, f_t, flag_u)
68  rhs = rhs * f_t
69  rhs0 = rhs
70 
71  ityp = fstrsolid%VELOCITY_ngrp_type(ig0)
72 
73  idofs = ityp/10
74  idofe = ityp - idofs*10
75 
76  is0 = hecmesh%node_group%grp_index(ig-1) + 1
77  ie0 = hecmesh%node_group%grp_index(ig )
78 
79  do ik = is0, ie0
80  in = hecmesh%node_group%grp_item(ik)
81  do idof = idofs, idofe
82 
83  if( present(iter) ) then ! increment
84  if( iter>1 ) then
85  rhs = 0.d0
86  else
87  rhs = &
88  + b2*fstrdynamic%VEL (ndof*in-(ndof-idof),1) &
89  + b3*fstrdynamic%ACC (ndof*in-(ndof-idof),1) &
90  + b4*rhs0
91  endif
92  else
93  rhs = fstrdynamic%DISP(ndof*in-(ndof-idof),1) &
94  + b2*fstrdynamic%VEL (ndof*in-(ndof-idof),1) &
95  + b3*fstrdynamic%ACC (ndof*in-(ndof-idof),1) &
96  + b4*rhs0
97  endif
98  if(present(conmat)) then
99  call hecmw_mat_ass_bc(hecmat, in, idof, rhs, conmat)
100  else
101  call hecmw_mat_ass_bc(hecmat, in, idof, rhs)
102  endif
103  if( fstr_is_contact_active() .and. fstrparam%contact_algo == kcaslagrange &
104  .and. fstrparam%nlgeom .and. fstrdynamic%idx_resp == 1 ) then
105  if(present(conmat)) then
106  call fstr_mat_ass_bc_contact(conmat,fstrmat,in,idof,rhs)
107  else
108  call fstr_mat_ass_bc_contact(hecmat,fstrmat,in,idof,rhs)
109  endif
110  endif
111 
112  !for output reaction force
113  fstrsolid%REACTION(ndof*(in-1)+idof) = fstrsolid%QFORCE(ndof*(in-1)+idof)
114  enddo
115  enddo
116  enddo
117  !C
118  !C-- end of implicit dynamic analysis
119  !C
120 
121  !C=============================C
122  !C-- explicit dynamic analysis
123  !C=============================C
124  else if( fstrdynamic%idx_eqa == 11 ) then
125  !C
126  do ig0 = 1, fstrsolid%VELOCITY_ngrp_tot
127  ig = fstrsolid%VELOCITY_ngrp_ID(ig0)
128  rhs = fstrsolid%VELOCITY_ngrp_val(ig0)
129 
130  call table_dyn(hecmesh, fstrsolid, fstrdynamic, ig0, f_t, flag_u)
131  rhs = rhs * f_t
132  rhs0 = rhs
133 
134  ityp = fstrsolid%VELOCITY_ngrp_type(ig0)
135 
136  is0 = hecmesh%node_group%grp_index(ig-1) + 1
137  ie0 = hecmesh%node_group%grp_index(ig )
138  idofs = ityp/10
139  idofe = ityp - idofs*10
140 
141  do ik = is0, ie0
142  in = hecmesh%node_group%grp_item(ik)
143  do idof = idofs, idofe
144  rhs = fstrdynamic%DISP(ndof*in-(ndof-idof),3) &
145  + c1*rhs0
146  hecmat%B (ndof*in-(ndof-idof)) = rhs
147  fstrdynamic%VEC1(ndof*in-(ndof-idof)) = 1.0d0
148 
149  !for output reaction force
150  fstrsolid%REACTION(ndof*(in-1)+idof) = fstrsolid%QFORCE(ndof*(in-1)+idof)
151  end do
152  enddo
153  enddo
154  !C
155  !C-- end of explicit dynamic analysis
156  !C
157  end if
158  !
159  return
160  end subroutine dynamic_mat_ass_bc_vl
161 
162 
163  !C***
165  !C***
166  subroutine dynamic_bc_init_vl(hecMESH, hecMAT, fstrSOLID ,fstrDYNAMIC)
167  use m_fstr
168  use m_table_dyn
169 
170  implicit none
171  type(hecmwst_matrix) :: hecmat
172  type(hecmwst_local_mesh) :: hecMESH
173  type(fstr_solid) :: fstrSOLID
174  type(fstr_dynamic) :: fstrDYNAMIC
175 
176  integer(kind=kint) :: NDOF, ig0, ig, ityp, iS0, iE0, ik, in, idofS, idofE, idof
177 
178  integer(kind=kint) :: flag_u
179  real(kind=kreal) :: rhs, f_t
180 
181  if( fstrsolid%VELOCITY_type == kbctransit )return
182 
183  flag_u = 2
184  ndof = hecmat%NDOF
185 
186  do ig0 = 1, fstrsolid%VELOCITY_ngrp_tot
187  ig = fstrsolid%VELOCITY_ngrp_ID(ig0)
188  rhs = fstrsolid%VELOCITY_ngrp_val(ig0)
189 
190 
191  call table_dyn(hecmesh, fstrsolid, fstrdynamic, ig0, f_t, flag_u)
192  rhs = rhs * f_t
193 
194  ityp = fstrsolid%VELOCITY_ngrp_type(ig0)
195 
196  is0 = hecmesh%node_group%grp_index(ig-1) + 1
197  ie0 = hecmesh%node_group%grp_index(ig )
198  idofs = ityp/10
199  idofe = ityp - idofs*10
200 
201  do ik = is0, ie0
202  in = hecmesh%node_group%grp_item(ik)
203 
204  do idof = idofs, idofe
205  fstrdynamic%VEL (ndof*in-(ndof-idof),1) = rhs
206  end do
207  enddo
208  enddo
209 
210  return
211  end subroutine dynamic_bc_init_vl
212 
213  subroutine dynamic_explicit_ass_vl(hecMESH, hecMAT, fstrSOLID, fstrDYNAMIC, iter)
214  use m_fstr
215  use m_table_dyn
218  use mcontact
219 
220  implicit none
221  type(hecmwst_matrix) :: hecmat
222  type(hecmwst_local_mesh) :: hecMESH
223  type(fstr_solid) :: fstrSOLID
224  type(fstr_dynamic) :: fstrDYNAMIC
225  integer, optional :: iter
226 
227  integer(kind=kint) :: ig0, ig, ityp, NDOF, iS0, iE0, ik, in, idofS, idofE, idof
228  integer(kind=kint) :: dyn_step, flag_u
229  real(kind=kreal) :: b2, b3, b4, c1
230  real(kind=kreal) :: rhs, rhs0, f_t
231 
232  if( fstrsolid%VELOCITY_type == kbcinitial )return
233 
234  dyn_step = fstrdynamic%i_step
235  flag_u = 2
236 
237  if(dabs(fstrdynamic%ganma) .lt. 1.0e-20) then
238  if( hecmesh%my_rank == 0 ) then
239  write(imsg,*) 'stop due to fstrDYNAMIC%ganma = 0'
240  end if
241  call hecmw_abort( hecmw_comm_get_comm())
242  end if
243 
244  b2 = fstrdynamic%t_delta &
245  *(fstrdynamic%ganma-fstrdynamic%beta)/fstrdynamic%ganma
246  b3 = fstrdynamic%t_delta**2 &
247  *(fstrdynamic%ganma-2.0*fstrdynamic%beta) &
248  /(2.0*fstrdynamic%ganma)
249  b4 = fstrdynamic%t_delta*fstrdynamic%beta/fstrdynamic%ganma
250  c1 = 2.0*fstrdynamic%t_delta
251 
252  ndof = hecmat%NDOF
253 
254  !C
255  do ig0 = 1, fstrsolid%VELOCITY_ngrp_tot
256  ig = fstrsolid%VELOCITY_ngrp_ID(ig0)
257  rhs = fstrsolid%VELOCITY_ngrp_val(ig0)
258 
259  call table_dyn(hecmesh, fstrsolid, fstrdynamic, ig0, f_t, flag_u)
260  rhs = rhs * f_t
261  rhs0 = rhs
262 
263  ityp = fstrsolid%VELOCITY_ngrp_type(ig0)
264 
265  is0 = hecmesh%node_group%grp_index(ig-1) + 1
266  ie0 = hecmesh%node_group%grp_index(ig )
267  idofs = ityp/10
268  idofe = ityp - idofs*10
269 
270  do ik = is0, ie0
271  in = hecmesh%node_group%grp_item(ik)
272  do idof = idofs, idofe
273  rhs = fstrdynamic%DISP(ndof*in-(ndof-idof),3) &
274  + c1*rhs0
275  hecmat%B(ndof*in-(ndof-idof)) = rhs* fstrdynamic%VEC1(ndof*in-(ndof-idof))
276  ! fstrDYNAMIC%VEC1(NDOF*in-(NDOF-idof)) = 1.0d0
277 
278  !for output reaction force
279  fstrsolid%REACTION(ndof*(in-1)+idof) = fstrsolid%QFORCE(ndof*(in-1)+idof)
280  end do
281  enddo
282  enddo
283 
284  end subroutine dynamic_explicit_ass_vl
285 
286 end module m_dynamic_mat_ass_bc_vl
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 velocity boundary condition in dynamic analysis.
subroutine dynamic_bc_init_vl(hecMESH, hecMAT, fstrSOLID, fstrDYNAMIC)
This function sets initial condition of velocity.
subroutine dynamic_explicit_ass_vl(hecMESH, hecMAT, fstrSOLID, fstrDYNAMIC, iter)
subroutine dynamic_mat_ass_bc_vl(hecMESH, hecMAT, fstrSOLID, fstrDYNAMIC, fstrPARAM, fstrMAT, iter, conMAT)
This subrouitne set velocity boundary condition in dynamic analysis.
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint), parameter imsg
Definition: m_fstr.f90:94
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