FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_solver_scaling_33.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 use hecmw_util
10 implicit none
11
12 private
13 real(kind=kreal), private, allocatable :: scale(:)
14
17
18contains
19
20 subroutine hecmw_solver_scaling_fw_33(hecMESH, hecMAT, COMMtime)
21 implicit none
22 type (hecmwst_local_mesh), intent(in) :: hecmesh
23 type (hecmwst_matrix), intent(inout) :: hecmat
24 real(kind=kreal), intent(inout) :: commtime
25 integer(kind=kint) :: n, np, ndof
26 real(kind=kreal), pointer :: d(:), al(:), au(:), b(:)
27 integer(kind=kint), pointer :: inl(:), ial(:), inu(:), iau(:)
28 integer(kind=kint) :: i, k, ip1, ip2, ip3, iq1, iq2, iq3
29 integer(kind=kint) :: isl, iel, isu, ieu, inod
30 real(kind=kreal) :: start_time, end_time
31
32 if (hecmw_mat_get_scaling(hecmat).eq.0) return
33
34 n = hecmat%N
35 np = hecmat%NP
36 ndof = hecmat%NDOF
37 d => hecmat%D
38 al => hecmat%AL
39 au => hecmat%AU
40 inl => hecmat%indexL
41 ial => hecmat%itemL
42 inu => hecmat%indexU
43 iau => hecmat%itemU
44 b => hecmat%B
45
46 allocate(scale(ndof*np))
47
48 do i= 1, n
49 scale(3*i-2)= 1.d0/dsqrt(dabs(d(9*i-8)))
50 scale(3*i-1)= 1.d0/dsqrt(dabs(d(9*i-4)))
51 scale(3*i )= 1.d0/dsqrt(dabs(d(9*i )))
52 enddo
53
54 start_time= hecmw_wtime()
55 call hecmw_update_3_r (hecmesh, scale, hecmesh%n_node)
56 end_time= hecmw_wtime()
57 commtime = commtime + end_time - start_time
58
59 do i= 1, np
60 ip1= 3*i-2
61 ip2= 3*i-1
62 ip3= 3*i
63 d(9*i-8)= d(9*i-8)*scale(ip1)*scale(ip1)
64 d(9*i-7)= d(9*i-7)*scale(ip1)*scale(ip2)
65 d(9*i-6)= d(9*i-6)*scale(ip1)*scale(ip3)
66 d(9*i-5)= d(9*i-5)*scale(ip2)*scale(ip1)
67 d(9*i-4)= d(9*i-4)*scale(ip2)*scale(ip2)
68 d(9*i-3)= d(9*i-3)*scale(ip2)*scale(ip3)
69 d(9*i-2)= d(9*i-2)*scale(ip3)*scale(ip1)
70 d(9*i-1)= d(9*i-1)*scale(ip3)*scale(ip2)
71 d(9*i )= d(9*i )*scale(ip3)*scale(ip3)
72
73 isl= inl(i-1) + 1
74 iel= inl(i )
75 !*voption indep (IAL,AL,SCALE)
76 do k= isl, iel
77 inod= ial(k)
78 iq1= 3*inod - 2
79 iq2= 3*inod - 1
80 iq3= 3*inod
81 al(9*k-8)= al(9*k-8)*scale(ip1)*scale(iq1)
82 al(9*k-7)= al(9*k-7)*scale(ip1)*scale(iq2)
83 al(9*k-6)= al(9*k-6)*scale(ip1)*scale(iq3)
84 al(9*k-5)= al(9*k-5)*scale(ip2)*scale(iq1)
85 al(9*k-4)= al(9*k-4)*scale(ip2)*scale(iq2)
86 al(9*k-3)= al(9*k-3)*scale(ip2)*scale(iq3)
87 al(9*k-2)= al(9*k-2)*scale(ip3)*scale(iq1)
88 al(9*k-1)= al(9*k-1)*scale(ip3)*scale(iq2)
89 al(9*k )= al(9*k )*scale(ip3)*scale(iq3)
90 enddo
91
92 isu= inu(i-1) + 1
93 ieu= inu(i )
94 !*voption indep (IAU,AU,SCALE)
95 do k= isu, ieu
96 inod= iau(k)
97 iq1= 3*inod - 2
98 iq2= 3*inod - 1
99 iq3= 3*inod
100 au(9*k-8)= au(9*k-8)*scale(ip1)*scale(iq1)
101 au(9*k-7)= au(9*k-7)*scale(ip1)*scale(iq2)
102 au(9*k-6)= au(9*k-6)*scale(ip1)*scale(iq3)
103 au(9*k-5)= au(9*k-5)*scale(ip2)*scale(iq1)
104 au(9*k-4)= au(9*k-4)*scale(ip2)*scale(iq2)
105 au(9*k-3)= au(9*k-3)*scale(ip2)*scale(iq3)
106 au(9*k-2)= au(9*k-2)*scale(ip3)*scale(iq1)
107 au(9*k-1)= au(9*k-1)*scale(ip3)*scale(iq2)
108 au(9*k )= au(9*k )*scale(ip3)*scale(iq3)
109 enddo
110 enddo
111 !*voption indep (B,SCALE)
112 do i= 1, n
113 b(3*i-2)= b(3*i-2) * scale(3*i-2)
114 b(3*i-1)= b(3*i-1) * scale(3*i-1)
115 b(3*i )= b(3*i ) * scale(3*i )
116 enddo
117 end subroutine hecmw_solver_scaling_fw_33
118
119 subroutine hecmw_solver_scaling_bk_33(hecMAT)
120 use hecmw_util
121 implicit none
122 type (hecmwst_matrix), intent(inout) :: hecmat
123 integer(kind=kint) :: n, np, ndof
124 real(kind=kreal), pointer :: d(:), al(:), au(:), b(:), x(:)
125 integer(kind=kint), pointer :: inl(:), ial(:), inu(:), iau(:)
126 integer(kind=kint) :: i, k, ip1, ip2, ip3, iq1, iq2, iq3
127 integer(kind=kint) :: isl, iel, isu, ieu, inod
128
129 if (hecmw_mat_get_scaling(hecmat).eq.0) return
130
131 n = hecmat%N
132 np = hecmat%NP
133 ndof = hecmat%NDOF
134 d => hecmat%D
135 al => hecmat%AL
136 au => hecmat%AU
137 inl => hecmat%indexL
138 ial => hecmat%itemL
139 inu => hecmat%indexU
140 iau => hecmat%itemU
141 b => hecmat%B
142 x => hecmat%X
143
144 !*voption indep (X,B,SCALE)
145 do i= 1, n
146 x(3*i-2)= x(3*i-2) * scale(3*i-2)
147 x(3*i-1)= x(3*i-1) * scale(3*i-1)
148 x(3*i )= x(3*i ) * scale(3*i )
149 b(3*i-2)= b(3*i-2) / scale(3*i-2)
150 b(3*i-1)= b(3*i-1) / scale(3*i-1)
151 b(3*i )= b(3*i ) / scale(3*i )
152 enddo
153
154 do i= 1, np
155 ip1= 3*i-2
156 ip2= 3*i-1
157 ip3= 3*i
158 d(9*i-8)= d(9*i-8)/(scale(ip1)*scale(ip1))
159 d(9*i-7)= d(9*i-7)/(scale(ip1)*scale(ip2))
160 d(9*i-6)= d(9*i-6)/(scale(ip1)*scale(ip3))
161 d(9*i-5)= d(9*i-5)/(scale(ip2)*scale(ip1))
162 d(9*i-4)= d(9*i-4)/(scale(ip2)*scale(ip2))
163 d(9*i-3)= d(9*i-3)/(scale(ip2)*scale(ip3))
164 d(9*i-2)= d(9*i-2)/(scale(ip3)*scale(ip1))
165 d(9*i-1)= d(9*i-1)/(scale(ip3)*scale(ip2))
166 d(9*i )= d(9*i )/(scale(ip3)*scale(ip3))
167
168 isl= inl(i-1) + 1
169 iel= inl(i )
170 !*voption indep (IAL,AL,SCALE)
171 do k= isl, iel
172 inod= ial(k)
173 iq1= 3*inod - 2
174 iq2= 3*inod - 1
175 iq3= 3*inod
176 al(9*k-8)= al(9*k-8)/(scale(ip1)*scale(iq1))
177 al(9*k-7)= al(9*k-7)/(scale(ip1)*scale(iq2))
178 al(9*k-6)= al(9*k-6)/(scale(ip1)*scale(iq3))
179 al(9*k-5)= al(9*k-5)/(scale(ip2)*scale(iq1))
180 al(9*k-4)= al(9*k-4)/(scale(ip2)*scale(iq2))
181 al(9*k-3)= al(9*k-3)/(scale(ip2)*scale(iq3))
182 al(9*k-2)= al(9*k-2)/(scale(ip3)*scale(iq1))
183 al(9*k-1)= al(9*k-1)/(scale(ip3)*scale(iq2))
184 al(9*k )= al(9*k )/(scale(ip3)*scale(iq3))
185 enddo
186
187 isu= inu(i-1) + 1
188 ieu= inu(i )
189 !*voption indep (IAU,AU,SCALE)
190 do k= isu, ieu
191 inod= iau(k)
192 iq1= 3*inod - 2
193 iq2= 3*inod - 1
194 iq3= 3*inod
195 au(9*k-8)= au(9*k-8)/(scale(ip1)*scale(iq1))
196 au(9*k-7)= au(9*k-7)/(scale(ip1)*scale(iq2))
197 au(9*k-6)= au(9*k-6)/(scale(ip1)*scale(iq3))
198 au(9*k-5)= au(9*k-5)/(scale(ip2)*scale(iq1))
199 au(9*k-4)= au(9*k-4)/(scale(ip2)*scale(iq2))
200 au(9*k-3)= au(9*k-3)/(scale(ip2)*scale(iq3))
201 au(9*k-2)= au(9*k-2)/(scale(ip3)*scale(iq1))
202 au(9*k-1)= au(9*k-1)/(scale(ip3)*scale(iq2))
203 au(9*k )= au(9*k )/(scale(ip3)*scale(iq3))
204 enddo
205 enddo
206
207 deallocate(scale)
208 end subroutine hecmw_solver_scaling_bk_33
209
integer(kind=kint) function, public hecmw_mat_get_scaling(hecmat)
subroutine, public hecmw_solver_scaling_bk_33(hecmat)
subroutine, public hecmw_solver_scaling_fw_33(hecmesh, hecmat, 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_3_r(hecmesh, val, n)