FrontISTR  5.7.0
Large-scale structural analysis program with finit element method
hecmw_solver_scaling_44.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 
18 contains
19 
20  subroutine hecmw_solver_scaling_fw_44(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, ip4, iq1, iq2, iq3, iq4
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(4*i-3)= 1.d0/dsqrt(dabs(d(16*i-15)))
50  scale(4*i-2)= 1.d0/dsqrt(dabs(d(16*i-10)))
51  scale(4*i-1)= 1.d0/dsqrt(dabs(d(16*i- 5)))
52  scale(4*i )= 1.d0/dsqrt(dabs(d(16*i )))
53  enddo
54 
55  start_time= hecmw_wtime()
56  call hecmw_update_r (hecmesh, scale, hecmesh%n_node, ndof)
57  end_time= hecmw_wtime()
58  commtime = commtime + end_time - start_time
59 
60  do i= 1, np
61  ip1= 4*i-3
62  ip2= 4*i-2
63  ip3= 4*i-1
64  ip4= 4*i
65  d(16*i-15)= d(16*i-15)*scale(ip1)*scale(ip1)
66  d(16*i-14)= d(16*i-14)*scale(ip1)*scale(ip2)
67  d(16*i-13)= d(16*i-13)*scale(ip1)*scale(ip3)
68  d(16*i-12)= d(16*i-12)*scale(ip1)*scale(ip4)
69  d(16*i-11)= d(16*i-11)*scale(ip2)*scale(ip1)
70  d(16*i-10)= d(16*i-10)*scale(ip2)*scale(ip2)
71  d(16*i- 9)= d(16*i- 9)*scale(ip2)*scale(ip3)
72  d(16*i- 8)= d(16*i- 8)*scale(ip2)*scale(ip4)
73  d(16*i- 7)= d(16*i- 7)*scale(ip3)*scale(ip1)
74  d(16*i- 6)= d(16*i- 6)*scale(ip3)*scale(ip2)
75  d(16*i- 5)= d(16*i- 5)*scale(ip3)*scale(ip3)
76  d(16*i- 4)= d(16*i- 4)*scale(ip3)*scale(ip4)
77  d(16*i- 3)= d(16*i- 3)*scale(ip4)*scale(ip1)
78  d(16*i- 2)= d(16*i- 2)*scale(ip4)*scale(ip2)
79  d(16*i- 1)= d(16*i- 1)*scale(ip4)*scale(ip3)
80  d(16*i- 0)= d(16*i- 0)*scale(ip4)*scale(ip4)
81 
82  isl= inl(i-1) + 1
83  iel= inl(i )
84  !*voption indep (IAL,AL,SCALE)
85  do k= isl, iel
86  inod= ial(k)
87  iq1= 4*inod - 3
88  iq2= 4*inod - 2
89  iq3= 4*inod - 1
90  iq4= 4*inod
91  al(16*k-15)= al(16*k-15)*scale(ip1)*scale(iq1)
92  al(16*k-14)= al(16*k-14)*scale(ip1)*scale(iq2)
93  al(16*k-13)= al(16*k-13)*scale(ip1)*scale(iq3)
94  al(16*k-12)= al(16*k-12)*scale(ip1)*scale(iq4)
95  al(16*k-11)= al(16*k-11)*scale(ip2)*scale(iq1)
96  al(16*k-10)= al(16*k-10)*scale(ip2)*scale(iq2)
97  al(16*k- 9)= al(16*k- 9)*scale(ip2)*scale(iq3)
98  al(16*k- 8)= al(16*k- 8)*scale(ip2)*scale(iq4)
99  al(16*k- 7)= al(16*k- 7)*scale(ip3)*scale(iq1)
100  al(16*k- 6)= al(16*k- 6)*scale(ip3)*scale(iq2)
101  al(16*k- 5)= al(16*k- 5)*scale(ip3)*scale(iq3)
102  al(16*k- 4)= al(16*k- 4)*scale(ip3)*scale(iq4)
103  al(16*k- 3)= al(16*k- 3)*scale(ip4)*scale(iq1)
104  al(16*k- 2)= al(16*k- 2)*scale(ip4)*scale(iq2)
105  al(16*k- 1)= al(16*k- 1)*scale(ip4)*scale(iq3)
106  al(16*k- 0)= al(16*k- 0)*scale(ip4)*scale(iq4)
107  enddo
108 
109  isu= inu(i-1) + 1
110  ieu= inu(i )
111  !*voption indep (IAU,AU,SCALE)
112  do k= isu, ieu
113  inod= iau(k)
114  iq1= 4*inod - 3
115  iq2= 4*inod - 2
116  iq3= 4*inod - 1
117  iq4= 4*inod
118  au(16*k-15)= au(16*k-15)*scale(ip1)*scale(iq1)
119  au(16*k-14)= au(16*k-14)*scale(ip1)*scale(iq2)
120  au(16*k-13)= au(16*k-13)*scale(ip1)*scale(iq3)
121  au(16*k-12)= au(16*k-12)*scale(ip1)*scale(iq4)
122  au(16*k-11)= au(16*k-11)*scale(ip2)*scale(iq1)
123  au(16*k-10)= au(16*k-10)*scale(ip2)*scale(iq2)
124  au(16*k- 9)= au(16*k- 9)*scale(ip2)*scale(iq3)
125  au(16*k- 8)= au(16*k- 8)*scale(ip2)*scale(iq4)
126  au(16*k- 7)= au(16*k- 7)*scale(ip3)*scale(iq1)
127  au(16*k- 6)= au(16*k- 6)*scale(ip3)*scale(iq2)
128  au(16*k- 5)= au(16*k- 5)*scale(ip3)*scale(iq3)
129  au(16*k- 4)= au(16*k- 4)*scale(ip3)*scale(iq4)
130  au(16*k- 3)= au(16*k- 3)*scale(ip4)*scale(iq1)
131  au(16*k- 2)= au(16*k- 2)*scale(ip4)*scale(iq2)
132  au(16*k- 1)= au(16*k- 1)*scale(ip4)*scale(iq3)
133  au(16*k- 0)= au(16*k- 0)*scale(ip4)*scale(iq4)
134  enddo
135  enddo
136  !*voption indep (B,SCALE)
137  do i= 1, n
138  b(4*i-3)= b(4*i-3) * scale(4*i-3)
139  b(4*i-2)= b(4*i-2) * scale(4*i-2)
140  b(4*i-1)= b(4*i-1) * scale(4*i-1)
141  b(4*i )= b(4*i ) * scale(4*i )
142  enddo
143  end subroutine hecmw_solver_scaling_fw_44
144 
145  subroutine hecmw_solver_scaling_bk_44(hecMAT)
147  implicit none
148  type (hecmwst_matrix), intent(inout) :: hecmat
149  integer(kind=kint) :: n, np, ndof
150  real(kind=kreal), pointer :: d(:), al(:), au(:), b(:), x(:)
151  integer(kind=kint), pointer :: inl(:), ial(:), inu(:), iau(:)
152  integer(kind=kint) :: i, k, ip1, ip2, ip3, ip4, iq1, iq2, iq3, iq4
153  integer(kind=kint) :: isl, iel, isu, ieu, inod
154 
155  if (hecmw_mat_get_scaling(hecmat).eq.0) return
156 
157  n = hecmat%N
158  np = hecmat%NP
159  ndof = hecmat%NDOF
160  d => hecmat%D
161  al => hecmat%AL
162  au => hecmat%AU
163  inl => hecmat%indexL
164  ial => hecmat%itemL
165  inu => hecmat%indexU
166  iau => hecmat%itemU
167  b => hecmat%B
168  x => hecmat%X
169 
170  !*voption indep (X,B,SCALE)
171  do i= 1, n
172  x(4*i-3)= x(4*i-3) * scale(4*i-3)
173  x(4*i-2)= x(4*i-2) * scale(4*i-2)
174  x(4*i-1)= x(4*i-1) * scale(4*i-1)
175  x(4*i )= x(4*i ) * scale(4*i )
176  b(4*i-3)= b(4*i-3) / scale(4*i-3)
177  b(4*i-2)= b(4*i-2) / scale(4*i-2)
178  b(4*i-1)= b(4*i-1) / scale(4*i-1)
179  b(4*i )= b(4*i ) / scale(4*i )
180  enddo
181 
182  do i= 1, np
183  ip1= 4*i-3
184  ip2= 4*i-2
185  ip3= 4*i-1
186  ip4= 4*i
187  d(16*i-15)= d(16*i-15)/(scale(ip1)*scale(ip1))
188  d(16*i-14)= d(16*i-14)/(scale(ip1)*scale(ip2))
189  d(16*i-13)= d(16*i-13)/(scale(ip1)*scale(ip3))
190  d(16*i-12)= d(16*i-12)/(scale(ip1)*scale(ip4))
191  d(16*i-11)= d(16*i-11)/(scale(ip2)*scale(ip1))
192  d(16*i-10)= d(16*i-10)/(scale(ip2)*scale(ip2))
193  d(16*i- 9)= d(16*i- 9)/(scale(ip2)*scale(ip3))
194  d(16*i- 8)= d(16*i- 8)/(scale(ip2)*scale(ip4))
195  d(16*i- 7)= d(16*i- 7)/(scale(ip3)*scale(ip1))
196  d(16*i- 6)= d(16*i- 6)/(scale(ip3)*scale(ip2))
197  d(16*i- 5)= d(16*i- 5)/(scale(ip3)*scale(ip3))
198  d(16*i- 4)= d(16*i- 4)/(scale(ip3)*scale(ip4))
199  d(16*i- 3)= d(16*i- 3)/(scale(ip4)*scale(ip1))
200  d(16*i- 2)= d(16*i- 2)/(scale(ip4)*scale(ip2))
201  d(16*i- 1)= d(16*i- 1)/(scale(ip4)*scale(ip3))
202  d(16*i- 0)= d(16*i- 0)/(scale(ip4)*scale(ip4))
203  isl= inl(i-1) + 1
204  iel= inl(i )
205  !*voption indep (IAL,AL,SCALE)
206  do k= isl, iel
207  inod= ial(k)
208  iq1= 4*inod - 3
209  iq2= 4*inod - 2
210  iq3= 4*inod - 1
211  iq4= 4*inod
212  al(16*k-15)= al(16*k-15)/(scale(ip1)*scale(iq1))
213  al(16*k-14)= al(16*k-14)/(scale(ip1)*scale(iq2))
214  al(16*k-13)= al(16*k-13)/(scale(ip1)*scale(iq3))
215  al(16*k-12)= al(16*k-12)/(scale(ip1)*scale(iq4))
216  al(16*k-11)= al(16*k-11)/(scale(ip2)*scale(iq1))
217  al(16*k-10)= al(16*k-10)/(scale(ip2)*scale(iq2))
218  al(16*k- 9)= al(16*k- 9)/(scale(ip2)*scale(iq3))
219  al(16*k- 8)= al(16*k- 8)/(scale(ip2)*scale(iq4))
220  al(16*k- 7)= al(16*k- 7)/(scale(ip3)*scale(iq1))
221  al(16*k- 6)= al(16*k- 6)/(scale(ip3)*scale(iq2))
222  al(16*k- 5)= al(16*k- 5)/(scale(ip3)*scale(iq3))
223  al(16*k- 4)= al(16*k- 4)/(scale(ip3)*scale(iq4))
224  al(16*k- 3)= al(16*k- 3)/(scale(ip4)*scale(iq1))
225  al(16*k- 2)= al(16*k- 2)/(scale(ip4)*scale(iq2))
226  al(16*k- 1)= al(16*k- 1)/(scale(ip4)*scale(iq3))
227  al(16*k- 0)= al(16*k- 0)/(scale(ip4)*scale(iq4))
228  enddo
229 
230  isu= inu(i-1) + 1
231  ieu= inu(i )
232  !*voption indep (IAU,AU,SCALE)
233  do k= isu, ieu
234  inod= iau(k)
235  iq1= 4*inod - 3
236  iq2= 4*inod - 2
237  iq3= 4*inod - 1
238  iq4= 4*inod
239  au(16*k-15)= au(16*k-15)/(scale(ip1)*scale(iq1))
240  au(16*k-14)= au(16*k-14)/(scale(ip1)*scale(iq2))
241  au(16*k-13)= au(16*k-13)/(scale(ip1)*scale(iq3))
242  au(16*k-12)= au(16*k-12)/(scale(ip1)*scale(iq4))
243  au(16*k-11)= au(16*k-11)/(scale(ip2)*scale(iq1))
244  au(16*k-10)= au(16*k-10)/(scale(ip2)*scale(iq2))
245  au(16*k- 9)= au(16*k- 9)/(scale(ip2)*scale(iq3))
246  au(16*k- 8)= au(16*k- 8)/(scale(ip2)*scale(iq4))
247  au(16*k- 7)= au(16*k- 7)/(scale(ip3)*scale(iq1))
248  au(16*k- 6)= au(16*k- 6)/(scale(ip3)*scale(iq2))
249  au(16*k- 5)= au(16*k- 5)/(scale(ip3)*scale(iq3))
250  au(16*k- 4)= au(16*k- 4)/(scale(ip3)*scale(iq4))
251  au(16*k- 3)= au(16*k- 3)/(scale(ip4)*scale(iq1))
252  au(16*k- 2)= au(16*k- 2)/(scale(ip4)*scale(iq2))
253  au(16*k- 1)= au(16*k- 1)/(scale(ip4)*scale(iq3))
254  au(16*k- 0)= au(16*k- 0)/(scale(ip4)*scale(iq4))
255  enddo
256  enddo
257 
258  deallocate(scale)
259  end subroutine hecmw_solver_scaling_bk_44
260 
261 end module hecmw_solver_scaling_44
hecmw_solver_scaling_44::hecmw_solver_scaling_bk_44
subroutine, public hecmw_solver_scaling_bk_44(hecMAT)
Definition: hecmw_solver_scaling_44.f90:146
hecmw_util::hecmw_wtime
real(kind=kreal) function hecmw_wtime()
Definition: hecmw_util_f.F90:549
hecmw_util
I/O and Utility.
Definition: hecmw_util_f.F90:7
hecmw_util::hecmwst_local_mesh
Definition: hecmw_util_f.F90:234
m_hecmw_comm_f
Definition: hecmw_comm_f.F90:6
hecmw_util::kreal
integer(kind=4), parameter kreal
Definition: hecmw_util_f.F90:16
hecmw_solver_scaling_44::hecmw_solver_scaling_fw_44
subroutine, public hecmw_solver_scaling_fw_44(hecMESH, hecMAT, COMMtime)
Definition: hecmw_solver_scaling_44.f90:21
hecmw_solver_scaling_44
Definition: hecmw_solver_scaling_44.f90:6
hecmw_matrix_misc
Definition: hecmw_matrix_misc.f90:6
hecmw_matrix_misc::hecmw_mat_get_scaling
integer(kind=kint) function, public hecmw_mat_get_scaling(hecMAT)
Definition: hecmw_matrix_misc.f90:394
m_hecmw_comm_f::hecmw_update_r
subroutine hecmw_update_r(hecMESH, val, n, m)
Definition: hecmw_comm_f.F90:683
hecmw_util::hecmwst_matrix
Definition: hecmw_util_f.F90:444