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