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