18 private :: fstr_update_ndforce_mpc
33 integer(kind=kint),
intent(in) :: cstep
34 type(hecmwst_local_mesh),
intent(in) :: hecmesh
35 type(hecmwst_matrix),
intent(inout) :: hecmat
37 type(hecmwst_matrix),
intent(inout),
optional :: conmat
39 integer(kind=kint) :: ndof, idof
40 real(kind=kreal) :: factor
42 factor = fstrsolid%factor(2)
45 do idof=1, hecmesh%n_node* hecmesh%n_dof
46 hecmat%B(idof)=fstrsolid%GL(idof)-fstrsolid%QFORCE(idof)
56 call fstr_update_ndforce_mpc( hecmesh, hecmat%B )
64 call hecmw_update_3_r(hecmesh,hecmat%B,hecmesh%n_node)
65 else if( ndof==2 )
then
66 call hecmw_update_2_r(hecmesh,hecmat%B,hecmesh%n_node)
67 else if( ndof==6 )
then
68 call hecmw_update_m_r(hecmesh,hecmat%B,hecmesh%n_node,6)
73 subroutine fstr_update_ndforce_mpc( hecMESH, B )
75 type(hecmwst_local_mesh),
intent(in) :: hecmesh
76 real(kind=kreal),
intent(inout) :: b(:)
78 integer(kind=kint) ndof, ig0, is0, ie0, ik, in, idof
79 real(kind=kreal) :: rhs, lambda
82 outer:
do ig0=1,hecmesh%mpc%n_mpc
83 is0= hecmesh%mpc%mpc_index(ig0-1)+1
84 ie0= hecmesh%mpc%mpc_index(ig0)
86 if (hecmesh%mpc%mpc_dof(ik) > ndof) cycle outer
89 in = hecmesh%mpc%mpc_item(is0)
90 idof = hecmesh%mpc%mpc_dof(is0)
91 rhs = hecmesh%mpc%mpc_val(is0)
92 lambda = b(ndof*(in-1)+idof)/rhs
95 in = hecmesh%mpc%mpc_item(ik)
96 idof = hecmesh%mpc%mpc_dof(ik)
97 rhs = hecmesh%mpc%mpc_val(ik)
98 b(ndof*(in-1)+idof) = b(ndof*(in-1)+idof) - rhs*lambda
101 end subroutine fstr_update_ndforce_mpc
105 integer(kind=kint),
intent(in) :: cstep
106 type(hecmwst_local_mesh),
intent(in) :: hecmesh
108 real(kind=kreal),
intent(inout) :: b(:)
110 integer(kind=kint) ndof, ig0, ig, ityp, is0, ie0, ik, in, idof1, idof2, idof
111 integer(kind=kint) :: grpid
112 real(kind=kreal) :: rhs
115 fstrsolid%REACTION = 0.d0
117 do ig0= 1, fstrsolid%BOUNDARY_ngrp_tot
118 grpid = fstrsolid%BOUNDARY_ngrp_GRPID(ig0)
120 ig= fstrsolid%BOUNDARY_ngrp_ID(ig0)
121 rhs= fstrsolid%BOUNDARY_ngrp_val(ig0)
122 ityp= fstrsolid%BOUNDARY_ngrp_type(ig0)
123 is0= hecmesh%node_group%grp_index(ig-1) + 1
124 ie0= hecmesh%node_group%grp_index(ig )
126 in = hecmesh%node_group%grp_item(ik)
128 idof2 = ityp - idof1*10
129 if( fstrsolid%BOUNDARY_ngrp_rotID(ig0) > 0 )
then
134 b( ndof*(in-1) + idof ) = 0.d0
136 fstrsolid%REACTION(ndof*(in-1)+idof) = fstrsolid%QFORCE(ndof*(in-1)+idof)
145 real(kind=kreal),
intent(in) :: force(:)
146 type(hecmwst_local_mesh),
intent(in) :: hecmesh
155 real(kind=kreal),
intent(in) :: force(:), displacement(:)
156 type(hecmwst_local_mesh),
intent(in) :: hecmesh
159 call hecmw_innerproduct_r(hecmesh, ndof, force, displacement,
fstr_get_energy)
166 type(hecmwst_local_mesh),
intent(in) :: hecmesh
167 type(hecmwst_matrix),
intent(in) :: hecmat
170 character(len=13) :: flag
171 real(kind=kreal) :: tmp1, tmp2, bi
172 integer :: i, i0, ndof
173 if( flag==
'residualForce' )
then
175 call hecmw_innerproduct_r(hecmesh,ndof,hecmat%B,hecmat%B,tmp1)
177 i0 = hecmesh%n_node*ndof
178 do i=1,fstrmat%num_lagrange
182 call hecmw_allreduce_r1(hecmesh,tmp2,hecmw_sum)
184 elseif( flag==
' force' )
then
194 type(hecmwst_matrix),
intent(in) :: hecmat
196 type(hecmwst_matrix),
intent(in) :: conmat
197 type(hecmwst_local_mesh),
intent(in) :: hecmesh
199 real(kind=kreal) :: rhsb
200 integer(kind=kint) :: i,j,n,i0,ndof,n_loc,nndof
201 integer(kind=kint) :: offset, pid, lid
202 integer(kind=kint),
allocatable :: displs(:)
203 real(kind=kreal),
allocatable :: rhs_con_all(:), rhs_con(:)
206 nndof = hecmat%N*ndof
207 n_loc = nndof + fstrmat%num_lagrange
208 allocate(displs(0:
nprocs))
211 call hecmw_allreduce_i(hecmesh, displs,
nprocs+1, hecmw_sum)
213 displs(i) = displs(i-1) + displs(i)
217 allocate(rhs_con_all(n))
219 rhs_con_all(i) = 0.0d0
221 do i= hecmat%N+1,hecmat%NP
223 pid = hecmesh%node_ID(i*2)
224 lid = hecmesh%node_ID(i*2-1)
225 i0 = displs(pid) + (lid-1)*ndof
226 if((i0 < 0.or.i0 > n))
then
229 if(conmat%b((i-1)*ndof+j) /= 0.0d0)
then
230 print *,
myrank,
'i0',i,
'conMAT%b',conmat%b((i-1)*ndof+j)
235 rhs_con_all(i0+1:i0+ndof) = conmat%b((i-1)*ndof+1:i*ndof)
239 call hecmw_allreduce_r(hecmesh, rhs_con_all, n, hecmw_sum)
240 allocate(rhs_con(n_loc))
242 rhs_con(i) = rhs_con_all(offset+i) + hecmat%B(i) + conmat%B(i)
244 deallocate(rhs_con_all)
246 do i=1,fstrmat%num_lagrange
247 rhs_con(nndof+i) = conmat%B(i0+i)
250 rhsb = dot_product(rhs_con, rhs_con)
251 call hecmw_allreduce_r1(hecmesh, rhsb, hecmw_sum)
260 type(hecmwst_matrix),
intent(in) :: hecmat
262 type(hecmwst_local_mesh),
intent(in) :: hecmesh
263 real(kind=kreal) :: rhsx
264 integer(kind=kint) :: nndof, npndof, i
266 nndof = hecmat%N * hecmat%NDOF
267 npndof = hecmat%NP * hecmat%NDOF
270 rhsx = rhsx + hecmat%X(i) ** 2
272 do i=1,fstrmat%num_lagrange
273 rhsx = rhsx + hecmat%X(npndof+i) ** 2
275 call hecmw_allreduce_r1(hecmesh, rhsx, hecmw_sum)
This module provides function to calcualte residual of nodal force.
subroutine, public fstr_update_ndforce_spc(cstep, hecMESH, fstrSOLID, B)
subroutine, public fstr_update_ndforce(cstep, hecMESH, hecMAT, fstrSOLID, conMAT)
real(kind=kreal) function, public fstr_get_x_norm_contact(hecMAT, fstrMAT, hecMESH)
real(kind=kreal) function, public fstr_get_residual(force, hecMESH)
Calculate magnitude of a real vector.
real(kind=kreal) function, public fstr_get_norm_para_contact(hecMAT, fstrMAT, conMAT, hecMESH)
real(kind=kreal) function fstr_get_energy(force, displacement, hecMESH)
Calculate magnitude of a real vector.
real(kind=kreal) function, public fstr_get_norm_contact(flag, hecMESH, hecMAT, fstrSOLID, fstrMAT)
Calculate square norm.
This module provides functions to deal with spring force.
subroutine fstr_update_ndforce_spring(cstep, hecMESH, fstrSOLID, B)
This module defined coomon data and basic structures for analysis.
integer(kind=kint) myrank
PARALLEL EXECUTION.
integer(kind=kint) nprocs
logical function fstr_isboundaryactive(fstrSOLID, nbc, cstep)
This subroutine read in used-defined loading tangent.
subroutine uresidual(cstep, factor, residual)
This subroutine take consider of user-defined external loading.