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