FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_ML_helper_nn_f.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 subroutine hecmw_ml_getrow_nn(id, n_requested_rows, requested_rows, &
7  allocated_space, cols, values, row_lengths, ierr)
8  use hecmw_util
9  use hecmw_mat_id
10  implicit none
11  integer(kind=kint), intent(in) :: id
12  integer(kind=kint), intent(in) :: n_requested_rows
13  integer(kind=kint), intent(in) :: requested_rows(n_requested_rows)
14  integer(kind=kint), intent(in) :: allocated_space
15  integer(kind=kint), intent(out) :: cols(allocated_space)
16  real(kind=kreal), intent(out) :: values(allocated_space)
17  integer(kind=kint), intent(out) :: row_lengths(n_requested_rows)
18  integer(kind=kint), intent(out) :: ierr
19  type(hecmwst_matrix), pointer :: hecMAT
20  type(hecmwst_local_mesh), pointer :: hecMESH
21  integer(kind=kint) :: m, i, row, inod, idof, nl, nd, nu, js, je, j, jj, jdof, start, ndof
22  ierr = 0
23  call hecmw_mat_id_get(id, hecmat, hecmesh)
24  ndof = hecmat%NDOF
25  m = 1
26  do i = 1, n_requested_rows
27  row = requested_rows(i) + 1 ! '+1' for Fortran-numbering
28  inod = (row-1)/ndof + 1
29  idof = row - (inod-1)*ndof
30  nl = (hecmat%indexL(inod) - hecmat%indexL(inod-1)) * ndof
31  nd = ndof
32  nu = (hecmat%indexU(inod) - hecmat%indexU(inod-1)) * ndof
33  if (allocated_space < m + nl + nd + nu) return
34  start = m
35  js = hecmat%indexL(inod-1)+1
36  je = hecmat%indexL(inod)
37  do j = js, je
38  jj = hecmat%itemL(j)
39  do jdof = 1, ndof
40  cols(m) = (jj-1)*ndof + jdof - 1 ! '-1' for C-numbering
41  values(m) = hecmat%AL((j-1)*ndof*ndof + (idof-1)*ndof + jdof)
42  m = m+1
43  enddo
44  enddo
45  do jdof = 1, ndof
46  cols(m) = (inod-1)*ndof + jdof - 1 ! '-1' for C-numbering
47  values(m) = hecmat%D((inod-1)*ndof*ndof + (idof-1)*ndof + jdof)
48  m = m+1
49  enddo
50  js = hecmat%indexU(inod-1)+1
51  je = hecmat%indexU(inod)
52  do j = js, je
53  jj = hecmat%itemU(j)
54  do jdof = 1, ndof
55  cols(m) = (jj-1)*ndof + jdof - 1 ! '-1' for C-numbering
56  values(m) = hecmat%AU((j-1)*ndof*ndof + (idof-1)*ndof + jdof)
57  m = m+1
58  enddo
59  enddo
60  row_lengths(i) = m - start
61  enddo
62  ierr = 1
63 end subroutine hecmw_ml_getrow_nn
64 
65 subroutine hecmw_ml_matvec_nn(id, in_length, p, out_length, ap, ierr)
66  use hecmw_util
67  use hecmw_mat_id
69  implicit none
70  integer(kind=kint), intent(in) :: id
71  integer(kind=kint), intent(in) :: in_length
72  real(kind=kreal), intent(in) :: p(in_length)
73  integer(kind=kint), intent(in) :: out_length
74  real(kind=kreal), intent(out) :: ap(out_length)
75  integer(kind=kint), intent(out) :: ierr
76  type(hecmwst_matrix), pointer :: hecMAT
77  type(hecmwst_local_mesh), pointer :: hecMESH
78  real(kind=kreal), allocatable :: w(:)
79  integer(kind=kint) :: i
80  call hecmw_mat_id_get(id, hecmat, hecmesh)
81  allocate(w(hecmat%NP*hecmat%NDOF))
82  do i = 1, hecmat%N*hecmat%NDOF
83  w(i) = p(i)
84  enddo
85  call hecmw_matvec(hecmesh, hecmat, w, ap)
86  deallocate(w)
87  ierr = 0
88 end subroutine hecmw_ml_matvec_nn
89 
90 subroutine hecmw_ml_comm_nn(id, x, ierr)
91  use hecmw_util
92  use hecmw_mat_id
93  use m_hecmw_comm_f
94  implicit none
95  integer(kind=kint), intent(in) :: id
96  real(kind=kreal), intent(inout) :: x(*)
97  integer(kind=kint), intent(out) :: ierr
98  type(hecmwst_matrix), pointer :: hecMAT
99  type(hecmwst_local_mesh), pointer :: hecMESH
100  call hecmw_mat_id_get(id, hecmat, hecmesh)
101  call hecmw_update_m_r (hecmesh, x, hecmat%NP, hecmat%NDOF)
102  ierr = 0
103 end subroutine hecmw_ml_comm_nn
104 
106  use hecmw_util
107  use hecmw_mat_id
109  implicit none
110  integer(kind=kint), intent(in) :: id
111  integer(kind=kint), intent(out) :: ierr
112  type(hecmwst_matrix), pointer :: hecMAT
113  type(hecmwst_local_mesh), pointer :: hecMESH
114  call hecmw_mat_id_get(id, hecmat, hecmesh)
115  call hecmw_precond_diag_nn_setup(hecmat)
116  ierr = 0
117 end subroutine hecmw_ml_smoother_diag_setup_nn
118 
119 subroutine hecmw_ml_smoother_diag_apply_nn(id, x_length, x, rhs_length, rhs, ierr)
120  use hecmw_util
121  use hecmw_mat_id
125  implicit none
126  integer(kind=kint), intent(in) :: id
127  integer(kind=kint), intent(in) :: x_length
128  real(kind=kreal), intent(inout) :: x(x_length)
129  integer(kind=kint), intent(in) :: rhs_length
130  real(kind=kreal), intent(in) :: rhs(rhs_length)
131  integer(kind=kint), intent(out) :: ierr
132  type(hecmwst_matrix), pointer :: hecMAT
133  type(hecmwst_local_mesh), pointer :: hecMESH
134 
135  real(kind=kreal), allocatable :: resid(:)
136  integer(kind=kint) :: i
137  real(kind=kreal) :: commtime
138  integer(kind=kint) :: num_sweeps, i_sweep
139 
140  call hecmw_mat_id_get(id, hecmat, hecmesh)
141  num_sweeps = hecmw_mat_get_solver_opt6(hecmat)
142  allocate(resid(hecmat%NP * hecmat%NDOF))
143  do i_sweep = 1, num_sweeps
144  ! {resid} = {rhs} - [A] {x}
145  call hecmw_matresid_nn(hecmesh, hecmat, x, rhs, resid, commtime)
146  ! {delta_x} = [M]^-1 {resid}
147  call hecmw_precond_diag_nn_apply(resid, hecmat%NDOF)
148  ! {x} = {x} + {delta_x}
149  do i=1,x_length
150  x(i) = x(i) + resid(i)
151  enddo
152  enddo
153  deallocate(resid)
154  ierr = 0
155 end subroutine hecmw_ml_smoother_diag_apply_nn
156 
158  use hecmw_util
159  use hecmw_mat_id
161  implicit none
162  integer(kind=kint), intent(in) :: id
163  integer(kind=kint), intent(out) :: ierr
164  type(hecmwst_matrix), pointer :: hecMAT
165  type(hecmwst_local_mesh), pointer :: hecMESH
166  call hecmw_mat_id_get(id, hecmat, hecmesh)
168  ierr = 0
169 end subroutine hecmw_ml_smoother_diag_clear_nn
170 
172  use hecmw_util
173  use hecmw_mat_id
175  implicit none
176  integer(kind=kint), intent(in) :: id
177  integer(kind=kint), intent(out) :: ierr
178  type(hecmwst_matrix), pointer :: hecMAT
179  type(hecmwst_local_mesh), pointer :: hecMESH
180  call hecmw_mat_id_get(id, hecmat, hecmesh)
181  call hecmw_precond_ssor_nn_setup(hecmat)
182  ierr = 0
183 end subroutine hecmw_ml_smoother_ssor_setup_nn
184 
185 subroutine hecmw_ml_smoother_ssor_apply_nn(id, x_length, x, rhs_length, rhs, ierr)
186  use hecmw_util
187  use hecmw_mat_id
191  implicit none
192  integer(kind=kint), intent(in) :: id
193  integer(kind=kint), intent(in) :: x_length
194  real(kind=kreal), intent(inout) :: x(x_length)
195  integer(kind=kint), intent(in) :: rhs_length
196  real(kind=kreal), intent(in) :: rhs(rhs_length)
197  integer(kind=kint), intent(out) :: ierr
198  type(hecmwst_matrix), pointer :: hecMAT
199  type(hecmwst_local_mesh), pointer :: hecMESH
200 
201  real(kind=kreal), allocatable :: resid(:)
202  integer(kind=kint) :: i
203  real(kind=kreal) :: commtime
204  integer(kind=kint) :: num_sweeps, i_sweep
205 
206  call hecmw_mat_id_get(id, hecmat, hecmesh)
207  num_sweeps = hecmw_mat_get_solver_opt6(hecmat)
208  allocate(resid(hecmat%NP * hecmat%NDOF))
209  do i_sweep = 1, num_sweeps
210  ! {resid} = {rhs} - [A] {x}
211  call hecmw_matresid_nn(hecmesh, hecmat, x, rhs, resid, commtime)
212  ! {delta_x} = [M]^-1 {resid}
213  call hecmw_precond_ssor_nn_apply(resid, hecmat%NDOF)
214  ! {x} = {x} + {delta_x}
215  do i=1,x_length
216  x(i) = x(i) + resid(i)
217  enddo
218  enddo
219  deallocate(resid)
220  ierr = 0
221 end subroutine hecmw_ml_smoother_ssor_apply_nn
222 
224  use hecmw_util
225  use hecmw_mat_id
227  implicit none
228  integer(kind=kint), intent(in) :: id
229  integer(kind=kint), intent(out) :: ierr
230  type(hecmwst_matrix), pointer :: hecMAT
231  type(hecmwst_local_mesh), pointer :: hecMESH
232  call hecmw_mat_id_get(id, hecmat, hecmesh)
233  call hecmw_precond_ssor_nn_clear(hecmat)
234  ierr = 0
235 end subroutine hecmw_ml_smoother_ssor_clear_nn
subroutine hecmw_ml_smoother_diag_apply_nn(id, x_length, x, rhs_length, rhs, ierr)
subroutine hecmw_ml_smoother_ssor_apply_nn(id, x_length, x, rhs_length, rhs, ierr)
subroutine hecmw_ml_comm_nn(id, x, ierr)
subroutine hecmw_ml_getrow_nn(id, n_requested_rows, requested_rows, allocated_space, cols, values, row_lengths, ierr)
subroutine hecmw_ml_matvec_nn(id, in_length, p, out_length, ap, ierr)
subroutine hecmw_ml_smoother_ssor_setup_nn(id, ierr)
subroutine hecmw_ml_smoother_diag_clear_nn(id, ierr)
subroutine hecmw_ml_smoother_diag_setup_nn(id, ierr)
subroutine hecmw_ml_smoother_ssor_clear_nn(id, ierr)
subroutine, public hecmw_mat_id_get(id, hecMAT, hecMESH)
integer(kind=kint) function, public hecmw_mat_get_solver_opt6(hecMAT)
subroutine, public hecmw_precond_diag_nn_setup(hecMAT)
subroutine, public hecmw_precond_diag_nn_clear()
subroutine, public hecmw_precond_diag_nn_apply(WW, NDOF)
subroutine, public hecmw_precond_ssor_nn_setup(hecMAT)
subroutine, public hecmw_precond_ssor_nn_clear(hecMAT)
subroutine, public hecmw_precond_ssor_nn_apply(ZP, NDOF)
subroutine, public hecmw_matresid_nn(hecMESH, hecMAT, X, B, R, COMMtime)
subroutine, public hecmw_matvec(hecMESH, hecMAT, X, Y, COMMtime)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
subroutine hecmw_update_m_r(hecMESH, val, n, m)