FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_jadm_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!-------------------------------------------------------------------------------
7
9 use hecmw_util
11 implicit none
12
13 private
14
15 public :: hecmw_jad_init_44
16 public :: hecmw_jad_finalize_44
18 public :: hecmw_jad_matvec_44
19
20 !C---------------------- AU&AL
21 real(kind=kreal), allocatable :: ajad(:)
22 integer(kind=kint), allocatable :: JAJAD(:)
23 integer(kind=kint), allocatable :: JADORD(:)
24 integer(kind=kint), allocatable :: IAJAD(:)
25 integer(kind=kint) :: MJAD
26 real(kind=kreal), allocatable :: wp1(:), wp2(:), wp3(:), wp4(:)
27 integer(kind=kint) :: INITIALIZED = 0
28
29contains
30
31 subroutine hecmw_jad_init_44(hecMAT)
32 type(hecmwst_matrix) :: hecmat
33 allocate(wp1(hecmat%NP), wp2(hecmat%NP), wp3(hecmat%NP), wp4(hecmat%NP))
34 allocate(ajad((hecmat%NPL+hecmat%NPU)*16))
35 allocate(jajad(hecmat%NPL+hecmat%NPU))
36 allocate(jadord(hecmat%NP))
37 allocate(iajad(hecmat%NP+1))
38 call repack(hecmat%N, hecmat, mjad, ajad, jajad, iajad, jadord)
39 initialized = 1
40 end subroutine hecmw_jad_init_44
41
43 deallocate(ajad)
44 deallocate(jajad)
45 deallocate(jadord)
46 deallocate(iajad)
47 deallocate(wp1,wp2,wp3,wp4)
48 initialized = 0
49 end subroutine hecmw_jad_finalize_44
50
52 integer(kind=kint) :: hecmw_jad_is_initialized_44
53 hecmw_jad_is_initialized_44 = initialized
55
56 subroutine hecmw_jad_matvec_44(hecMESH, hecMAT, X, Y, COMMtime)
57 type(hecmwst_local_mesh), intent(in) :: hecmesh
58 type(hecmwst_matrix), intent(in), target :: hecmat
59 real(kind=kreal), intent(in) :: x(:)
60 real(kind=kreal), intent(out) :: y(:)
61 real(kind=kreal), intent(inout) :: commtime
62 real(kind=kreal) :: start_time, end_time
63 real(kind=kreal), pointer :: d(:)
64 integer(kind=kint) :: i
65 real(kind=kreal) :: x1, x2, x3, x4
66
67 start_time= hecmw_wtime()
68 call hecmw_update_4_r (hecmesh, x, hecmat%NP)
69 end_time= hecmw_wtime()
70 commtime = commtime + end_time - start_time
71
72 d => hecmat%D
73
74 !$OMP PARALLEL PRIVATE(i)
75 !$OMP DO
76 do i= 1, hecmat%N
77 x1= x(4*i-3)
78 x2= x(4*i-2)
79 x3= x(4*i-1)
80 x4= x(4*i )
81 y(4*i-3)= d(16*i-15)*x1 + d(16*i-14)*x2 + d(16*i-13)*x3 + d(16*i-12)*x4
82 y(4*i-2)= d(16*i-11)*x1 + d(16*i-10)*x2 + d(16*i- 9)*x3 + d(16*i- 8)*x4
83 y(4*i-1)= d(16*i- 7)*x1 + d(16*i- 6)*x2 + d(16*i- 5)*x3 + d(16*i- 4)*x4
84 y(4*i )= d(16*i- 3)*x1 + d(16*i- 2)*x2 + d(16*i- 1)*x3 + d(16*i- 0)*x4
85 enddo
86 !$OMP END DO
87 !$OMP END PARALLEL
88 call matjad(hecmat%N, mjad, iajad, jajad, ajad, jadord, x, y, wp1, wp2, wp3, wp4)
89 end subroutine hecmw_jad_matvec_44
90
91 subroutine repack(N, hecMAT, MJAD, AJAD, JAJAD, IAJAD, JADORD)
92 use hecmw_util
93 !C---------------------------------
94 type (hecmwst_matrix) :: hecmat
95 !C----------------------
96 integer(kind = kint) :: n, mjad
97 real(kind = kreal), dimension(*) :: ajad
98 integer(kind = kint), dimension(*) :: jajad
99 integer(kind = kint), dimension(*) :: iajad
100 integer(kind = kint), dimension(*) :: jadord
101
102 integer(kind = kint) :: ijad, maxnz, minnz
103 integer(kind = kint) :: i, j, js, je, in, jc
104 integer(kind = kint), allocatable :: len(:), lenz(:), jadreord(:)
105
106 allocate(len(n))
107 allocate(jadreord(n))
108 do i=1,n
109 len(i)= hecmat%indexL(i) - hecmat%indexL(i-1) &
110 & + hecmat%indexU(i) - hecmat%indexU(i-1)
111 end do
112 maxnz=maxval(len(1:n))
113 minnz=minval(len(1:n))
114 mjad =maxnz
115 allocate(lenz(0:mjad))
116 lenz = 0
117 do i=1,n
118 lenz(len(i))=lenz(len(i))+1
119 enddo
120 do i=maxnz-1,minnz,-1
121 lenz(i)=lenz(i)+lenz(i+1)
122 enddo
123 do i=1,n
124 jadord(i)=lenz(len(i))
125 lenz(len(i))=lenz(len(i))-1
126 enddo
127 do i=1,n
128 jadreord(jadord(i))=i
129 enddo
130 do i=1,n
131 lenz(len(jadreord(i)))=i
132 enddo
133 do i=maxnz-1,1,-1
134 lenz(i)=max(lenz(i+1),lenz(i))
135 enddo
136 iajad(1)=1
137 do i=1,maxnz
138 iajad(i+1)=iajad(i)+lenz(i)
139 enddo
140 len=0
141 do i= 1, n
142 ijad=jadord(i)
143 js= hecmat%indexL(i-1) + 1
144 je= hecmat%indexL(i )
145 do j=js,je
146 in = hecmat%itemL(j)
147 len(ijad)=len(ijad)+1
148 jc=iajad(len(ijad))+ijad-1
149 ajad(jc*16-15:jc*16) = hecmat%AL(16*j-15:16*j)
150 jajad(jc) = in
151 end do
152 end do
153 do i= 1, n
154 ijad=jadord(i)
155 js= hecmat%indexU(i-1) + 1
156 je= hecmat%indexU(i )
157 do j=js,je
158 in = hecmat%itemU(j)
159 len(ijad)=len(ijad)+1
160 jc=iajad(len(ijad))+ijad-1
161 ajad(jc*16-15:jc*16) = hecmat%AU(16*j-15:16*j)
162 jajad(jc) = in
163 end do
164 end do
165 deallocate(len)
166 deallocate(jadreord)
167 deallocate(lenz)
168 end subroutine repack
169
170 subroutine matjad(N, MJAD, IAJAD, JAJAD, AJAD, JADORD, X, Y, W1, W2, W3, W4)
171 use hecmw_util
172 integer(kind=kint) :: n, mjad
173 integer(kind=kint) :: iajad(*), jajad(*), jadord(*)
174 real(kind=kreal) :: ajad(*), x(*), y(*), w1(*), w2(*), w3(*), w4(*)
175
176 integer(kind=kint) :: i, k, nz, ixx
177 real(kind=kreal) :: x1, x2, x3, x4
178
179 !$OMP PARALLEL PRIVATE(I,K,X1,X2,X3,X4,IXX)
180 !$OMP DO
181 do i=1,n
182 w1(i)=0.d0
183 w2(i)=0.d0
184 w3(i)=0.d0
185 w4(i)=0.d0
186 enddo
187 !$OMP END DO
188
189 do nz=1,mjad
190 !$OMP DO
191 do k=iajad(nz),iajad(nz+1)-1
192 x1=x(jajad(k)*4-3)
193 x2=x(jajad(k)*4-2)
194 x3=x(jajad(k)*4-1)
195 x4=x(jajad(k)*4 )
196 ixx = k-iajad(nz)+1
197 w1(ixx)=w1(ixx) + ajad(k*16-15)*x1 + ajad(k*16-14)*x2 + ajad(k*16-13)*x3 + ajad(k*16-12)*x4
198 w2(ixx)=w2(ixx) + ajad(k*16-11)*x1 + ajad(k*16-10)*x2 + ajad(k*16- 9)*x3 + ajad(k*16- 8)*x4
199 w3(ixx)=w3(ixx) + ajad(k*16- 7)*x1 + ajad(k*16- 6)*x2 + ajad(k*16- 5)*x3 + ajad(k*16- 4)*x4
200 w4(ixx)=w4(ixx) + ajad(k*16- 3)*x1 + ajad(k*16- 2)*x2 + ajad(k*16- 1)*x3 + ajad(k*16- 0)*x4
201 enddo
202 !$OMP END DO
203 enddo
204
205 !$OMP DO
206 do i=1,n
207 y(4*i-3)=y(4*i-3)+w1(jadord(i))
208 y(4*i-2)=y(4*i-2)+w2(jadord(i))
209 y(4*i-1)=y(4*i-1)+w3(jadord(i))
210 y(4*i )=y(4*i )+w4(jadord(i))
211 enddo
212 !$OMP END DO
213 !$OMP END PARALLEL
214 end subroutine matjad
215
216end module hecmw_jad_type_44
Jagged Diagonal Matrix storage for vector processors. Original code was provided by JAMSTEC.
subroutine, public hecmw_jad_init_44(hecmat)
subroutine, public hecmw_jad_finalize_44()
subroutine, public hecmw_jad_matvec_44(hecmesh, hecmat, x, y, commtime)
integer(kind=kint) function, public hecmw_jad_is_initialized_44()
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
real(kind=kreal) function hecmw_wtime()
subroutine hecmw_update_4_r(hecmesh, val, n)