FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
sparse_matrix.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 use hecmw_util
9 implicit none
10
11 private
12
15
19
20 public :: sparse_matrix
22 public :: sparse_matrix_init
24 public :: sparse_matrix_dump
27 public :: sparse_matrix_is_sym
28
29 integer(kind=kint), parameter :: sparse_matrix_type_csr=1
30 integer(kind=kint), parameter :: sparse_matrix_type_coo=2
31
32 integer(kind=kint), parameter :: sparse_matrix_symtype_asym=0
33 integer(kind=kint), parameter :: sparse_matrix_symtype_spd=1
34 integer(kind=kint), parameter :: sparse_matrix_symtype_sym=2
35
37 integer(kind=kint) :: type
38 integer(kind=kint) :: symtype
39 integer(kind=kint) :: n, n_loc
40 integer(kind=kint) :: nz
41 integer(kind=kint), pointer :: irn(:) => null()
42 integer(kind=kint), pointer :: jcn(:) => null()
43 real(kind=kreal), pointer :: a(:) => null()
44 real(kind=kreal), pointer :: rhs(:)
45 integer(kind=kint) :: offset
46 integer(kind=kint), pointer :: n_counts(:) => null()
47 integer(kind=kint), pointer :: displs(:) => null()
48 integer(kind=kint), pointer :: conv_ext(:) => null()
49 integer(kind=kint) :: iterlog
50 integer(kind=kint) :: timelog
51 logical :: is_initialized = .false.
52 end type sparse_matrix
53
54contains
55
56 !!! public subroutines
57
58 subroutine sparse_matrix_set_type(spMAT, type, symtype)
59 type (sparse_matrix), intent(inout) :: spmat
60 integer(kind=kint), intent(in) :: type
61 integer(kind=kint), intent(in) :: symtype
62 if (type == sparse_matrix_type_csr .or. &
63 type == sparse_matrix_type_coo) then
64 spmat%type = type
65 else
66 write(*,*) 'ERROR: unknown sparse matrix type'
67 stop
68 endif
69 if (symtype == sparse_matrix_symtype_asym .or. &
70 symtype == sparse_matrix_symtype_spd .or. &
71 symtype == sparse_matrix_symtype_sym) then
72 spmat%symtype = symtype
73 else
74 write(*,*) 'ERROR: unknown symmetry type for sparse matrix'
75 stop
76 endif
77 end subroutine sparse_matrix_set_type
78
79 subroutine sparse_matrix_init(spMAT, N_loc, NZ)
80 type (sparse_matrix), intent(inout) :: spmat
81 integer(kind=kint), intent(in) :: n_loc
82 integer(kind=kint), intent(in) :: nz
83 if (spmat%is_initialized) then
84 !write(*,*) "WARNING: sparse_matrix_init_prof: already initialized; freeing"
85 call sparse_matrix_finalize(spmat)
86 endif
87 call sparse_matrix_set_dimension(spmat, n_loc)
88 call sparse_matrix_set_offsets(spmat)
89 call sparse_matrix_allocate_arrays(spmat, nz)
90 spmat%is_initialized = .true.
91 end subroutine sparse_matrix_init
92
93 subroutine sparse_matrix_finalize(spMAT)
94 type (sparse_matrix), intent(inout) :: spmat
95 call sparse_matrix_free_arrays(spmat)
96 call sparse_matrix_clear(spmat)
97 spmat%is_initialized = .false.
98 end subroutine sparse_matrix_finalize
99
100 subroutine sparse_matrix_dump(spMAT)
101 type(sparse_matrix), intent(in) :: spmat
102 integer(kind=kint),save :: num_call
103 character(len=128) :: fname
104 integer(kind=kint) :: i,myrank
105 myrank=hecmw_comm_get_rank()
106 write(fname,'(A,I0,A,I0)') 'sparse-matrix-',num_call,'.coo.',myrank
107 num_call = num_call + 1
108 write(*,*) 'Dumping sparse matrix to ', fname
109 open(91,file=fname)
110 if (spmat%type == sparse_matrix_type_csr) then
111 if (spmat%symtype == sparse_matrix_symtype_asym) &
112 write(91,*) '%SPARSE MATRIX CSR ASYM'
113 if (spmat%symtype == sparse_matrix_symtype_spd) &
114 write(91,*) '%SPARSE MATRIX CSR SPD'
115 if (spmat%symtype == sparse_matrix_symtype_sym) &
116 write(91,*) '%SPARSE MATRIX CSR SYM'
117 write(91,*) spmat%N_loc, spmat%N_loc, spmat%NZ
118 do i=1,spmat%N_loc+1
119 write(91,*) spmat%IRN(i)
120 enddo
121 do i=1,spmat%NZ
122 write(91,*) spmat%JCN(i), spmat%A(i)
123 enddo
124 elseif (spmat%type == sparse_matrix_type_coo) then
125 if (spmat%symtype == sparse_matrix_symtype_asym) &
126 write(91,*) '%SPARSE MATRIX COO ASYM'
127 if (spmat%symtype == sparse_matrix_symtype_spd) &
128 write(91,*) '%SPARSE MATRIX COO SPD'
129 if (spmat%symtype == sparse_matrix_symtype_sym) &
130 write(91,*) '%SPARSE MATRIX COO SYM'
131 write(91,*) spmat%N_loc, spmat%N_loc, spmat%NZ
132 do i=1,spmat%NZ
133 write(91,*) spmat%IRN(i), spmat%JCN(i), spmat%A(i)
134 enddo
135 endif
136 close(91)
137 end subroutine sparse_matrix_dump
138
139 subroutine sparse_matrix_gather_rhs(spMAT, rhs_all)
140 type (sparse_matrix), intent(in) :: spmat
141 real(kind=kreal), intent(out) :: rhs_all(:)
142 integer(kind=kint) :: ierr,i
143 if (hecmw_comm_get_size() == 1) then
144 do i=1,spmat%N_loc
145 rhs_all(i) = spmat%rhs(i)
146 enddo
147 else
148 call hecmw_gatherv_real(spmat%rhs, spmat%N_loc, &
149 rhs_all, spmat%N_COUNTS, spmat%DISPLS, &
151 endif
152 end subroutine sparse_matrix_gather_rhs
153
154 subroutine sparse_matrix_scatter_rhs(spMAT, rhs_all)
155 type (sparse_matrix), intent(inout) :: spmat
156 real(kind=kreal), intent(in) :: rhs_all(:)
157 integer(kind=kint) :: ierr,i
158 if (hecmw_comm_get_size() == 1) then
159 do i=1,spmat%N_loc
160 spmat%rhs(i) = rhs_all(i)
161 enddo
162 else
163 call hecmw_scatterv_real( &
164 rhs_all, spmat%N_COUNTS, spmat%DISPLS, &
165 spmat%rhs, spmat%N_loc, &
167 endif
168 end subroutine sparse_matrix_scatter_rhs
169
170 function sparse_matrix_is_sym(spMAT)
171 type(sparse_matrix), intent(inout) :: spmat
172 logical :: sparse_matrix_is_sym
173 sparse_matrix_is_sym = (spmat%symtype > 0)
174 end function sparse_matrix_is_sym
175
176 !!! private subroutines
177
178 subroutine sparse_matrix_set_dimension(spMAT, N_loc)
179 type (sparse_matrix), intent(inout) :: spmat
180 integer(kind=kint), intent(in) :: n_loc
181 integer(kind=kint) :: ierr
182 spmat%N_loc = n_loc
183 if (hecmw_comm_get_size() == 1) then
184 spmat%N = spmat%N_loc
185 else
186 call hecmw_allreduce_int_1(spmat%N_loc, spmat%N, hecmw_sum, &
188 endif
189 end subroutine sparse_matrix_set_dimension
190
191 subroutine sparse_matrix_set_offsets(spMAT)
192 type (sparse_matrix), intent(inout) :: spmat
193 integer(kind=kint) :: i,ierr,nprocs,myrank
194 nprocs=hecmw_comm_get_size()
195 myrank=hecmw_comm_get_rank()
196 ! OFFSET
197 !if (myrank == 0) then
198 if (associated(spmat%N_COUNTS)) deallocate(spmat%N_COUNTS)
199 if (associated(spmat%DISPLS)) deallocate(spmat%DISPLS)
200 allocate(spmat%N_COUNTS(nprocs), spmat%DISPLS(nprocs), stat=ierr)
201 if (ierr /= 0) then
202 write(*,*) " Allocation error, spMAT%N_COUNTS, spMAT%DISPLS"
204 endif
205 !endif
206 if (nprocs > 1) then
207 call hecmw_allgather_int_1(spmat%N_loc, &
208 spmat%N_COUNTS, hecmw_comm_get_comm())
209 endif
210 spmat%DISPLS(1)=0
211 do i=1,nprocs-1
212 spmat%DISPLS(i+1)=spmat%DISPLS(i)+spmat%N_COUNTS(i)
213 enddo
214 spmat%offset = spmat%DISPLS(myrank+1)
215 end subroutine sparse_matrix_set_offsets
216
217 subroutine sparse_matrix_allocate_arrays(spMAT, NZ)
218 type(sparse_matrix), intent(inout) :: spmat
219 integer(kind=kint), intent(in) :: nz
220 integer(kind=kint) :: n_loc
221 integer(kind=kint) :: ierr
222 if (associated(spmat%IRN)) deallocate(spmat%IRN)
223 if (associated(spmat%JCN)) deallocate(spmat%JCN)
224 if (associated(spmat%A)) deallocate(spmat%A)
225 ierr = -1
226 n_loc=spmat%N_loc
227 if (spmat%type == sparse_matrix_type_csr) then
228 allocate(spmat%IRN(n_loc+1), spmat%JCN(nz), spmat%A(nz), stat=ierr)
229 elseif (spmat%type == sparse_matrix_type_coo) then
230 allocate(spmat%IRN(nz), spmat%JCN(nz), spmat%A(nz), stat=ierr)
231 endif
232 if (ierr /= 0) then
233 write(*,*) " Allocation error, spMAT%IRN, spMAT%JCN, spMAT%A"
235 endif
236 spmat%NZ = nz
237 end subroutine sparse_matrix_allocate_arrays
238
239 subroutine sparse_matrix_free_arrays(spMAT)
240 type (sparse_matrix), intent(inout) :: spmat
241 if (associated(spmat%IRN)) deallocate(spmat%IRN)
242 if (associated(spmat%JCN)) deallocate(spmat%JCN)
243 if (associated(spmat%A)) deallocate(spmat%A)
244 if (associated(spmat%N_COUNTS)) deallocate(spmat%N_COUNTS)
245 if (associated(spmat%DISPLS)) deallocate(spmat%DISPLS)
246 if (associated(spmat%conv_ext)) deallocate(spmat%conv_ext)
247 end subroutine sparse_matrix_free_arrays
248
249 subroutine sparse_matrix_clear(spMAT)
250 type(sparse_matrix), intent(inout) :: spmat
251 spmat%type = 0
252 spmat%symtype = 0
253 spmat%N = 0
254 spmat%N_loc = 0
255 spmat%NZ = 0
256 spmat%IRN => null()
257 spmat%JCN => null()
258 spmat%A => null()
259 spmat%offset = 0
260 spmat%N_COUNTS => null()
261 spmat%DISPLS => null()
262 spmat%conv_ext => null()
263 spmat%iterlog = 0
264 spmat%timelog = 0
265 end subroutine sparse_matrix_clear
266
267end module m_sparse_matrix
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=kint), parameter hecmw_sum
integer(kind=kint) function hecmw_comm_get_size()
integer(kind=kint) function hecmw_comm_get_comm()
integer(kind=4), parameter kreal
integer(kind=kint) function hecmw_comm_get_rank()
subroutine hecmw_abort(comm)
subroutine hecmw_gatherv_real(sbuf, sc, rbuf, rcs, disp, root, comm)
subroutine hecmw_scatterv_real(sbuf, scs, disp, rbuf, rc, root, comm)
subroutine hecmw_allgather_int_1(sval, rbuf, comm)
subroutine hecmw_allreduce_int_1(sval, rval, op, comm)
This module provides DOF based sparse matrix data structure (CSR and COO)
integer(kind=kint), parameter, public sparse_matrix_type_coo
integer(kind=kint), parameter, public sparse_matrix_type_csr
subroutine, public sparse_matrix_set_type(spmat, type, symtype)
integer(kind=kint), parameter, public sparse_matrix_symtype_sym
subroutine, public sparse_matrix_gather_rhs(spmat, rhs_all)
subroutine, public sparse_matrix_finalize(spmat)
subroutine, public sparse_matrix_init(spmat, n_loc, nz)
integer(kind=kint), parameter, public sparse_matrix_symtype_asym
integer(kind=kint), parameter, public sparse_matrix_symtype_spd
subroutine, public sparse_matrix_scatter_rhs(spmat, rhs_all)
subroutine, public sparse_matrix_dump(spmat)
logical function, public sparse_matrix_is_sym(spmat)