FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_ML_helper_33_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
6subroutine hecmw_ml_getrow_33(id, n_requested_rows, requested_rows, &
7 allocated_space, cols, values, row_lengths, ierr)
8 use hecmw_util
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
22 ierr = 0
23 call hecmw_mat_id_get(id, hecmat, hecmesh)
24 m = 1
25 do i = 1, n_requested_rows
26 row = requested_rows(i) + 1 ! '+1' for Fortran-numbering
27 inod = (row-1)/3 + 1
28 idof = row - (inod-1)*3
29 nl = (hecmat%indexL(inod) - hecmat%indexL(inod-1)) * 3
30 nd = 3
31 nu = (hecmat%indexU(inod) - hecmat%indexU(inod-1)) * 3
32 if (allocated_space < m + nl + nd + nu) return
33 start = m
34 js = hecmat%indexL(inod-1)+1
35 je = hecmat%indexL(inod)
36 do j = js, je
37 jj = hecmat%itemL(j)
38 do jdof = 1, 3
39 cols(m) = (jj-1)*3 + jdof - 1 ! '-1' for C-numbering
40 values(m) = hecmat%AL((j-1)*9 + (idof-1)*3 + jdof)
41 m = m+1
42 enddo
43 enddo
44 do jdof = 1, 3
45 cols(m) = (inod-1)*3 + jdof - 1 ! '-1' for C-numbering
46 values(m) = hecmat%D((inod-1)*9 + (idof-1)*3 + jdof)
47 m = m+1
48 enddo
49 js = hecmat%indexU(inod-1)+1
50 je = hecmat%indexU(inod)
51 do j = js, je
52 jj = hecmat%itemU(j)
53 do jdof = 1, 3
54 cols(m) = (jj-1)*3 + jdof - 1 ! '-1' for C-numbering
55 values(m) = hecmat%AU((j-1)*9 + (idof-1)*3 + jdof)
56 m = m+1
57 enddo
58 enddo
59 row_lengths(i) = m - start
60 enddo
61 ierr = 1
62end subroutine hecmw_ml_getrow_33
63
64subroutine hecmw_ml_matvec_33(id, in_length, p, out_length, ap, ierr)
65 use hecmw_util
66 use hecmw_mat_id
68 implicit none
69 integer(kind=kint), intent(in) :: id
70 integer(kind=kint), intent(in) :: in_length
71 real(kind=kreal), intent(in) :: p(in_length)
72 integer(kind=kint), intent(in) :: out_length
73 real(kind=kreal), intent(out) :: ap(out_length)
74 integer(kind=kint), intent(out) :: ierr
75 type(hecmwst_matrix), pointer :: hecMAT
76 type(hecmwst_local_mesh), pointer :: hecMESH
77 real(kind=kreal), allocatable :: w(:)
78 integer(kind=kint) :: i
79 call hecmw_mat_id_get(id, hecmat, hecmesh)
80 allocate(w(hecmat%NP*hecmat%NDOF))
81 do i = 1, hecmat%N*hecmat%NDOF
82 w(i) = p(i)
83 enddo
84 call hecmw_matvec(hecmesh, hecmat, w, ap)
85 deallocate(w)
86 ierr = 0
87end subroutine hecmw_ml_matvec_33
88
89subroutine hecmw_ml_comm_33(id, x, ierr)
90 use hecmw_util
91 use hecmw_mat_id
93 implicit none
94 integer(kind=kint), intent(in) :: id
95 real(kind=kreal), intent(inout) :: x(*)
96 integer(kind=kint), intent(out) :: ierr
97 type(hecmwst_matrix), pointer :: hecMAT
98 type(hecmwst_local_mesh), pointer :: hecMESH
99 call hecmw_mat_id_get(id, hecmat, hecmesh)
100 call hecmw_update_3_r (hecmesh, x, hecmesh%n_node)
101 ierr = 0
102end subroutine hecmw_ml_comm_33
103
105 use hecmw_util
106 use hecmw_mat_id
108 implicit none
109 integer(kind=kint), intent(in) :: id
110 integer(kind=kint), intent(out) :: ierr
111 type(hecmwst_matrix), pointer :: hecMAT
112 type(hecmwst_local_mesh), pointer :: hecMESH
113 call hecmw_mat_id_get(id, hecmat, hecmesh)
114 call hecmw_precond_diag_33_setup(hecmat)
115 ierr = 0
117
118subroutine hecmw_ml_smoother_diag_apply_33(id, x_length, x, rhs_length, rhs, ierr)
119 use hecmw_util
120 use hecmw_mat_id
124 implicit none
125 integer(kind=kint), intent(in) :: id
126 integer(kind=kint), intent(in) :: x_length
127 real(kind=kreal), intent(inout) :: x(x_length)
128 integer(kind=kint), intent(in) :: rhs_length
129 real(kind=kreal), intent(in) :: rhs(rhs_length)
130 integer(kind=kint), intent(out) :: ierr
131 type(hecmwst_matrix), pointer :: hecMAT
132 type(hecmwst_local_mesh), pointer :: hecMESH
133
134 real(kind=kreal), allocatable :: resid(:)
135 integer(kind=kint) :: i
136 real(kind=kreal) :: commtime
137 integer(kind=kint) :: num_sweeps, i_sweep
138
139 call hecmw_mat_id_get(id, hecmat, hecmesh)
140 num_sweeps = hecmw_mat_get_solver_opt6(hecmat)
141 allocate(resid(hecmat%NP * hecmat%NDOF))
142 do i_sweep = 1, num_sweeps
143 ! {resid} = {rhs} - [A] {x}
144 call hecmw_matresid_33(hecmesh, hecmat, x, rhs, resid, commtime)
145 ! {delta_x} = [M]^-1 {resid}
147 ! {x} = {x} + {delta_x}
148 do i=1,x_length
149 x(i) = x(i) + resid(i)
150 enddo
151 enddo
152 deallocate(resid)
153 ierr = 0
155
157 use hecmw_util
158 use hecmw_mat_id
160 implicit none
161 integer(kind=kint), intent(in) :: id
162 integer(kind=kint), intent(out) :: ierr
163 type(hecmwst_matrix), pointer :: hecMAT
164 type(hecmwst_local_mesh), pointer :: hecMESH
165 call hecmw_mat_id_get(id, hecmat, hecmesh)
167 ierr = 0
169
171 use hecmw_util
172 use hecmw_mat_id
174 implicit none
175 integer(kind=kint), intent(in) :: id
176 integer(kind=kint), intent(out) :: ierr
177 type(hecmwst_matrix), pointer :: hecMAT
178 type(hecmwst_local_mesh), pointer :: hecMESH
179 call hecmw_mat_id_get(id, hecmat, hecmesh)
180 call hecmw_precond_ssor_33_setup(hecmat)
181 ierr = 0
183
184subroutine hecmw_ml_smoother_ssor_apply_33(id, x_length, x, rhs_length, rhs, ierr)
185 use hecmw_util
186 use hecmw_mat_id
190 implicit none
191 integer(kind=kint), intent(in) :: id
192 integer(kind=kint), intent(in) :: x_length
193 real(kind=kreal), intent(inout) :: x(x_length)
194 integer(kind=kint), intent(in) :: rhs_length
195 real(kind=kreal), intent(in) :: rhs(rhs_length)
196 integer(kind=kint), intent(out) :: ierr
197 type(hecmwst_matrix), pointer :: hecMAT
198 type(hecmwst_local_mesh), pointer :: hecMESH
199
200 real(kind=kreal), allocatable :: resid(:)
201 integer(kind=kint) :: i
202 real(kind=kreal) :: commtime
203 integer(kind=kint) :: num_sweeps, i_sweep
204
205 call hecmw_mat_id_get(id, hecmat, hecmesh)
206 num_sweeps = hecmw_mat_get_solver_opt6(hecmat)
207 allocate(resid(hecmat%NP * hecmat%NDOF))
208 do i_sweep = 1, num_sweeps
209 ! {resid} = {rhs} - [A] {x}
210 call hecmw_matresid_33(hecmesh, hecmat, x, rhs, resid, commtime)
211 ! {delta_x} = [M]^-1 {resid}
213 ! {x} = {x} + {delta_x}
214 do i=1,x_length
215 x(i) = x(i) + resid(i)
216 enddo
217 enddo
218 deallocate(resid)
219 ierr = 0
221
223 use hecmw_util
224 use hecmw_mat_id
226 implicit none
227 integer(kind=kint), intent(in) :: id
228 integer(kind=kint), intent(out) :: ierr
229 type(hecmwst_matrix), pointer :: hecMAT
230 type(hecmwst_local_mesh), pointer :: hecMESH
231 call hecmw_mat_id_get(id, hecmat, hecmesh)
232 call hecmw_precond_ssor_33_clear(hecmat)
233 ierr = 0
subroutine hecmw_ml_smoother_ssor_apply_33(id, x_length, x, rhs_length, rhs, ierr)
subroutine hecmw_ml_getrow_33(id, n_requested_rows, requested_rows, allocated_space, cols, values, row_lengths, ierr)
subroutine hecmw_ml_smoother_diag_setup_33(id, ierr)
subroutine hecmw_ml_smoother_ssor_clear_33(id, ierr)
subroutine hecmw_ml_smoother_diag_clear_33(id, ierr)
subroutine hecmw_ml_smoother_diag_apply_33(id, x_length, x, rhs_length, rhs, ierr)
subroutine hecmw_ml_comm_33(id, x, ierr)
subroutine hecmw_ml_matvec_33(id, in_length, p, out_length, ap, ierr)
subroutine hecmw_ml_smoother_ssor_setup_33(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_33_clear()
subroutine, public hecmw_precond_diag_33_setup(hecmat)
subroutine, public hecmw_precond_diag_33_apply(ww)
subroutine, public hecmw_precond_ssor_33_apply(zp)
subroutine, public hecmw_precond_ssor_33_setup(hecmat)
subroutine, public hecmw_precond_ssor_33_clear(hecmat)
subroutine, public hecmw_matresid_33(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_3_r(hecmesh, val, n)