FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_solver_las.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 
7  use hecmw_util
14 
15  implicit none
16 
17  private
18 
19  public :: hecmw_matvec
20  public :: hecmw_matvec_set_async
21  public :: hecmw_matvec_unset_async
22  public :: hecmw_matresid
23  public :: hecmw_rel_resid_l2
24  public :: hecmw_tvec
25  public :: hecmw_ttvec
26  public :: hecmw_ttmattvec
27  public :: hecmw_matvec_clear_timer
28  public :: hecmw_matvec_get_timer
29  public :: hecmw_mat_diag_sr
30  public :: hecmw_mat_add
31  public :: hecmw_mat_multiple
32 
33  real(kind=kreal), save :: time_ax = 0.d0
34 
35 contains
36 
37  !C
38  !C***
39  !C*** hecmw_matvec
40  !C***
41  !C
42  subroutine hecmw_matvec (hecMESH, hecMAT, X, Y, COMMtime)
43  use hecmw_util
44 
45  implicit none
46  type (hecmwst_local_mesh), intent(in) :: hecmesh
47  type (hecmwst_matrix), intent(in), target :: hecmat
48  real(kind=kreal), intent(in) :: x(:)
49  real(kind=kreal), intent(out) :: y(:)
50  real(kind=kreal), intent(inout), optional :: commtime
51  select case(hecmat%NDOF)
52  case (3)
53  call hecmw_matvec_33(hecmesh, hecmat, x, y, time_ax, commtime)
54  case (4)
55  call hecmw_matvec_44(hecmesh, hecmat, x, y, time_ax,commtime)
56  case (6)
57  call hecmw_matvec_66(hecmesh, hecmat, x, y, time_ax,commtime)
58  case default
59  call hecmw_matvec_nn(hecmesh, hecmat, x, y, time_ax, commtime)
60  end select
61 
62  end subroutine hecmw_matvec
63 
64  !C
65  !C***
66  !C*** hecmw_matvec_set_async
67  !C***
68  !C
69  subroutine hecmw_matvec_set_async (hecMAT)
70  use hecmw_util
71  implicit none
72  type (hecmwst_matrix), intent(in) :: hecmat
73 
74  end subroutine hecmw_matvec_set_async
75 
76  !C
77  !C***
78  !C*** hecmw_matvec_unset_async
79  !C***
80  !C
82  implicit none
83  end subroutine hecmw_matvec_unset_async
84 
85  !C
86  !C***
87  !C*** hecmw_matresid
88  !C***
89  !C
90  subroutine hecmw_matresid (hecMESH, hecMAT, X, B, R, COMMtime)
91  use hecmw_util
94  implicit none
95  type (hecmwst_local_mesh), intent(in) :: hecmesh
96  type (hecmwst_matrix), intent(in) :: hecmat
97  real(kind=kreal), intent(in) :: x(:), b(:)
98  real(kind=kreal), intent(out) :: r(:)
99  real(kind=kreal), intent(inout), optional :: commtime
100 
101  select case(hecmat%NDOF)
102  case (3)
103  call hecmw_matresid_33(hecmesh, hecmat, x, b, r, commtime)
104  case default
105  call hecmw_matresid_nn(hecmesh, hecmat, x, b, r, commtime)
106  end select
107  end subroutine hecmw_matresid
108 
109  !C
110  !C***
111  !C*** hecmw_rel_resid_L2
112  !C***
113  !C
114  function hecmw_rel_resid_l2 (hecMESH, hecMAT, COMMtime)
115  use hecmw_util
117 
118  implicit none
119  real(kind=kreal) :: hecmw_rel_resid_l2
120  type ( hecmwst_local_mesh ), intent(in) :: hecmesh
121  type ( hecmwst_matrix ), intent(in) :: hecmat
122  real(kind=kreal), intent(inout), optional :: commtime
123 
124  real(kind=kreal), allocatable :: r(:)
125  real(kind=kreal) :: bnorm2, rnorm2
126  real(kind=kreal) :: tcomm
127 
128  allocate(r(hecmat%NDOF*hecmat%NP))
129 
130  tcomm = 0.d0
131  call hecmw_innerproduct_r(hecmesh, hecmat%NDOF, &
132  hecmat%B, hecmat%B, bnorm2, tcomm)
133  if (bnorm2 == 0.d0) then
134  bnorm2 = 1.d0
135  endif
136  call hecmw_matresid(hecmesh, hecmat, hecmat%X, hecmat%B, r, tcomm)
137  call hecmw_innerproduct_r(hecmesh, hecmat%NDOF, r, r, rnorm2, tcomm)
138  hecmw_rel_resid_l2 = sqrt(rnorm2 / bnorm2)
139 
140  if (present(commtime)) commtime = commtime + tcomm
141 
142  deallocate(r)
143  end function hecmw_rel_resid_l2
144 
145  !C
146  !C***
147  !C*** hecmw_Tvec
148  !C***
149  !C
150  subroutine hecmw_tvec (hecMESH, ndof, X, Y, COMMtime)
151  use hecmw_util
154 
155  implicit none
156  type (hecmwst_local_mesh), intent(in) :: hecmesh
157  integer(kind=kint), intent(in) :: ndof
158  real(kind=kreal), intent(in) :: x(:)
159  real(kind=kreal), intent(out) :: y(:)
160  real(kind=kreal), intent(inout) :: commtime
161 
162  select case(ndof)
163  case (3)
164  call hecmw_tvec_33(hecmesh, x, y, commtime)
165  case default
166  call hecmw_tvec_nn(hecmesh, ndof, x, y, commtime)
167  end select
168 
169  end subroutine hecmw_tvec
170 
171  !C
172  !C***
173  !C*** hecmw_Ttvec
174  !C***
175  !C
176  subroutine hecmw_ttvec (hecMESH, ndof, X, Y, COMMtime)
177  use hecmw_util
180  implicit none
181  type (hecmwst_local_mesh), intent(in) :: hecmesh
182  integer(kind=kint), intent(in) :: ndof
183  real(kind=kreal), intent(in) :: x(:)
184  real(kind=kreal), intent(out) :: y(:)
185  real(kind=kreal), intent(inout) :: commtime
186 
187  select case(ndof)
188  case (3)
189  call hecmw_ttvec_33(hecmesh, x, y, commtime)
190  case default
191  call hecmw_ttvec_nn(hecmesh, ndof, x, y, commtime)
192  end select
193 
194  end subroutine hecmw_ttvec
195 
196  !C
197  !C***
198  !C*** hecmw_TtmatTvec
199  !C***
200  !C
201  subroutine hecmw_ttmattvec (hecMESH, hecMAT, X, Y, W, COMMtime)
202  use hecmw_util
203  implicit none
204  type (hecmwst_local_mesh), intent(in) :: hecmesh
205  type (hecmwst_matrix), intent(in) :: hecmat
206  real(kind=kreal), intent(in) :: x(:)
207  real(kind=kreal), intent(out) :: y(:), w(:)
208  real(kind=kreal), intent(inout) :: commtime
209 
210  ! call hecmw_Tvec(hecMESH, X, Y, COMMtime)
211  ! call hecmw_matvec(hecMESH, hecMAT, Y, W, COMMtime)
212  ! call hecmw_Ttvec(hecMESH, W, Y, COMMtime)
213  select case(hecmesh%n_dof)
214  case (3)
215  call hecmw_ttmattvec_33 (hecmesh, hecmat, x, y, w, commtime)
216  case default
217  call hecmw_ttmattvec_nn (hecmesh, hecmat, x, y, w, commtime)
218  end select
219 
220  end subroutine hecmw_ttmattvec
221 
222  !C
223  !C***
224  !C*** hecmw_matvec_clear_timer
225  !C***
226  !C
228  implicit none
229  time_ax = 0.d0
230  end subroutine hecmw_matvec_clear_timer
231 
232  !C
233  !C***
234  !C*** hecmw_matvec_get_timer
235  !C***
236  !C
238  implicit none
239  real(kind=kreal) :: hecmw_matvec_get_timer
240  hecmw_matvec_get_timer = time_ax
241  end function hecmw_matvec_get_timer
242 
243  !C
244  !C***
245  !C*** hecmw_mat_diag_sr
246  !C***
247  !C
248  subroutine hecmw_mat_diag_sr(hecMESH, hecMAT, COMMtime)
249  use hecmw_util
252  implicit none
253  type (hecmwst_local_mesh), intent(in) :: hecmesh
254  type (hecmwst_matrix), intent(inout), target :: hecmat
255  real(kind=kreal), intent(inout), optional :: commtime
256 
257  select case(hecmesh%n_dof)
258  case (3)
259  call hecmw_mat_diag_sr_33(hecmesh, hecmat, commtime)
260  case default
261  call hecmw_mat_diag_sr_nn(hecmesh, hecmat, commtime)
262  end select
263 
264  end subroutine hecmw_mat_diag_sr
265 
266  subroutine hecmw_mat_add(hecMAT1, hecMAT2, hecMAT3)
267  use hecmw_util
270  implicit none
271  type (hecmwst_matrix), intent(inout), target :: hecmat1, hecmat2, hecmat3
272 
273  call hecmw_mat_add_nn(hecmat1, hecmat2, hecmat3)
274  end subroutine hecmw_mat_add
275 
276  subroutine hecmw_mat_multiple(hecMAT, alpha)
277  use hecmw_util
280  implicit none
281  type (hecmwst_matrix), intent(inout), target :: hecmat
282  real(kind=kreal), intent(in) :: alpha
283 
284  call hecmw_mat_multiple_nn(hecmat, alpha)
285  end subroutine hecmw_mat_multiple
286 end module hecmw_solver_las
subroutine, public hecmw_matvec_33(hecMESH, hecMAT, X, Y, time_Ax, COMMtime)
subroutine, public hecmw_mat_diag_sr_33(hecMESH, hecMAT, COMMtime)
subroutine, public hecmw_matresid_33(hecMESH, hecMAT, X, B, R, COMMtime)
subroutine, public hecmw_ttmattvec_33(hecMESH, hecMAT, X, Y, W, COMMtime)
subroutine, public hecmw_ttvec_33(hecMESH, X, Y, COMMtime)
subroutine, public hecmw_tvec_33(hecMESH, X, Y, COMMtime)
subroutine, public hecmw_matvec_44(hecMESH, hecMAT, X, Y, time_Ax, COMMtime)
subroutine, public hecmw_matvec_66(hecMESH, hecMAT, X, Y, time_Ax, COMMtime)
subroutine, public hecmw_mat_diag_sr_nn(hecMESH, hecMAT, COMMtime)
subroutine, public hecmw_mat_multiple_nn(hecMAT, alpha)
subroutine, public hecmw_ttmattvec_nn(hecMESH, hecMAT, X, Y, W, COMMtime)
subroutine, public hecmw_matresid_nn(hecMESH, hecMAT, X, B, R, COMMtime)
subroutine, public hecmw_ttvec_nn(hecMESH, ndof, X, Y, COMMtime)
subroutine, public hecmw_tvec_nn(hecMESH, ndof, X, Y, COMMtime)
subroutine, public hecmw_matvec_nn(hecMESH, hecMAT, X, Y, time_Ax, COMMtime)
subroutine, public hecmw_mat_add_nn(hecMAT1, hecMAT2, hecMAT3)
subroutine, public hecmw_mat_add(hecMAT1, hecMAT2, hecMAT3)
subroutine, public hecmw_matvec_clear_timer
subroutine, public hecmw_matvec_set_async(hecMAT)
subroutine, public hecmw_ttmattvec(hecMESH, hecMAT, X, Y, W, COMMtime)
subroutine, public hecmw_ttvec(hecMESH, ndof, X, Y, COMMtime)
real(kind=kreal) function, public hecmw_rel_resid_l2(hecMESH, hecMAT, COMMtime)
subroutine, public hecmw_tvec(hecMESH, ndof, X, Y, COMMtime)
subroutine, public hecmw_mat_multiple(hecMAT, alpha)
subroutine, public hecmw_matresid(hecMESH, hecMAT, X, B, R, COMMtime)
subroutine, public hecmw_mat_diag_sr(hecMESH, hecMAT, COMMtime)
real(kind=kreal) function, public hecmw_matvec_get_timer()
subroutine, public hecmw_matvec(hecMESH, hecMAT, X, Y, COMMtime)
subroutine, public hecmw_matvec_unset_async
subroutine hecmw_innerproduct_r(hecMESH, ndof, X, Y, sum, COMMtime)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal