FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_solver_las_11.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!-------------------------------------------------------------------------------
5
7
8 private
9
10 public :: hecmw_matvec_11
11 public :: hecmw_matresid_11
12 public :: hecmw_rel_resid_l2_11
13
14contains
15
16 !C
17 !C***
18 !C*** hecmw_matvec_11
19 !C***
20 !C
21 subroutine hecmw_matvec_11 (hecMESH, hecMAT, X, Y, n, COMMtime)
22 use hecmw_util
24
25 implicit none
26 integer(kind=kint):: n
27 real(kind=kreal), dimension(n) :: x, y
28 type (hecmwst_matrix) :: hecmat
29 type (hecmwst_local_mesh) :: hecmesh
30 real(kind=kreal), optional :: commtime
31
32 real(kind=kreal) :: start_time, end_time
33 integer(kind=kint) :: i, j, js, je, in
34 real(kind=kreal) :: yv
35
36 start_time= hecmw_wtime()
37 call hecmw_update_1_r (hecmesh, x, hecmat%NP)
38 end_time= hecmw_wtime()
39 if (present(commtime)) commtime = commtime + end_time - start_time
40
41 do i= 1, hecmat%N
42 yv= hecmat%D(i) * x(i)
43 js= hecmat%indexL(i-1) + 1
44 je= hecmat%indexL(i )
45 do j= js, je
46 in= hecmat%itemL(j)
47 yv= yv + hecmat%AL(j) * x(in)
48 enddo
49 js= hecmat%indexU(i-1) + 1
50 je= hecmat%indexU(i )
51 do j= js, je
52 in= hecmat%itemU(j)
53 yv= yv + hecmat%AU(j) * x(in)
54 enddo
55 y(i)= yv
56 enddo
57
58 end subroutine hecmw_matvec_11
59
60 !C
61 !C***
62 !C*** hecmw_matresid_11
63 !C***
64 !C
65 subroutine hecmw_matresid_11 (hecMESH, hecMAT, X, B, R, COMMtime)
66 use hecmw_util
67
68 implicit none
69 real(kind=kreal) :: x(:), b(:), r(:)
70 type (hecmwst_matrix) :: hecmat
71 type (hecmwst_local_mesh) :: hecmesh
72 real(kind=kreal), optional :: commtime
73
74 integer(kind=kint) :: i
75 real(kind=kreal) :: tcomm
76
77 tcomm = 0.d0
78 call hecmw_matvec_11 (hecmesh, hecmat, x, r, hecmat%NP, tcomm)
79 if (present(commtime)) commtime = commtime + tcomm
80
81 do i = 1, hecmat%N
82 r(i) = b(i) - r(i)
83 enddo
84
85 end subroutine hecmw_matresid_11
86
87 !C
88 !C***
89 !C*** hecmw_rel_resid_L2_11
90 !C***
91 !C
92 function hecmw_rel_resid_l2_11 (hecMESH, hecMAT, COMMtime)
93 use hecmw_util
95
96 implicit none
97 real(kind=kreal) :: hecmw_rel_resid_l2_11
98 type ( hecmwst_local_mesh ), intent(in) :: hecmesh
99 type ( hecmwst_matrix ), intent(in) :: hecmat
100 real(kind=kreal), optional :: commtime
101
102 real(kind=kreal) :: r(hecmat%NDOF*hecmat%NP)
103 real(kind=kreal) :: bnorm2, rnorm2
104 real(kind=kreal) :: tcomm
105
106 tcomm = 0.d0
107 call hecmw_innerproduct_r(hecmesh, hecmat%NDOF, hecmat%B, hecmat%B, bnorm2, tcomm)
108 if (bnorm2 == 0.d0) then
109 bnorm2 = 1.d0
110 endif
111 call hecmw_matresid_11(hecmesh, hecmat, hecmat%X, hecmat%B, r, tcomm)
112 call hecmw_innerproduct_r(hecmesh, hecmat%NDOF, r, r, rnorm2, tcomm)
113 if (present(commtime)) commtime = commtime + tcomm
114 hecmw_rel_resid_l2_11 = sqrt(rnorm2 / bnorm2)
115
116 end function hecmw_rel_resid_l2_11
117
118end module hecmw_solver_las_11
real(kind=kreal) function, public hecmw_rel_resid_l2_11(hecmesh, hecmat, commtime)
subroutine, public hecmw_matresid_11(hecmesh, hecmat, x, b, r, commtime)
subroutine, public hecmw_matvec_11(hecmesh, hecmat, x, y, n, commtime)
subroutine hecmw_innerproduct_r(hecmesh, ndof, x, y, sum, commtime)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
real(kind=kreal) function hecmw_wtime()
subroutine hecmw_update_1_r(hecmesh, val, n)