FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
fstr_AddBC.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  implicit none
9 contains
10 
12  !------------------------------------------------------------------------------------------*
13  subroutine fstr_addbc(cstep,hecMESH,hecMAT,fstrSOLID,fstrPARAM,fstrMAT,iter,conMAT)
14  !------------------------------------------------------------------------------------------*
15  use m_fstr
18  use mcontact
19  use m_static_lib_1d
20  use m_utilities
21  integer, intent(in) :: cstep
22  type(hecmwst_local_mesh) :: hecMESH
23  type(hecmwst_matrix) :: hecMAT
24  type(fstr_solid) :: fstrSOLID
25  type(fstr_param) :: fstrPARAM
26  type(fstrst_matrix_contact_lagrange) :: fstrMAT
27  integer(kind=kint) :: iter
28  type(hecmwst_matrix), optional :: conMAT
29 
30  integer(kind=kint) :: ig0, ig, ityp, idofS, idofE, idof, iS0, iE0, ik, in
31  real(kind=kreal) :: rhs, factor
32  integer(kind=kint) :: ndof, grpid
33 
34  !for rotation
35  integer(kind=kint) :: n_rot, rid
36  type(trotinfo) :: rinfo
37  real(kind=kreal) :: ccoord(3), cdiff(3), cdiff0(3)
38  real(kind=kreal) :: cdisp(3), cddisp(3)
39 
40  !
41  ndof = hecmat%NDOF
42  factor = fstrsolid%FACTOR(2)-fstrsolid%FACTOR(1)
43 
44  if( cstep<=fstrsolid%nstep_tot .and. fstrsolid%step_ctrl(cstep)%solution==stepvisco ) then
45  factor = 0.d0
46  if( fstrsolid%FACTOR(1) < 1.d-10 ) factor = 1.d0
47  endif
48  if( iter>1 ) factor=0.d0
49 
50  n_rot = fstrsolid%BOUNDARY_ngrp_rot
51  if( n_rot > 0 ) call fstr_rotinfo_init(n_rot, rinfo)
52 
53  ! ----- Prescibed displacement Boundary Conditions
54  do ig0 = 1, fstrsolid%BOUNDARY_ngrp_tot
55  grpid = fstrsolid%BOUNDARY_ngrp_GRPID(ig0)
56  if( .not. fstr_isboundaryactive( fstrsolid, grpid, cstep ) ) cycle
57  ig = fstrsolid%BOUNDARY_ngrp_ID(ig0)
58  rhs = fstrsolid%BOUNDARY_ngrp_val(ig0)
59  !
60  rhs= rhs*factor
61  !
62  ityp = fstrsolid%BOUNDARY_ngrp_type(ig0)
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  if( fstrsolid%BOUNDARY_ngrp_rotID(ig0) > 0 ) then ! setup rotation information
70  rid = fstrsolid%BOUNDARY_ngrp_rotID(ig0)
71  if( .not. rinfo%conds(rid)%active ) then
72  rinfo%conds(rid)%active = .true.
73  rinfo%conds(rid)%center_ngrp_id = fstrsolid%BOUNDARY_ngrp_centerID(ig0)
74  rinfo%conds(rid)%torque_ngrp_id = ig
75  endif
76  do idof=idofs,idofe
77  if( idof>ndof ) then
78  rinfo%conds(rid)%vec(idof-ndof) = rhs
79  else
80  rinfo%conds(rid)%vec(idof) = rhs
81  endif
82  enddo
83  cycle
84  endif
85 
86  !
87  do ik = is0, ie0
88  in = hecmesh%node_group%grp_item(ik)
89  !
90  do idof = idofs, idofe
91  if(present(conmat)) then
92  call hecmw_mat_ass_bc(hecmat, in, idof, rhs, conmat)
93  else
94  call hecmw_mat_ass_bc(hecmat, in, idof, rhs)
95  endif
96  if( fstr_is_contact_active() .and. fstrparam%solution_type == kststatic &
97  .and. fstrparam%contact_algo == kcaslagrange ) then
98  if(present(conmat)) then
99  call fstr_mat_ass_bc_contact(conmat,fstrmat,in,idof,rhs)
100  else
101  call fstr_mat_ass_bc_contact(hecmat,fstrmat,in,idof,rhs)
102  endif
103  endif
104 
105  enddo
106  enddo
107  enddo
108 
109  !Apply rotational boundary condition
110  do rid = 1, n_rot
111  if( .not. rinfo%conds(rid)%active ) cycle
112  cdiff = 0.d0
113  cdiff0 = 0.d0
114  cddisp = 0.d0
115 
116  if( factor > 0.d0 ) then
117  ig = rinfo%conds(rid)%center_ngrp_id
118  do idof = 1, ndof
119  ccoord(idof) = hecmw_ngrp_get_totalvalue(hecmesh, ig, ndof, idof, hecmesh%node)
120  cdisp(idof) = hecmw_ngrp_get_totalvalue(hecmesh, ig, ndof, idof, fstrsolid%unode)
121  cddisp(idof) = hecmw_ngrp_get_totalvalue(hecmesh, ig, ndof, idof, hecmat%B)
122  enddo
123  ccoord(1:ndof) = ccoord(1:ndof) + cdisp(1:ndof)
124  endif
125 
126  ig = rinfo%conds(rid)%torque_ngrp_id
127  is0 = hecmesh%node_group%grp_index(ig-1) + 1
128  ie0 = hecmesh%node_group%grp_index(ig )
129  do ik = is0, ie0
130  in = hecmesh%node_group%grp_item(ik)
131  if( factor > 0.d0 ) then
132  cdiff0(1:ndof) = hecmesh%node(ndof*(in-1)+1:ndof*in)+fstrsolid%unode(ndof*(in-1)+1:ndof*in)-ccoord(1:ndof)
133  cdiff(1:ndof) = cdiff0(1:ndof)
134  call rotate_3dvector_by_rodrigues_formula(rinfo%conds(rid)%vec(1:ndof),cdiff(1:ndof))
135  endif
136  do idof = 1, ndof
137  rhs = cdiff(idof)-cdiff0(idof)+cddisp(idof)
138  if(present(conmat)) then
139  call hecmw_mat_ass_bc(hecmat, in, idof, rhs, conmat)
140  else
141  call hecmw_mat_ass_bc(hecmat, in, idof, rhs)
142  endif
143  if( fstr_is_contact_active() .and. fstrparam%solution_type == kststatic &
144  .and. fstrparam%contact_algo == kcaslagrange ) then
145  if(present(conmat)) then
146  call fstr_mat_ass_bc_contact(conmat,fstrmat,in,idof,rhs)
147  else
148  call fstr_mat_ass_bc_contact(hecmat,fstrmat,in,idof,rhs)
149  endif
150  endif
151  enddo
152  enddo
153  enddo
154  if( n_rot > 0 ) call fstr_rotinfo_finalize(rinfo)
155 
156  !
157  ! ------ Truss element Diagonal Modification
158  call truss_diag_modify(hecmat,hecmesh)
159 
160  !
161  ! ------ Equation boundary conditions
162  do ig0=1,fstrsolid%n_fix_mpc
163  if( fstrsolid%mpc_const(ig0) == 0.d0 ) cycle
164  ! we need to confirm if it is active in curr step here
165  rhs = fstrsolid%mpc_const(ig0)*factor
166  hecmesh%mpc%mpc_const(ig0) = rhs
167  enddo
168 
169  end subroutine fstr_addbc
170 
171 end module m_fstr_addbc
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 provides a function to deal with prescribed displacement.
Definition: fstr_AddBC.f90:7
subroutine fstr_addbc(cstep, hecMESH, hecMAT, fstrSOLID, fstrPARAM, fstrMAT, iter, conMAT)
Add Essential Boundary Conditions.
Definition: fstr_AddBC.f90:14
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 kststatic
Definition: m_fstr.f90:37
logical function fstr_isboundaryactive(fstrSOLID, nbc, cstep)
Definition: m_fstr.f90:997
This module provide common functions of 3D truss elements.
subroutine truss_diag_modify(hecMAT, hecMESH)
This module provides aux functions.
Definition: utilities.f90:6
subroutine rotate_3dvector_by_rodrigues_formula(r, v)
Definition: utilities.f90:515
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...
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138