FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_mat_resid_contact.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 use m_fstr
8 implicit none
9
10 private
14
15contains
16
17 subroutine fstr_get_residual_contact(hecMESH, hecMAT, fstrMAT, r, rlag)
18 implicit none
19 type(hecmwst_local_mesh) :: hecmesh
20 type(hecmwst_matrix) :: hecmat
21 type(fstrst_matrix_contact_lagrange), intent(in) :: fstrmat
22 real(kind=kreal), intent(out) :: r(:)
23 real(kind=kreal), intent(out) :: rlag(:)
24 real(kind=kreal) :: tcomm, sum
25 integer(kind=kint) :: i, idof, ii, i0, ls, le, l, loc0, ll, j, loc, l0, k
26 integer(kind=kint) :: ndof, npndof
27 !! r = b - K x - Ulag lag
28 !! rlag = blag - Llag x
29 ndof=hecmat%NDOF
30 npndof=hecmat%NP*ndof
31 do i=1,npndof
32 r(i)=0.d0
33 enddo
34 ! r = b - K x
35 tcomm = 0.d0
36 call hecmw_matresid(hecmesh, hecmat, hecmat%X, hecmat%B, r, tcomm)
37 if (fstrmat%num_lagrange > 0) then
38 ! r = r - Ulag lag
39 do i=1,hecmat%N
40 i0=(i-1)*ndof
41 do idof=1,ndof
42 ii=i0+idof
43 ls=fstrmat%indexU_lagrange(i-1)+1
44 le=fstrmat%indexU_lagrange(i)
45 sum=0.d0
46 do l=ls,le
47 loc=(l-1)*ndof+idof
48 ll=fstrmat%itemU_lagrange(l)
49 !sum=sum+fstrMAT%AU_lagrange(loc)*fstrMAT%Lagrange(ll)
50 sum=sum+fstrmat%AU_lagrange(loc)*hecmat%X(npndof+ll)
51 enddo
52 r(ii)=r(ii)-sum
53 enddo
54 enddo
55 ! rlag = blag - Llag x
56 do i=1,fstrmat%num_lagrange
57 ls=fstrmat%indexL_lagrange(i-1)+1
58 le=fstrmat%indexL_lagrange(i)
59 sum=0.d0
60 do l=ls,le
61 loc0=(l-1)*ndof
62 j=fstrmat%itemL_lagrange(l)
63 l0=(j-1)*ndof
64 do k=1,ndof
65 loc=loc0+k
66 ll=l0+k
67 sum=sum+fstrmat%AL_lagrange(loc)*hecmat%X(ll)
68 enddo
69 enddo
70 rlag(i)=hecmat%B(npndof+i)-sum
71 enddo
72 end if
73 end subroutine fstr_get_residual_contact
74
76 function fstr_get_rel_resid_l2_contact(hecMESH, hecMAT, fstrMAT)
77 implicit none
78 real(kind=kreal) :: fstr_get_rel_resid_l2_contact
79 type(hecmwst_local_mesh) :: hecmesh
80 type(hecmwst_matrix) :: hecmat
81 type(fstrst_matrix_contact_lagrange), intent(in) :: fstrmat
82 real(kind=kreal), allocatable :: r(:)
83 real(kind=kreal), allocatable :: rlag(:)
84 real(kind=kreal) :: bnorm2, rnorm2
85 real(kind=kreal) :: rlagnorm2
86 real(kind=kreal) :: tcomm
87 integer(kind=kint) :: i
88 allocate(r(hecmat%NDOF*hecmat%NP),rlag(fstrmat%num_lagrange))
89 ! |b|^2
90 call hecmw_innerproduct_r(hecmesh, hecmat%NDOF, hecmat%B, hecmat%B, bnorm2, tcomm)
91 call fstr_get_residual_contact(hecmesh, hecmat, fstrmat, r, rlag)
92 ! |r|^2
93 call hecmw_innerproduct_r(hecmesh, hecmat%NDOF, r, r, rnorm2, tcomm)
94 ! |rlag|^2
95 rlagnorm2=0.d0
96 do i=1,fstrmat%num_lagrange
97 rlagnorm2=rlagnorm2+rlag(i)*rlag(i)
98 enddo
99 deallocate(r,rlag)
100 call hecmw_allreduce_r1(hecmesh, rlagnorm2, hecmw_sum)
101 ! |r_total|^2 = |r|^2 + |rlag|^2
102 if (fstrmat%num_lagrange > 0) rnorm2=rnorm2+rlagnorm2
103 ! |r_total| / |b|
104 fstr_get_rel_resid_l2_contact = sqrt(rnorm2 / bnorm2)
106
108 function fstr_get_resid_max_contact(hecMESH, hecMAT, fstrMAT)
109 implicit none
110 real(kind=kreal) :: fstr_get_resid_max_contact
111 type(hecmwst_local_mesh) :: hecmesh
112 type(hecmwst_matrix) :: hecmat
113 type(fstrst_matrix_contact_lagrange), intent(in) :: fstrmat
114 real(kind=kreal), allocatable :: r(:)
115 real(kind=kreal), allocatable :: rlag(:)
116 real(kind=kreal) :: rmax, rlagmax
117 real(kind=kreal) :: tcomm
118 allocate(r(hecmat%NDOF*hecmat%NP),rlag(fstrmat%num_lagrange))
119 call fstr_get_residual_contact(hecmesh, hecmat, fstrmat, r, rlag)
120 rmax = maxval(dabs(r))
121 if (fstrmat%num_lagrange > 0) then
122 rlagmax = maxval(dabs(rlag))
123 if (rlagmax > rmax) rmax = rlagmax
124 endif
125 deallocate(r,rlag)
126 call hecmw_allreduce_r1(hecmesh, rmax, hecmw_max)
128 end function fstr_get_resid_max_contact
129
This module provides functions of reconstructing.
subroutine, public fstr_get_residual_contact(hecmesh, hecmat, fstrmat, r, rlag)
real(kind=kreal) function, public fstr_get_rel_resid_l2_contact(hecmesh, hecmat, fstrmat)
This function calculates relative L2 residual.
real(kind=kreal) function, public fstr_get_resid_max_contact(hecmesh, hecmat, fstrmat)
This function calculates maximum residual.
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
Structure for Lagrange multiplier-related part of stiffness matrix (Lagrange multiplier-related matri...