FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_solver_las_nn.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
8 implicit none
9
10 private
11
12 public :: hecmw_matvec_nn
15 public :: hecmw_matresid_nn
16 public :: hecmw_rel_resid_l2_nn
17 public :: hecmw_tvec_nn
18 public :: hecmw_ttvec_nn
19 public :: hecmw_ttmattvec_nn
20 public :: hecmw_mat_diag_sr_nn
21 public :: hecmw_mat_add_nn
22 public :: hecmw_mat_multiple_nn
23
24 ! ! for communication hiding in matvec
25 ! integer(kind=kint), save, allocatable :: index_o(:), item_o(:)
26 ! real(kind=kreal), save, allocatable :: A_o(:)
27 logical, save :: async_matvec_flg = .false.
28
29contains
30
31 !C
32 !C***
33 !C*** hecmw_matvec_nn
34 !C***
35 !C
36 subroutine hecmw_matvec_nn (hecMESH, hecMAT, X, Y, time_Ax, COMMtime)
37 use hecmw_util
39 implicit none
40 type (hecmwst_local_mesh), intent(in) :: hecmesh
41 type (hecmwst_matrix), intent(in), target :: hecmat
42 real(kind=kreal), intent(in) :: x(:)
43 real(kind=kreal), intent(out) :: y(:)
44 real(kind=kreal), intent(inout) :: time_ax
45 real(kind=kreal), intent(inout), optional :: commtime
46
47 real(kind=kreal) :: tcomm
48 real(kind=kreal), allocatable :: wk(:)
49
50 tcomm = 0.d0
51
52 if (hecmw_mat_get_flag_mpcmatvec(hecmat) /= 0) then
53 allocate(wk(hecmat%NP * hecmat%NDOF))
54 call hecmw_ttmattvec_nn(hecmesh, hecmat, x, y, wk, tcomm)
55 deallocate(wk)
56 else
57 call hecmw_matvec_nn_inner(hecmesh, hecmat, x, y, tcomm)
58 endif
59
60 if (present(commtime)) commtime = commtime + tcomm
61 end subroutine hecmw_matvec_nn
62
63 !C
64 !C***
65 !C*** hecmw_matvec_nn_set_async
66 !C***
67 !C
68 subroutine hecmw_matvec_nn_set_async (hecMAT)
69 use hecmw_util
70 implicit none
71 type (hecmwst_matrix), intent(in) :: hecmat
72 ! integer(kind=kint) :: i, j, jS, jE, idx, in
73
74 ! allocate(index_o(0:hecMAT%N))
75 ! index_o(0) = 0
76 ! do i = 1, hecMAT%N
77 ! jS= hecMAT%indexU(i-1) + 1
78 ! jE= hecMAT%indexU(i )
79 ! idx = index_o(i-1)
80 ! do j= jS, jE
81 ! in = hecMAT%itemU(j)
82 ! if (in <= hecMAT%N) cycle
83 ! idx = idx + 1
84 ! enddo
85 ! index_o(i) = idx
86 ! enddo
87 ! allocate(item_o(idx))
88 ! allocate(A_o(idx*9))
89 ! do i = 1, hecMAT%N
90 ! jS= hecMAT%indexU(i-1) + 1
91 ! jE= hecMAT%indexU(i )
92 ! idx = index_o(i-1)
93 ! do j= jS, jE
94 ! in = hecMAT%itemU(j)
95 ! if (in <= hecMAT%N) cycle
96 ! idx = idx + 1
97 ! item_o(idx) = hecMAT%itemU(j) - hecMAT%N
98 ! A_o(9*idx-8:9*idx) = hecMAT%AU(9*j-8:9*j)
99 ! enddo
100 ! enddo
101 ! async_matvec_flg = .true.
102 end subroutine hecmw_matvec_nn_set_async
103
104 !C
105 !C***
106 !C*** hecmw_matvec_nn_unset_async
107 !C***
108 !C
110 implicit none
111 ! if (allocated(index_o)) deallocate(index_o)
112 ! if (allocated(item_o)) deallocate(item_o)
113 ! if (allocated(A_o)) deallocate(A_o)
114 ! async_matvec_flg = .false.
115 end subroutine hecmw_matvec_nn_unset_async
116
117 !C
118 !C***
119 !C*** hecmw_matvec_nn_inner ( private subroutine )
120 !C***
121 !C
122 subroutine hecmw_matvec_nn_inner (hecMESH, hecMAT, X, Y, time_Ax, COMMtime)
123 use hecmw_util
129 !$ use omp_lib
130
131 implicit none
132 type (hecmwst_local_mesh), intent(in) :: hecmesh
133 type (hecmwst_matrix), intent(in), target :: hecmat
134 real(kind=kreal), intent(in) :: x(:)
135 real(kind=kreal), intent(out) :: y(:)
136 real(kind=kreal), intent(inout) :: time_ax
137 real(kind=kreal), intent(inout), optional :: commtime
138
139 real(kind=kreal) :: start_time, end_time, tcomm
140 integer(kind=kint) :: i, j, k, l, js, je, in
141 real(kind=kreal) :: yv(hecmat%NDOF), xv(hecmat%NDOF)
142
143 integer(kind=kint) :: n, np, ndof, ndof2
144 integer(kind=kint), pointer :: indexl(:), iteml(:), indexu(:), itemu(:)
145 real(kind=kreal), pointer :: al(:), au(:), d(:)
146
147 ! added for turning >>>
148 integer, parameter :: numofblockperthread = 100
149 logical, save :: isfirst = .true.
150 integer, save :: numofthread = 1
151 integer, save, allocatable :: startpos(:), endpos(:)
152 integer(kind=kint), save :: sectorcachesize0, sectorcachesize1
153 integer(kind=kint) :: threadnum, blocknum, numofblock
154 integer(kind=kint) :: numofelement, elementcount, blockindex
155 real(kind=kreal) :: numofelementperblock
156 ! <<< added for turning
157
158 if (hecmw_jad_is_initialized().ne.0) then
159 tcomm = 0.d0
160 start_time = hecmw_wtime()
161 call hecmw_jad_matvec(hecmesh, hecmat, x, y, tcomm)
162 end_time = hecmw_wtime()
163 time_ax = time_ax + end_time - start_time - tcomm
164 if (present(commtime)) commtime = commtime + tcomm
165 else
166
167 n = hecmat%N
168 np = hecmat%NP
169 indexl => hecmat%indexL
170 indexu => hecmat%indexU
171 iteml => hecmat%itemL
172 itemu => hecmat%itemU
173 al => hecmat%AL
174 au => hecmat%AU
175 d => hecmat%D
176 ndof = hecmat%NDOF
177 ndof2 = ndof*ndof
178
179 ! added for turning >>>
180 if (.not. isfirst) then
181 numofblock = numofthread * numofblockperthread
182 if (endpos(numofblock-1) .ne. n-1) then
183 deallocate(startpos, endpos)
184 isfirst = .true.
185 endif
186 endif
187 if (isfirst) then
188 !$ numOfThread = omp_get_max_threads()
189 numofblock = numofthread * numofblockperthread
190 allocate (startpos(0 : numofblock - 1), endpos(0 : numofblock - 1))
191 numofelement = n + indexl(n) + indexu(n)
192 numofelementperblock = dble(numofelement) / numofblock
193 blocknum = 0
194 elementcount = 0
195 startpos(blocknum) = 1
196 do i= 1, n
197 elementcount = elementcount + 1
198 elementcount = elementcount + (indexl(i) - indexl(i-1))
199 elementcount = elementcount + (indexu(i) - indexu(i-1))
200 if (elementcount > (blocknum + 1) * numofelementperblock) then
201 endpos(blocknum) = i
202 ! write(9000+hecMESH%my_rank,*) mod(blockNum, numOfThread), &
203 ! startPos(blockNum), endPos(blockNum)
204 blocknum = blocknum + 1
205 startpos(blocknum) = i + 1
206 if (blocknum == (numofblock - 1)) exit
207 endif
208 enddo
209 endpos(blocknum) = n
210 ! write(9000+hecMESH%my_rank,*) mod(blockNum, numOfThread), &
211 ! startPos(blockNum), endPos(blockNum)
212 ! for irregular data
213 do i= blocknum+1, numofblock-1
214 startpos(i) = n
215 endpos(i) = n-1
216 ! write(9000+hecMESH%my_rank,*) mod(i, numOfThread), &
217 ! startPos(i), endPos(i)
218 end do
219
220 call hecmw_tuning_fx_calc_sector_cache(np, ndof, &
221 sectorcachesize0, sectorcachesize1)
222
223 isfirst = .false.
224 endif
225 ! <<< added for turning
226
227 start_time= hecmw_wtime()
228 ! if (async_matvec_flg) then
229 ! call hecmw_update_3_R_async (hecMESH, X, NP, ireq)
230 ! else
231 call hecmw_update_m_r (hecmesh, x, np, ndof)
232 ! endif
233 end_time= hecmw_wtime()
234 if (present(commtime)) commtime = commtime + end_time - start_time
235
236 start_time = hecmw_wtime()
237
238 !call fapp_start("loopInMatvec33", 1, 0)
239 !call start_collection("loopInMatvec33")
240
241 !OCL CACHE_SECTOR_SIZE(sectorCacheSize0,sectorCacheSize1)
242 !OCL CACHE_SUBSECTOR_ASSIGN(X)
243
244 !$OMP PARALLEL DEFAULT(NONE) &
245 !$OMP&PRIVATE(i,XV,YV,jS,jE,j,k,l,in,threadNum,blockNum,blockIndex) &
246 !$OMP&SHARED(D,AL,AU,indexL,itemL,indexU,itemU,X,Y,startPos,endPos,numOfThread,N,NDOF,NDOF2,async_matvec_flg)
247 threadnum = 0
248 !$ threadNum = omp_get_thread_num()
249 do blocknum = 0 , numofblockperthread - 1
250 blockindex = blocknum * numofthread + threadnum
251 do i = startpos(blockindex), endpos(blockindex)
252 do k=1,ndof
253 xv(k) = x(ndof*(i-1)+k)
254 end do
255 yv(:)=0.0d0
256 do k=1,ndof
257 do l=1,ndof
258 yv(k)=yv(k)+d(ndof2*(i-1)+(k-1)*ndof+l)*xv(l)
259 end do
260 end do
261 js= indexl(i-1) + 1
262 je= indexl(i )
263 do j= js, je
264 in = iteml(j)
265 do k=1,ndof
266 xv(k) = x(ndof*(in-1)+k)
267 end do
268 do k=1,ndof
269 do l=1,ndof
270 yv(k)=yv(k)+al(ndof2*(j-1)+(k-1)*ndof+l)*xv(l)
271 end do
272 end do
273 enddo
274 js= indexu(i-1) + 1
275 je= indexu(i )
276 do j= js, je
277 in = itemu(j)
278 ! if (async_matvec_flg .and. in > N) cycle
279 do k=1,ndof
280 xv(k) = x(ndof*(in-1)+k)
281 end do
282 do k=1,ndof
283 do l=1,ndof
284 yv(k)=yv(k)+au(ndof2*(j-1)+(k-1)*ndof+l)*xv(l)
285 end do
286 end do
287 enddo
288 do k=1,ndof
289 y(ndof*(i-1)+k) = yv(k)
290 end do
291 enddo
292 enddo
293 !$OMP END PARALLEL
294
295 !OCL END_CACHE_SUBSECTOR
296 !OCL END_CACHE_SECTOR_SIZE
297
298 !call stop_collection("loopInMatvec33")
299 !call fapp_stop("loopInMatvec33", 1, 0)
300
301 end_time = hecmw_wtime()
302 time_ax = time_ax + end_time - start_time
303
304
305
306 endif
307
308 if (hecmat%cmat%n_val > 0) then
309 call hecmw_cmat_multvec_add( hecmat%cmat, x, y, np * ndof )
310 end if
311 end subroutine hecmw_matvec_nn_inner
312
313
314
315 !C
316 !C***
317 !C*** hecmw_matresid_nn
318 !C***
319 !C
320 subroutine hecmw_matresid_nn (hecMESH, hecMAT, X, B, R, COMMtime)
321 use hecmw_util
322 implicit none
323 type (hecmwst_local_mesh), intent(in) :: hecmesh
324 type (hecmwst_matrix), intent(in) :: hecmat
325 real(kind=kreal), intent(in) :: x(:), b(:)
326 real(kind=kreal), intent(out) :: r(:)
327 real(kind=kreal), intent(inout), optional :: commtime
328
329 integer(kind=kint) :: i
330 real(kind=kreal) :: tcomm
331
332 tcomm = 0.d0
333 call hecmw_matvec_nn (hecmesh, hecmat, x, r, tcomm)
334 if (present(commtime)) commtime = commtime + tcomm
335 !$omp parallel default(none),private(i),shared(hecMAT,R,B)
336 !$omp do
337 do i = 1, hecmat%N * hecmat%NDOF
338 r(i) = b(i) - r(i)
339 enddo
340 !$omp end do
341 !$omp end parallel
342 end subroutine hecmw_matresid_nn
343
344 !C
345 !C***
346 !C*** hecmw_rel_resid_L2_nn
347 !C***
348 !C
349 function hecmw_rel_resid_l2_nn (hecMESH, hecMAT, COMMtime)
350 use hecmw_util
352 implicit none
353 real(kind=kreal) :: hecmw_rel_resid_l2_nn
354 type ( hecmwst_local_mesh ), intent(in) :: hecmesh
355 type ( hecmwst_matrix ), intent(in) :: hecmat
356 real(kind=kreal), intent(inout), optional :: commtime
357
358 real(kind=kreal), allocatable :: r(:)
359 real(kind=kreal) :: bnorm2, rnorm2
360 real(kind=kreal) :: tcomm
361
362 allocate(r(hecmat%NDOF*hecmat%NP))
363
364 tcomm = 0.d0
365 call hecmw_innerproduct_r(hecmesh, hecmat%NDOF, &
366 hecmat%B, hecmat%B, bnorm2, tcomm)
367 if (bnorm2 == 0.d0) then
368 bnorm2 = 1.d0
369 endif
370 call hecmw_matresid_nn(hecmesh, hecmat, hecmat%X, hecmat%B, r, tcomm)
371 call hecmw_innerproduct_r(hecmesh, hecmat%NDOF, r, r, rnorm2, tcomm)
372 hecmw_rel_resid_l2_nn = sqrt(rnorm2 / bnorm2)
373
374 if (present(commtime)) commtime = commtime + tcomm
375
376 deallocate(r)
377 end function hecmw_rel_resid_l2_nn
378
379 !C
380 !C***
381 !C*** hecmw_Tvec_nn
382 !C***
383 !C
384 subroutine hecmw_tvec_nn (hecMESH, ndof, X, Y, COMMtime)
385 use hecmw_util
387 implicit none
388 type (hecmwst_local_mesh), intent(in) :: hecmesh
389 integer(kind=kint), intent(in) :: ndof
390 real(kind=kreal), intent(in) :: x(:)
391 real(kind=kreal), intent(out) :: y(:)
392 real(kind=kreal), intent(inout) :: commtime
393
394 real(kind=kreal) :: start_time, end_time
395 integer(kind=kint) :: i, j, jj, k, kk
396
397 start_time= hecmw_wtime()
398 call hecmw_update_m_r (hecmesh, x, hecmesh%n_node, ndof)
399 end_time= hecmw_wtime()
400 commtime = commtime + end_time - start_time
401
402 !$omp parallel default(none),private(i,k,kk,j,jj),shared(hecMESH,X,Y),firstprivate(ndof)
403 !$omp do
404 do i= 1, hecmesh%nn_internal * ndof
405 y(i)= x(i)
406 enddo
407 !$omp end do
408
409 !$omp do
410 outer: do i= 1, hecmesh%mpc%n_mpc
411 do j= hecmesh%mpc%mpc_index(i-1) + 1, hecmesh%mpc%mpc_index(i)
412 if (hecmesh%mpc%mpc_dof(j) > ndof) cycle outer
413 enddo
414 k = hecmesh%mpc%mpc_index(i-1) + 1
415 kk = ndof * (hecmesh%mpc%mpc_item(k) - 1) + hecmesh%mpc%mpc_dof(k)
416 y(kk) = 0.d0
417 do j= hecmesh%mpc%mpc_index(i-1) + 2, hecmesh%mpc%mpc_index(i)
418 jj = ndof * (hecmesh%mpc%mpc_item(j) - 1) + hecmesh%mpc%mpc_dof(j)
419 y(kk) = y(kk) - hecmesh%mpc%mpc_val(j) * x(jj)
420 enddo
421 enddo outer
422 !$omp end do
423 !$omp end parallel
424
425 end subroutine hecmw_tvec_nn
426
427 !C
428 !C***
429 !C*** hecmw_Ttvec_nn
430 !C***
431 !C
432 subroutine hecmw_ttvec_nn (hecMESH, ndof, X, Y, COMMtime)
433 use hecmw_util
435 implicit none
436 type (hecmwst_local_mesh), intent(in) :: hecmesh
437 integer(kind=kint), intent(in) :: ndof
438 real(kind=kreal), intent(in) :: x(:)
439 real(kind=kreal), intent(out) :: y(:)
440 real(kind=kreal), intent(inout) :: commtime
441
442 real(kind=kreal) :: start_time, end_time
443 integer(kind=kint) :: i, j, jj, k, kk
444
445 start_time= hecmw_wtime()
446 call hecmw_update_m_r (hecmesh, x, hecmesh%n_node,ndof)
447 end_time= hecmw_wtime()
448 commtime = commtime + end_time - start_time
449
450 !$omp parallel default(none),private(i,k,kk,j,jj),shared(hecMESH,X,Y),firstprivate(ndof)
451 !$omp do
452 do i= 1, hecmesh%nn_internal * ndof
453 y(i)= x(i)
454 enddo
455 !$omp end do
456
457 !$omp do
458 outer: do i= 1, hecmesh%mpc%n_mpc
459 do j= hecmesh%mpc%mpc_index(i-1) + 1, hecmesh%mpc%mpc_index(i)
460 if (hecmesh%mpc%mpc_dof(j) > ndof) cycle outer
461 enddo
462 k = hecmesh%mpc%mpc_index(i-1) + 1
463 kk = ndof * (hecmesh%mpc%mpc_item(k) - 1) + hecmesh%mpc%mpc_dof(k)
464 y(kk) = 0.d0
465 do j= hecmesh%mpc%mpc_index(i-1) + 2, hecmesh%mpc%mpc_index(i)
466 jj = ndof * (hecmesh%mpc%mpc_item(j) - 1) + hecmesh%mpc%mpc_dof(j)
467 y(jj) = y(jj) - hecmesh%mpc%mpc_val(j) * x(kk)
468 enddo
469 enddo outer
470 !$omp end do
471 !$omp end parallel
472
473
474 end subroutine hecmw_ttvec_nn
475
476 !C
477 !C***
478 !C*** hecmw_TtmatTvec_nn
479 !C***
480 !C
481 subroutine hecmw_ttmattvec_nn (hecMESH, hecMAT, X, Y, W, COMMtime)
482 use hecmw_util
483 implicit none
484 type (hecmwst_local_mesh), intent(in) :: hecmesh
485 type (hecmwst_matrix), intent(in) :: hecmat
486 real(kind=kreal), intent(in) :: x(:)
487 real(kind=kreal), intent(out) :: y(:), w(:)
488 real(kind=kreal), intent(inout) :: commtime
489
490 call hecmw_tvec_nn(hecmesh, hecmat%NDOF, x, y, commtime)
491 call hecmw_matvec_nn_inner(hecmesh, hecmat, y, w, commtime)
492 call hecmw_ttvec_nn(hecmesh, hecmat%NDOF, w, y, commtime)
493
494 end subroutine hecmw_ttmattvec_nn
495
496 !C
497 !C***
498 !C*** hecmw_mat_diag_sr_nn
499 !C***
500 !C
501 subroutine hecmw_mat_diag_sr_nn(hecMESH, hecMAT, COMMtime)
502 use hecmw_util
504 implicit none
505 type (hecmwst_local_mesh), intent(in) :: hecmesh
506 type (hecmwst_matrix), intent(inout), target :: hecmat
507 real(kind=kreal), intent(inout), optional :: commtime
508 real(kind=kreal), allocatable :: w(:,:)
509 real(kind=kreal), pointer :: d(:)
510 integer(kind=kint) :: ip, ndof, i, j
511 real(kind=kreal) :: start_time, end_time
512 ndof = hecmat%NDOF
513 allocate(w(ndof*hecmat%NP,ndof))
514 d => hecmat%D
515 do ip= 1, hecmat%N
516 do i=1,ndof
517 do j=1,ndof
518 w(ndof*(ip-1)+i,j) = d(ndof*ndof*(ip-1)+(i-1)*ndof+j)
519 end do
520 end do
521 enddo
522 start_time= hecmw_wtime()
523 do i=1,ndof
524 call hecmw_update_m_r (hecmesh, w(:,i), hecmat%NP, ndof)
525 end do
526 end_time= hecmw_wtime()
527 if (present(commtime)) commtime = commtime + end_time - start_time
528 do ip= hecmat%N+1, hecmat%NP
529 do i=1,ndof
530 do j=1,ndof
531 d(ndof*ndof*(ip-1)+(i-1)*ndof+j) = w(ndof*(ip-1)+i,j)
532 end do
533 end do
534 enddo
535 deallocate(w)
536 end subroutine hecmw_mat_diag_sr_nn
537
538 subroutine hecmw_mat_add_nn(hecMAT1, hecMAT2, hecMAT3)
539 use hecmw_util
540 implicit none
541 type (hecmwst_matrix) :: hecmat1, hecmat2, hecmat3
542 integer(kind=kint) :: i
543
544 do i = 1, hecmat1%NP*hecmat1%NDOF*hecmat1%NDOF
545 hecmat3%D(i) = hecmat1%D(i) + hecmat2%D(i)
546 enddo
547
548 do i = 1, hecmat1%NPU*hecmat1%NDOF*hecmat1%NDOF
549 hecmat3%AU(i) = hecmat1%AU(i) + hecmat2%AU(i)
550 enddo
551
552 do i = 1, hecmat1%NPL*hecmat1%NDOF*hecmat1%NDOF
553 hecmat3%AL(i) = hecmat1%AL(i) + hecmat2%AL(i)
554 enddo
555 end subroutine hecmw_mat_add_nn
556
557 subroutine hecmw_mat_multiple_nn(hecMAT, alpha)
558 use hecmw_util
559 implicit none
560 type (hecmwst_matrix) :: hecmat
561 real(kind=kreal), intent(in) :: alpha
562 integer(kind=kint) :: i
563
564 do i = 1, hecmat%NP*hecmat%NDOF*hecmat%NDOF
565 hecmat%D(i) = alpha*hecmat%D(i)
566 enddo
567
568 do i = 1, hecmat%NPU*hecmat%NDOF*hecmat%NDOF
569 hecmat%AU(i) = alpha*hecmat%AU(i)
570 enddo
571
572 do i = 1, hecmat%NPL*hecmat%NDOF*hecmat%NDOF
573 hecmat%AL(i) = alpha*hecmat%AL(i)
574 enddo
575 end subroutine hecmw_mat_multiple_nn
576end module hecmw_solver_las_nn
Jagged Diagonal Matrix storage for vector processors. Original code was provided by JAMSTEC.
Definition: hecmw_jadm.f90:8
integer(kind=kint) function, public hecmw_jad_is_initialized()
Definition: hecmw_jadm.f90:56
subroutine, public hecmw_jad_matvec(hecmesh, hecmat, x, y, commtime)
Definition: hecmw_jadm.f90:61
subroutine, public hecmw_cmat_multvec_add(cmat, x, y, len)
integer(kind=kint) function, public hecmw_mat_get_flag_mpcmatvec(hecmat)
subroutine, public hecmw_mat_diag_sr_nn(hecmesh, hecmat, commtime)
subroutine, public hecmw_ttvec_nn(hecmesh, ndof, x, y, commtime)
subroutine, public hecmw_matvec_nn_set_async(hecmat)
subroutine, public hecmw_mat_add_nn(hecmat1, hecmat2, hecmat3)
subroutine, public hecmw_matvec_nn_unset_async
subroutine, public hecmw_matresid_nn(hecmesh, hecmat, x, b, r, commtime)
real(kind=kreal) function, public hecmw_rel_resid_l2_nn(hecmesh, hecmat, commtime)
subroutine, public hecmw_tvec_nn(hecmesh, ndof, x, y, commtime)
subroutine, public hecmw_ttmattvec_nn(hecmesh, hecmat, x, y, w, commtime)
subroutine, public hecmw_matvec_nn(hecmesh, hecmat, x, y, time_ax, commtime)
subroutine, public hecmw_mat_multiple_nn(hecmat, alpha)
subroutine hecmw_innerproduct_r(hecmesh, ndof, x, y, sum, commtime)
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
real(kind=kreal) function hecmw_wtime()
subroutine hecmw_update_m_r(hecmesh, val, n, m)