FrontISTR  6.0.0
Large-scale structural analysis program with finit element method
hecmw_precond_SSOR_66.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 
6 !C
7 !C***
8 !C*** module hecmw_precond_SSOR_66
9 !C***
10 !C
11 module hecmw_precond_ssor_66
12  use hecmw_util
13  use hecmw_matrix_misc
14  use m_hecmw_matrix_ordering_cm
15  use m_hecmw_matrix_ordering_mc
16  use hecmw_matrix_reorder
17 #ifndef _OPENACC
18  !$ use omp_lib
19 #endif
20 
21  private
22 
23  public:: hecmw_precond_ssor_66_setup
24  public:: hecmw_precond_ssor_66_apply
25  public:: hecmw_precond_ssor_66_clear
26 
27  integer(kind=kint) :: N
28  real(kind=kreal), pointer :: d(:) => null()
29  real(kind=kreal), pointer :: al(:) => null()
30  real(kind=kreal), pointer :: au(:) => null()
31  integer(kind=kint), pointer :: indexL(:) => null()
32  integer(kind=kint), pointer :: indexU(:) => null()
33  integer(kind=kint), pointer :: itemL(:) => null()
34  integer(kind=kint), pointer :: itemU(:) => null()
35  real(kind=kreal), pointer :: alu(:) => null()
36 
37  integer(kind=kint) :: NColor
38  integer(kind=kint), pointer :: COLORindex(:) => null()
39  integer(kind=kint), pointer :: perm(:) => null()
40  integer(kind=kint), pointer :: iperm(:) => null()
41 
42  logical, save :: isFirst = .true.
43 
44  logical, save :: INITIALIZED = .false.
45 
46  ! relaxation parameter of
47  ! M = 1/(2-OMEGA) (D/OMEGA + L) (D/OMEGA)^-1 (D/OMEGA + U)
48  ! OMEGA = 1 is the classic (D+L) D^-1 (D+U). _apply is handed only ZP, so setup
49  ! leaves the value here for it.
50  real(kind=kreal), save :: omega = 1.d0
51 
52  ! for tuning
53  integer(kind=kint), parameter :: numOfBlockPerThread = 100
54  integer(kind=kint), save :: numOfThread = 1, numofblock
55  integer(kind=kint), save, allocatable :: icToBlockIndex(:)
56  integer(kind=kint), save, allocatable :: blockIndexToColorIndex(:)
57  integer(kind=kint), save :: sectorCacheSize0, sectorCacheSize1
58 
59  integer(kind=kint), parameter :: SSOR_DEBUG = 0
60 
61 contains
62 
63  subroutine hecmw_precond_ssor_66_setup(hecMAT)
64  implicit none
65  type(hecmwst_matrix), intent(inout) :: hecmat
66  integer(kind=kint ) :: npl, npu
67  integer(kind=kint ) :: ncolor_in
68  real (kind=kreal) :: omega_inv
69  real (kind=kreal) :: alutmp(6,6), pw(6)
70  integer(kind=kint ) :: ii, i, j, k
71  integer(kind=kint ) :: nthreads = 1
72  integer(kind=kint ), allocatable :: perm_tmp(:)
73  real (kind=kreal) :: t0
74 
75  if (ssor_debug >= 1) then
76  t0 = hecmw_wtime()
77  write(*,*) 'DEBUG: SSOR setup start', hecmw_wtime()-t0
78  endif
79 
80  if (initialized) then
81  if (hecmat%Iarray(98) == 1) then ! need symbolic and numerical setup
82  call hecmw_precond_ssor_66_clear(hecmat)
83  else if (hecmat%Iarray(97) == 1) then ! need numerical setup only
84  call hecmw_precond_ssor_66_clear(hecmat) ! TEMPORARY
85  else
86  return
87  endif
88  endif
89 
90 #ifndef _OPENACC
91  !$ nthreads = omp_get_max_threads()
92 #endif
93 
94  n = hecmat%N
95  ! N = hecMAT%NP
96  ncolor_in = hecmw_mat_get_ncolor_in(hecmat)
97  omega = hecmw_mat_get_omega(hecmat)
98  omega_inv = 1.d0 / omega
99 
100 #ifdef _OPENACC
101  allocate(colorindex(0:n), perm_tmp(n), perm(n), iperm(n))
102  call hecmw_matrix_ordering_rcm(n, hecmat%indexL, hecmat%itemL, &
103  hecmat%indexU, hecmat%itemU, perm_tmp, iperm)
104  if (ssor_debug >= 1) write(*,*) 'DEBUG: RCM ordering done', hecmw_wtime()-t0
105  call hecmw_matrix_ordering_mc(n, hecmat%indexL, hecmat%itemL, &
106  hecmat%indexU, hecmat%itemU, perm_tmp, &
107  ncolor_in, ncolor, colorindex, perm, iperm)
108  if (ssor_debug >= 1) write(*,*) 'DEBUG: MC ordering done', hecmw_wtime()-t0
109  deallocate(perm_tmp)
110 
111 #else
112  if (nthreads == 1) then
113  ncolor = 1
114  allocate(colorindex(0:1), perm(n), iperm(n))
115  colorindex(0) = 0
116  colorindex(1) = n
117  do i=1,n
118  perm(i) = i
119  iperm(i) = i
120  end do
121  else
122  allocate(colorindex(0:n), perm_tmp(n), perm(n), iperm(n))
123  call hecmw_matrix_ordering_rcm(n, hecmat%indexL, hecmat%itemL, &
124  hecmat%indexU, hecmat%itemU, perm_tmp, iperm)
125  if (ssor_debug >= 1) write(*,*) 'DEBUG: RCM ordering done', hecmw_wtime()-t0
126  call hecmw_matrix_ordering_mc(n, hecmat%indexL, hecmat%itemL, &
127  hecmat%indexU, hecmat%itemU, perm_tmp, &
128  ncolor_in, ncolor, colorindex, perm, iperm)
129  if (ssor_debug >= 1) write(*,*) 'DEBUG: MC ordering done', hecmw_wtime()-t0
130  deallocate(perm_tmp)
131 
132  endif
133 #endif
134 
135  npl = 0
136  do i=1,n
137  do j=hecmat%indexU(i-1)+1,hecmat%indexU(i)
138  if( hecmat%itemU(j) > n ) exit
139  npl = npl + 1
140  enddo
141  enddo
142  npl = max(hecmat%indexL(n),npl)
143  npu = hecmat%indexU(n)
144  allocate(indexl(0:n), indexu(0:n), iteml(npl), itemu(npu))
145  call hecmw_matrix_reorder_profile(n, perm, iperm, &
146  hecmat%indexL, hecmat%indexU, hecmat%itemL, hecmat%itemU, &
147  indexl, indexu, iteml, itemu)
148  if (ssor_debug >= 1) write(*,*) 'DEBUG: reordering profile done', hecmw_wtime()-t0
149 
150 
151  allocate(d(36*n), al(36*npl), au(36*npu))
152  call hecmw_matrix_reorder_values(n, 6, perm, iperm, &
153  hecmat%indexL, hecmat%indexU, hecmat%itemL, hecmat%itemU, &
154  hecmat%AL, hecmat%AU, hecmat%D, &
155  indexl, indexu, iteml, itemu, al, au, d)
156  if (ssor_debug >= 1) write(*,*) 'DEBUG: reordering values done', hecmw_wtime()-t0
157 
158  call hecmw_matrix_reorder_renum_item(n, perm, indexl, iteml)
159  call hecmw_matrix_reorder_renum_item(n, perm, indexu, itemu)
160 
161  allocate(alu(36*n))
162  alu = 0.d0
163 
164  do ii= 1, 36*n
165  alu(ii) = d(ii)
166  enddo
167 
168 #ifdef _OPENACC
169  !$acc kernels
170  !$acc loop independent private(ALUtmp,PW)
171 #else
172  !$omp parallel default(none),private(ii,ALUtmp,k,i,j,PW),shared(N,ALU,OMEGA_INV)
173  !$omp do
174 #endif
175  do ii= 1, n
176  alutmp(1,1)= alu(36*ii-35) * omega_inv
177  alutmp(1,2)= alu(36*ii-34) * omega_inv
178  alutmp(1,3)= alu(36*ii-33) * omega_inv
179  alutmp(1,4)= alu(36*ii-32) * omega_inv
180  alutmp(1,5)= alu(36*ii-31) * omega_inv
181  alutmp(1,6)= alu(36*ii-30) * omega_inv
182 
183  alutmp(2,1)= alu(36*ii-29) * omega_inv
184  alutmp(2,2)= alu(36*ii-28) * omega_inv
185  alutmp(2,3)= alu(36*ii-27) * omega_inv
186  alutmp(2,4)= alu(36*ii-26) * omega_inv
187  alutmp(2,5)= alu(36*ii-25) * omega_inv
188  alutmp(2,6)= alu(36*ii-24) * omega_inv
189 
190  alutmp(3,1)= alu(36*ii-23) * omega_inv
191  alutmp(3,2)= alu(36*ii-22) * omega_inv
192  alutmp(3,3)= alu(36*ii-21) * omega_inv
193  alutmp(3,4)= alu(36*ii-20) * omega_inv
194  alutmp(3,5)= alu(36*ii-19) * omega_inv
195  alutmp(3,6)= alu(36*ii-18) * omega_inv
196 
197  alutmp(4,1)= alu(36*ii-17) * omega_inv
198  alutmp(4,2)= alu(36*ii-16) * omega_inv
199  alutmp(4,3)= alu(36*ii-15) * omega_inv
200  alutmp(4,4)= alu(36*ii-14) * omega_inv
201  alutmp(4,5)= alu(36*ii-13) * omega_inv
202  alutmp(4,6)= alu(36*ii-12) * omega_inv
203 
204  alutmp(5,1)= alu(36*ii-11) * omega_inv
205  alutmp(5,2)= alu(36*ii-10) * omega_inv
206  alutmp(5,3)= alu(36*ii-9 ) * omega_inv
207  alutmp(5,4)= alu(36*ii-8 ) * omega_inv
208  alutmp(5,5)= alu(36*ii-7 ) * omega_inv
209  alutmp(5,6)= alu(36*ii-6 ) * omega_inv
210 
211  alutmp(6,1)= alu(36*ii-5 ) * omega_inv
212  alutmp(6,2)= alu(36*ii-4 ) * omega_inv
213  alutmp(6,3)= alu(36*ii-3 ) * omega_inv
214  alutmp(6,4)= alu(36*ii-2 ) * omega_inv
215  alutmp(6,5)= alu(36*ii-1 ) * omega_inv
216  alutmp(6,6)= alu(36*ii ) * omega_inv
217 
218  do k= 1, 6
219  alutmp(k,k)= 1.d0/alutmp(k,k)
220  do i= k+1, 6
221  alutmp(i,k)= alutmp(i,k) * alutmp(k,k)
222  do j= k+1, 6
223  pw(j)= alutmp(i,j) - alutmp(i,k)*alutmp(k,j)
224  enddo
225  do j= k+1, 6
226  alutmp(i,j)= pw(j)
227  enddo
228  enddo
229  enddo
230 
231  alu(36*ii-35)= alutmp(1,1)
232  alu(36*ii-34)= alutmp(1,2)
233  alu(36*ii-33)= alutmp(1,3)
234  alu(36*ii-32)= alutmp(1,4)
235  alu(36*ii-31)= alutmp(1,5)
236  alu(36*ii-30)= alutmp(1,6)
237  alu(36*ii-29)= alutmp(2,1)
238  alu(36*ii-28)= alutmp(2,2)
239  alu(36*ii-27)= alutmp(2,3)
240  alu(36*ii-26)= alutmp(2,4)
241  alu(36*ii-25)= alutmp(2,5)
242  alu(36*ii-24)= alutmp(2,6)
243  alu(36*ii-23)= alutmp(3,1)
244  alu(36*ii-22)= alutmp(3,2)
245  alu(36*ii-21)= alutmp(3,3)
246  alu(36*ii-20)= alutmp(3,4)
247  alu(36*ii-19)= alutmp(3,5)
248  alu(36*ii-18)= alutmp(3,6)
249  alu(36*ii-17)= alutmp(4,1)
250  alu(36*ii-16)= alutmp(4,2)
251  alu(36*ii-15)= alutmp(4,3)
252  alu(36*ii-14)= alutmp(4,4)
253  alu(36*ii-13)= alutmp(4,5)
254  alu(36*ii-12)= alutmp(4,6)
255  alu(36*ii-11)= alutmp(5,1)
256  alu(36*ii-10)= alutmp(5,2)
257  alu(36*ii-9 )= alutmp(5,3)
258  alu(36*ii-8 )= alutmp(5,4)
259  alu(36*ii-7 )= alutmp(5,5)
260  alu(36*ii-6 )= alutmp(5,6)
261  alu(36*ii-5 )= alutmp(6,1)
262  alu(36*ii-4 )= alutmp(6,2)
263  alu(36*ii-3 )= alutmp(6,3)
264  alu(36*ii-2 )= alutmp(6,4)
265  alu(36*ii-1 )= alutmp(6,5)
266  alu(36*ii )= alutmp(6,6)
267  enddo
268 #ifdef _OPENACC
269  !$acc end kernels
270 #else
271  !$omp end do
272  !$omp end parallel
273 #endif
274 
275  isfirst = .true.
276 
277  initialized = .true.
278  hecmat%Iarray(98) = 0 ! symbolic setup done
279  hecmat%Iarray(97) = 0 ! numerical setup done
280 
281  if (ssor_debug >= 1) write(*,*) 'DEBUG: SSOR setup done', hecmw_wtime()-t0
282 
283  end subroutine hecmw_precond_ssor_66_setup
284 
285  subroutine setup_tuning_parameters
286  use hecmw_tuning_fx
287  implicit none
288  integer(kind=kint) :: blockindex, elementcount, numofelement, ii
289  real(kind=kreal) :: numofelementperblock
290  integer(kind=kint) :: ic, i
291  if (ssor_debug >= 1) write(*,*) 'DEBUG: setting up tuning parameters for SSOR'
292 #ifndef _OPENACC
293  !$ numOfThread = omp_get_max_threads()
294 #endif
295 
296  numofblock = numofthread * numofblockperthread
297  if (allocated(ictoblockindex)) deallocate(ictoblockindex)
298  if (allocated(blockindextocolorindex)) deallocate(blockindextocolorindex)
299  allocate (ictoblockindex(0:ncolor), &
300  blockindextocolorindex(0:numofblock + ncolor))
301  numofelement = n + indexl(n) + indexu(n)
302  numofelementperblock = dble(numofelement) / numofblock
303  blockindex = 0
304  ictoblockindex = -1
305  ictoblockindex(0) = 0
306  blockindextocolorindex = -1
307  blockindextocolorindex(0) = 0
308  do ic = 1, ncolor
309  elementcount = 0
310  ii = 1
311  do i = colorindex(ic-1)+1, colorindex(ic)
312  elementcount = elementcount + 1
313  elementcount = elementcount + (indexl(i) - indexl(i-1))
314  elementcount = elementcount + (indexu(i) - indexu(i-1))
315  if (elementcount > ii * numofelementperblock &
316  .or. i == colorindex(ic)) then
317  ii = ii + 1
318  blockindex = blockindex + 1
319  blockindextocolorindex(blockindex) = i
320  endif
321  enddo
322  ictoblockindex(ic) = blockindex
323  enddo
324  numofblock = blockindex
325 
326  call hecmw_tuning_fx_calc_sector_cache( n, 6, &
327  sectorcachesize0, sectorcachesize1 )
328  end subroutine setup_tuning_parameters
329 
330  subroutine hecmw_precond_ssor_66_apply(ZP)
331  implicit none
332  real(kind=kreal), intent(inout) :: zp(:)
333  integer(kind=kint) :: ic, i, iold, j, isl, iel, isu, ieu, k
334  real(kind=kreal) :: x1, x2, x3, x4, x5, x6
335  real(kind=kreal) :: sw1, sw2, sw3, sw4, sw5, sw6
336 
337  real(kind=kreal) :: omega_fac
338 
339  ! added for tuning >>>
340  integer(kind=kint) :: blockindex
341 
342 #ifndef _OPENACC
343  if (isfirst) then
344  call setup_tuning_parameters
345  isfirst = .false.
346  endif
347 #endif
348  ! <<< added for tuning
349 
350  omega_fac = 2.d0 - omega
351 
352 #ifndef _OPENACC
353  !call start_collection("loopInPrecond66")
354 
355  !OCL CACHE_SECTOR_SIZE(sectorCacheSize0,sectorCacheSize1)
356  !OCL CACHE_SUBSECTOR_ASSIGN(ZP)
357 
358  !$omp parallel default(none) &
359  !$omp&shared(NColor,indexL,itemL,indexU,itemU,AL,AU,D,ALU,perm,&
360  !$omp& ZP,icToBlockIndex,blockIndexToColorIndex,OMEGA_FAC) &
361  !$omp&private(SW1,SW2,SW3,SW4,SW5,SW6,X1,X2,X3,X4,X5,X6,ic,i,iold,isL,ieL,isU,ieU,j,k,blockIndex)
362 #endif
363 
364  !C-- FORWARD
365  do ic=1,ncolor
366 #ifdef _OPENACC
367  !$acc kernels
368  !$acc loop independent
369  do i = colorindex(ic-1)+1, colorindex(ic)
370 #else
371  !$omp do schedule (static, 1)
372  do blockindex = ictoblockindex(ic-1)+1, ictoblockindex(ic)
373  do i = blockindextocolorindex(blockindex-1)+1, &
374  blockindextocolorindex(blockindex)
375 #endif
376  ! do i = startPos(threadNum, ic), endPos(threadNum, ic)
377  iold = perm(i)
378  sw1= omega_fac * zp(6*iold-5)
379  sw2= omega_fac * zp(6*iold-4)
380  sw3= omega_fac * zp(6*iold-3)
381  sw4= omega_fac * zp(6*iold-2)
382  sw5= omega_fac * zp(6*iold-1)
383  sw6= omega_fac * zp(6*iold )
384  isl= indexl(i-1)+1
385  iel= indexl(i)
386  do j= isl, iel
387  !k= perm(itemL(j))
388  k= iteml(j)
389  x1= zp(6*k-5)
390  x2= zp(6*k-4)
391  x3= zp(6*k-3)
392  x4= zp(6*k-2)
393  x5= zp(6*k-1)
394  x6= zp(6*k )
395  sw1= sw1 -al(36*j-35)*x1 -al(36*j-34)*x2 -al(36*j-33)*x3 -al(36*j-32)*x4 -al(36*j-31)*x5 -al(36*j-30)*x6
396  sw2= sw2 -al(36*j-29)*x1 -al(36*j-28)*x2 -al(36*j-27)*x3 -al(36*j-26)*x4 -al(36*j-25)*x5 -al(36*j-24)*x6
397  sw3= sw3 -al(36*j-23)*x1 -al(36*j-22)*x2 -al(36*j-21)*x3 -al(36*j-20)*x4 -al(36*j-19)*x5 -al(36*j-18)*x6
398  sw4= sw4 -al(36*j-17)*x1 -al(36*j-16)*x2 -al(36*j-15)*x3 -al(36*j-14)*x4 -al(36*j-13)*x5 -al(36*j-12)*x6
399  sw5= sw5 -al(36*j-11)*x1 -al(36*j-10)*x2 -al(36*j-9 )*x3 -al(36*j-8 )*x4 -al(36*j-7 )*x5 -al(36*j-6 )*x6
400  sw6= sw6 -al(36*j-5 )*x1 -al(36*j-4 )*x2 -al(36*j-3 )*x3 -al(36*j-2 )*x4 -al(36*j-1 )*x5 -al(36*j )*x6
401  enddo ! j
402 
403  x1= sw1
404  x2= sw2
405  x3= sw3
406  x4= sw4
407  x5= sw5
408  x6= sw6
409  x2= x2 -alu(36*i-29)*x1
410  x3= x3 -alu(36*i-23)*x1 -alu(36*i-22)*x2
411  x4= x4 -alu(36*i-17)*x1 -alu(36*i-16)*x2 -alu(36*i-5)*x3
412  x5= x5 -alu(36*i-11)*x1 -alu(36*i-10)*x2 -alu(36*i-9)*x3 -alu(36*i-8)*x4
413  x6= x6 -alu(36*i-5 )*x1 -alu(36*i-4 )*x2 -alu(36*i-3)*x3 -alu(36*i-2)*x4 -alu(36*i-1)*x5
414  x6= alu(36*i )* x6
415  x5= alu(36*i-7 )*( x5 -alu(36*i-6)*x6 )
416  x4= alu(36*i-14)*( x4 -alu(36*i-12)*x6 -alu(36*i-13)*x5)
417  x3= alu(36*i-21)*( x3 -alu(36*i-18)*x6 -alu(36*i-19)*x5 -alu(36*i-20)*x4)
418  x2= alu(36*i-28)*( x2 -alu(36*i-24)*x6 -alu(36*i-25)*x5 -alu(36*i-26)*x4 -alu(36*i-27)*x3)
419  x1= alu(36*i-35)*( x1 -alu(36*i-30)*x6 -alu(36*i-31)*x5 -alu(36*i-32)*x4 -alu(36*i-33)*x3 -alu(36*i-34)*x2)
420  zp(6*iold-5)= x1
421  zp(6*iold-4)= x2
422  zp(6*iold-3)= x3
423  zp(6*iold-2)= x4
424  zp(6*iold-1)= x5
425  zp(6*iold )= x6
426 #ifdef _OPENACC
427  enddo
428  !$acc end kernels
429 #else
430  enddo ! i
431  enddo ! blockIndex
432  !$omp end do
433 #endif
434  enddo ! ic
435 
436  !C-- BACKWARD
437  do ic=ncolor, 1, -1
438 #ifdef _OPENACC
439  !$acc kernels
440  !$acc loop independent
441  do i = colorindex(ic-1)+1, colorindex(ic)
442 #else
443  !$omp do schedule (static, 1)
444  do blockindex = ictoblockindex(ic), ictoblockindex(ic-1)+1, -1
445  do i = blockindextocolorindex(blockindex), &
446  blockindextocolorindex(blockindex-1)+1, -1
447 #endif
448  ! do blockIndex = icToBlockIndex(ic-1)+1, icToBlockIndex(ic)
449  ! do i = blockIndexToColorIndex(blockIndex-1)+1, &
450  ! blockIndexToColorIndex(blockIndex)
451  ! do i = endPos(threadNum, ic), startPos(threadNum, ic), -1
452  sw1= 0.d0
453  sw2= 0.d0
454  sw3= 0.d0
455  sw4= 0.d0
456  sw5= 0.d0
457  sw6= 0.d0
458  isu= indexu(i-1) + 1
459  ieu= indexu(i)
460  do j= ieu, isu, -1
461  !k= perm(itemU(j))
462  k= itemu(j)
463  x1= zp(6*k-5)
464  x2= zp(6*k-4)
465  x3= zp(6*k-3)
466  x4= zp(6*k-2)
467  x5= zp(6*k-1)
468  x6= zp(6*k )
469  sw1= sw1 +au(36*j-35)*x1 +au(36*j-34)*x2 +au(36*j-33)*x3 +au(36*j-32)*x4 +au(36*j-31)*x5 +au(36*j-30)*x6
470  sw2= sw2 +au(36*j-29)*x1 +au(36*j-28)*x2 +au(36*j-27)*x3 +au(36*j-26)*x4 +au(36*j-25)*x5 +au(36*j-24)*x6
471  sw3= sw3 +au(36*j-23)*x1 +au(36*j-22)*x2 +au(36*j-21)*x3 +au(36*j-20)*x4 +au(36*j-19)*x5 +au(36*j-18)*x6
472  sw4= sw4 +au(36*j-17)*x1 +au(36*j-16)*x2 +au(36*j-15)*x3 +au(36*j-14)*x4 +au(36*j-13)*x5 +au(36*j-12)*x6
473  sw5= sw5 +au(36*j-11)*x1 +au(36*j-10)*x2 +au(36*j-9 )*x3 +au(36*j-8 )*x4 +au(36*j-7 )*x5 +au(36*j-6 )*x6
474  sw6= sw6 +au(36*j-5 )*x1 +au(36*j-4 )*x2 +au(36*j-3 )*x3 +au(36*j-2 )*x4 +au(36*j-1 )*x5 +au(36*j )*x6
475  enddo ! j
476 
477  x1= sw1
478  x2= sw2
479  x3= sw3
480  x4= sw4
481  x5= sw5
482  x6= sw6
483  x2= x2 -alu(36*i-29)*x1
484  x3= x3 -alu(36*i-23)*x1 -alu(36*i-22)*x2
485  x4= x4 -alu(36*i-17)*x1 -alu(36*i-16)*x2 -alu(36*i-5)*x3
486  x5= x5 -alu(36*i-11)*x1 -alu(36*i-10)*x2 -alu(36*i-9)*x3 -alu(36*i-8)*x4
487  x6= x6 -alu(36*i-5 )*x1 -alu(36*i-4 )*x2 -alu(36*i-3)*x3 -alu(36*i-2)*x4 -alu(36*i-1)*x5
488  x6= alu(36*i )* x6
489  x5= alu(36*i-7 )*( x5 -alu(36*i-6)*x6 )
490  x4= alu(36*i-14)*( x4 -alu(36*i-12)*x6 -alu(36*i-13)*x5)
491  x3= alu(36*i-21)*( x3 -alu(36*i-18)*x6 -alu(36*i-19)*x5 -alu(36*i-20)*x4)
492  x2= alu(36*i-28)*( x2 -alu(36*i-24)*x6 -alu(36*i-25)*x5 -alu(36*i-26)*x4 -alu(36*i-27)*x3)
493  x1= alu(36*i-35)*( x1 -alu(36*i-30)*x6 -alu(36*i-31)*x5 -alu(36*i-32)*x4 -alu(36*i-33)*x3 -alu(36*i-34)*x2)
494  iold = perm(i)
495  zp(6*iold-5)= zp(6*iold-5) -x1
496  zp(6*iold-4)= zp(6*iold-4) -x2
497  zp(6*iold-3)= zp(6*iold-3) -x3
498  zp(6*iold-2)= zp(6*iold-2) -x4
499  zp(6*iold-1)= zp(6*iold-1) -x5
500  zp(6*iold )= zp(6*iold ) -x6
501 #ifdef _OPENACC
502  enddo
503  !$acc end kernels
504 #else
505  enddo ! i
506  enddo ! blockIndex
507  !$omp end do
508 #endif
509  enddo ! ic
510 #ifndef _OPENACC
511  !$omp end parallel
512 
513  !OCL END_CACHE_SUBSECTOR
514  !OCL END_CACHE_SECTOR_SIZE
515 
516  !call stop_collection("loopInPrecond66")
517 #endif
518 
519  end subroutine hecmw_precond_ssor_66_apply
520 
521  subroutine hecmw_precond_ssor_66_clear(hecMAT)
522  implicit none
523  type(hecmwst_matrix), intent(inout) :: hecmat
524  if (associated(colorindex)) deallocate(colorindex)
525  if (associated(perm)) deallocate(perm)
526  if (associated(iperm)) deallocate(iperm)
527  if (associated(alu)) deallocate(alu)
528  if (associated(d)) deallocate(d)
529  if (associated(al)) deallocate(al)
530  if (associated(au)) deallocate(au)
531  if (associated(indexl)) deallocate(indexl)
532  if (associated(indexu)) deallocate(indexu)
533  if (associated(iteml)) deallocate(iteml)
534  if (associated(itemu)) deallocate(itemu)
535  nullify(colorindex)
536  nullify(perm)
537  nullify(iperm)
538  nullify(alu)
539  nullify(d)
540  nullify(al)
541  nullify(au)
542  nullify(indexl)
543  nullify(indexu)
544  nullify(iteml)
545  nullify(itemu)
546  initialized = .false.
547  end subroutine hecmw_precond_ssor_66_clear
548 
549 end module hecmw_precond_ssor_66
hecmw_matrix_misc
Definition: hecmw_matrix_misc.F90:6
hecmw_matrix_misc::hecmw_mat_get_omega
real(kind=kreal) function, public hecmw_mat_get_omega(hecMAT)
Definition: hecmw_matrix_misc.F90:832
hecmw_matrix_misc::hecmw_mat_get_ncolor_in
integer(kind=kint) function, public hecmw_mat_get_ncolor_in(hecMAT)
Definition: hecmw_matrix_misc.F90:644
hecmw_matrix_reorder
Definition: hecmw_matrix_reorder.f90:6
hecmw_matrix_reorder::hecmw_matrix_reorder_values
subroutine, public hecmw_matrix_reorder_values(N, NDOF, perm, iperm, indexL, indexU, itemL, itemU, AL, AU, D, indexLp, indexUp, itemLp, itemUp, ALp, AUp, Dp)
Definition: hecmw_matrix_reorder.f90:69
hecmw_matrix_reorder::hecmw_matrix_reorder_renum_item
subroutine, public hecmw_matrix_reorder_renum_item(N, perm, indexXp, itemXp)
Definition: hecmw_matrix_reorder.f90:91
hecmw_matrix_reorder::hecmw_matrix_reorder_profile
subroutine, public hecmw_matrix_reorder_profile(N, perm, iperm, indexL, indexU, itemL, itemU, indexLp, indexUp, itemLp, itemUp)
Definition: hecmw_matrix_reorder.f90:20
hecmw_precond_ssor_66
Definition: hecmw_precond_SSOR_66.F90:11
hecmw_precond_ssor_66::hecmw_precond_ssor_66_setup
subroutine, public hecmw_precond_ssor_66_setup(hecMAT)
Definition: hecmw_precond_SSOR_66.F90:64
hecmw_precond_ssor_66::hecmw_precond_ssor_66_clear
subroutine, public hecmw_precond_ssor_66_clear(hecMAT)
Definition: hecmw_precond_SSOR_66.F90:522
hecmw_precond_ssor_66::hecmw_precond_ssor_66_apply
subroutine, public hecmw_precond_ssor_66_apply(ZP)
Definition: hecmw_precond_SSOR_66.F90:331
hecmw_tuning_fx
Definition: hecmw_tuning_fx.f90:6
hecmw_tuning_fx::hecmw_tuning_fx_calc_sector_cache
subroutine, public hecmw_tuning_fx_calc_sector_cache(N, NDOF, sectorCacheSize0, sectorCacheSize1)
Definition: hecmw_tuning_fx.f90:25
hecmw_util
I/O and Utility.
Definition: hecmw_util_f.F90:7
hecmw_util::kreal
integer(kind=4), parameter kreal
Definition: hecmw_util_f.F90:16
hecmw_util::hecmw_wtime
real(kind=kreal) function hecmw_wtime()
Definition: hecmw_util_f.F90:635
m_hecmw_matrix_ordering_cm
Definition: hecmw_matrix_ordering_CM.f90:6
m_hecmw_matrix_ordering_cm::hecmw_matrix_ordering_rcm
subroutine, public hecmw_matrix_ordering_rcm(N, indexL, itemL, indexU, itemU, perm, iperm)
Definition: hecmw_matrix_ordering_CM.f90:67
m_hecmw_matrix_ordering_mc
Definition: hecmw_matrix_ordering_MC.f90:6
m_hecmw_matrix_ordering_mc::hecmw_matrix_ordering_mc
subroutine, public hecmw_matrix_ordering_mc(N, indexL, itemL, indexU, itemU, perm_cur, ncolor_in, ncolor_out, COLORindex, perm, iperm)
Definition: hecmw_matrix_ordering_MC.f90:19
hecmw_util::hecmwst_matrix
Definition: hecmw_util_f.F90:465