38 real(kind=
kreal),
intent(in) :: x(:)
39 real(kind=
kreal),
intent(out) :: y(:)
40 real(kind=
kreal),
intent(inout) :: time_ax
41 real(kind=
kreal),
intent(inout),
optional :: commtime
43 real(kind=
kreal) :: start_time, end_time, tcomm
44 integer(kind=kint) :: i, j, js, je, in
45 real(kind=
kreal) :: yv1, yv2, yv3, x1, x2, x3
46 real(kind=
kreal) :: yv4, yv5, yv6, x4, x5, x6
48 integer(kind=kint) :: n, np
49 integer(kind=kint),
pointer :: indexl(:), iteml(:), indexu(:), itemu(:)
50 real(kind=
kreal),
pointer :: al(:), au(:), d(:)
53 integer,
parameter :: numofblockperthread = 100
54 logical,
save :: isfirst = .true.
55 integer,
save :: numofthread = 1
56 integer,
save,
allocatable :: startpos(:), endpos(:)
57 integer(kind=kint),
save :: sectorcachesize0, sectorcachesize1
58 integer(kind=kint) :: threadnum, blocknum, numofblock
59 integer(kind=kint) :: numofelement, elementcount, blockindex
60 real(kind=
kreal) :: numofelementperblock
68 time_ax = time_ax + end_time - start_time - tcomm
69 if (
present(commtime)) commtime = commtime + tcomm
74 indexl => hecmat%indexL
75 indexu => hecmat%indexU
83 if (isfirst .eqv. .true.)
then
85 numofblock = numofthread * numofblockperthread
86 allocate (startpos(0 : numofblock - 1), endpos(0 : numofblock - 1))
87 numofelement = n + indexl(n) + indexu(n)
88 numofelementperblock = dble(numofelement) / numofblock
91 startpos(blocknum) = 1
93 elementcount = elementcount + 1
94 elementcount = elementcount + (indexl(i) - indexl(i-1))
95 elementcount = elementcount + (indexu(i) - indexu(i-1))
96 if (elementcount > (blocknum + 1) * numofelementperblock)
then
100 blocknum = blocknum + 1
101 startpos(blocknum) = i + 1
102 if (blocknum == (numofblock - 1))
exit
109 do i= blocknum+1, numofblock-1
117 sectorcachesize0, sectorcachesize1)
126 if (
present(commtime)) commtime = commtime + end_time - start_time
141 do blocknum = 0 , numofblockperthread - 1
142 blockindex = blocknum * numofthread + threadnum
143 do i = startpos(blockindex), endpos(blockindex)
150 yv1= d(36*i-35)*x1 + d(36*i-34)*x2 + d(36*i-33)*x3 + d(36*i-32)*x4 + d(36*i-31)*x5 + d(36*i-30)*x6
151 yv2= d(36*i-29)*x1 + d(36*i-28)*x2 + d(36*i-27)*x3 + d(36*i-26)*x4 + d(36*i-25)*x5 + d(36*i-24)*x6
152 yv3= d(36*i-23)*x1 + d(36*i-22)*x2 + d(36*i-21)*x3 + d(36*i-20)*x4 + d(36*i-19)*x5 + d(36*i-18)*x6
153 yv4= d(36*i-17)*x1 + d(36*i-16)*x2 + d(36*i-15)*x3 + d(36*i-14)*x4 + d(36*i-13)*x5 + d(36*i-12)*x6
154 yv5= d(36*i-11)*x1 + d(36*i-10)*x2 + d(36*i-9 )*x3 + d(36*i-8 )*x4 + d(36*i-7 )*x5 + d(36*i-6 )*x6
155 yv6= d(36*i-5 )*x1 + d(36*i-4 )*x2 + d(36*i-3 )*x3 + d(36*i-2 )*x4 + d(36*i-1 )*x5 + d(36*i )*x6
167 yv1= yv1 + 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
168 yv2= yv2 + 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
169 yv3= yv3 + 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
170 yv4= yv4 + 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
171 yv5= yv5 + 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
172 yv6= yv6 + 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
184 yv1= yv1 + 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
185 yv2= yv2 + 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
186 yv3= yv3 + 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
187 yv4= yv4 + 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
188 yv5= yv5 + 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
189 yv6= yv6 + 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
208 time_ax = time_ax + end_time - start_time
212 if (hecmat%cmat%n_val > 0)
then
228 real(kind=
kreal),
intent(in) :: x(:), b(:)
229 real(kind=
kreal),
intent(out) :: r(:)
230 real(kind=
kreal),
intent(inout),
optional :: commtime
232 integer(kind=kint) :: i
233 real(kind=
kreal) :: tcomm
237 if (
present(commtime)) commtime = commtime + tcomm
238 do i = 1, hecmat%N * 6
256 real(kind=
kreal),
intent(inout),
optional :: commtime
258 real(kind=
kreal),
allocatable :: r(:)
259 real(kind=
kreal) :: bnorm2, rnorm2
260 real(kind=
kreal) :: tcomm
262 allocate(r(hecmat%NDOF*hecmat%NP))
266 hecmat%B, hecmat%B, bnorm2, tcomm)
267 if (bnorm2 == 0.d0)
then
274 if (
present(commtime)) commtime = commtime + tcomm
289 real(kind=
kreal),
intent(in) :: x(:)
290 real(kind=
kreal),
intent(out) :: y(:)
291 real(kind=
kreal),
intent(inout) :: commtime
293 real(kind=
kreal) :: start_time, end_time
294 integer(kind=kint) :: i
299 commtime = commtime + end_time - start_time
301 do i= 1, hecmesh%nn_internal * hecmesh%n_dof
327 real(kind=
kreal),
intent(in) :: x(:)
328 real(kind=
kreal),
intent(out) :: y(:)
329 real(kind=
kreal),
intent(inout) :: commtime
331 real(kind=
kreal) :: start_time, end_time
332 integer(kind=kint) :: i
337 commtime = commtime + end_time - start_time
339 do i= 1, hecmesh%nn_internal * hecmesh%n_dof
365 real(kind=
kreal),
intent(in) :: x(:)
366 real(kind=
kreal),
intent(out) :: y(:), w(:)
367 real(kind=
kreal),
intent(inout) :: commtime
Jagged Diagonal Matrix storage for vector processors. Original code was provided by JAMSTEC.
subroutine, public hecmw_jad_matvec(hecMESH, hecMAT, X, Y, COMMtime)
integer(kind=kint) function, public hecmw_mat_get_usejad(hecMAT)
subroutine, public hecmw_matresid_66(hecMESH, hecMAT, X, B, R, COMMtime)
subroutine, public hecmw_tvec_66(hecMESH, X, Y, COMMtime)
subroutine, public hecmw_ttvec_66(hecMESH, X, Y, COMMtime)
real(kind=kreal) function, public hecmw_rel_resid_l2_66(hecMESH, hecMAT, COMMtime)
subroutine, public hecmw_matvec_66(hecMESH, hecMAT, X, Y, time_Ax, COMMtime)
subroutine, public hecmw_ttmattvec_66(hecMESH, hecMAT, X, Y, W, COMMtime)
subroutine hecmw_innerproduct_r(hecMESH, ndof, X, Y, sum, COMMtime)
subroutine, public hecmw_tuning_fx_calc_sector_cache(N, NDOF, sectorCacheSize0, sectorCacheSize1)
integer(kind=4), parameter kreal
real(kind=kreal) function hecmw_wtime()
subroutine hecmw_update_6_r(hecMESH, val, n)