FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_solver_misc.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 
8 contains
9 
10  !C
11  !C***
12  !C*** hecmw_innerProduct_I
13  !C***
14  !C
15  subroutine hecmw_innerproduct_i (hecMESH, ndof, X, Y, sum, COMMtime )
16  use hecmw_util
17  use m_hecmw_comm_f
18 
19  implicit none
20  type (hecmwST_local_mesh) :: hecMESH
21  integer(kind=kint) :: ndof
22  integer(kind=kint) :: X(:), Y(:)
23  integer(kind=kint) :: sum
24  real(kind=kreal), optional :: commtime
25 
26  integer(kind=kint) :: i
27  real(kind=kreal) :: start_time, end_time
28 
29  sum = 0
30  do i = 1, hecmesh%nn_internal * ndof
31  sum = sum + x(i)*y(i)
32  end do
33 
34  start_time= hecmw_wtime()
35  call hecmw_allreduce_i1 (hecmesh, sum, hecmw_sum)
36  end_time= hecmw_wtime()
37  if (present(commtime)) commtime = commtime + end_time - start_time
38 
39  end subroutine hecmw_innerproduct_i
40 
41  !C
42  !C***
43  !C*** hecmw_innerProduct_R
44  !C***
45  !C
46  subroutine hecmw_innerproduct_r (hecMESH, ndof, X, Y, sum, COMMtime )
47  use hecmw_util
48  use m_hecmw_comm_f
49 
50  implicit none
51  type (hecmwST_local_mesh) :: hecMESH
52  integer(kind=kint) :: ndof
53  real(kind=kreal) :: x(:), y(:)
54  real(kind=kreal) :: sum
55  real(kind=kreal), optional :: commtime
56 
57  integer(kind=kint) :: i
58  real(kind=kreal) :: start_time, end_time
59 
60  sum = 0.0d0
61  do i = 1, hecmesh%nn_internal * ndof
62  sum = sum + x(i)*y(i)
63  end do
64 
65  start_time= hecmw_wtime()
66  call hecmw_allreduce_r1 (hecmesh, sum, hecmw_sum)
67  end_time= hecmw_wtime()
68  if (present(commtime)) commtime = commtime + end_time - start_time
69 
70  end subroutine hecmw_innerproduct_r
71 
72  !C
73  !C***
74  !C*** hecmw_innerProduct_R_nocomm
75  !C***
76  !C
77  subroutine hecmw_innerproduct_r_nocomm (hecMESH, ndof, X, Y, sum)
78  use hecmw_util
79 
80  implicit none
81  type (hecmwST_local_mesh) :: hecMESH
82  integer(kind=kint) :: ndof
83  real(kind=kreal) :: x(:), y(:)
84  real(kind=kreal) :: sum
85 
86  integer(kind=kint) :: i
87 
88  sum = 0.0d0
89  do i = 1, hecmesh%nn_internal * ndof
90  sum = sum + x(i)*y(i)
91  end do
92 
93  end subroutine hecmw_innerproduct_r_nocomm
94 
95  !C
96  !C***
97  !C*** hecmw_axpy_R
98  !C***
99  !C
100  subroutine hecmw_axpy_r (hecMESH, ndof, alpha, X, Y)
101 
102  use hecmw_util
103 
104  implicit none
105  type (hecmwST_local_mesh) :: hecMESH
106  integer(kind=kint) :: ndof
107  real(kind=kreal) :: x(:), y(:)
108  real(kind=kreal) :: alpha
109 
110  integer(kind=kint) :: i
111 
112 
113  !$OMP PARALLEL
114  !$OMP DO
115  do i = 1, hecmesh%nn_internal * ndof
116  y(i) = alpha * x(i) + y(i)
117  end do
118  !$OMP END DO
119  !$OMP END PARALLEL
120 
121  end subroutine hecmw_axpy_r
122 
123  !C
124  !C***
125  !C*** hecmw_xpay_R
126  !C***
127  !C
128  subroutine hecmw_xpay_r (hecMESH, ndof, alpha, X, Y)
129 
130  use hecmw_util
131 
132  implicit none
133  type (hecmwST_local_mesh) :: hecMESH
134  integer(kind=kint) :: ndof
135  real(kind=kreal) :: x(:), y(:)
136  real(kind=kreal) :: alpha
137 
138  integer(kind=kint) :: i
139 
140  !$OMP PARALLEL
141  !$OMP DO
142  do i = 1, hecmesh%nn_internal * ndof
143  y(i) = x(i) + alpha * y(i)
144  end do
145  !$OMP END DO
146  !$OMP END PARALLEL
147 
148  end subroutine hecmw_xpay_r
149 
150  !C
151  !C***
152  !C*** hecmw_axpyz_R
153  !C***
154  !C
155  subroutine hecmw_axpyz_r (hecMESH, ndof, alpha, X, Y, Z)
156 
157  use hecmw_util
158 
159  implicit none
160  type (hecmwST_local_mesh) :: hecMESH
161  integer(kind=kint) :: ndof
162  real(kind=kreal) :: x(:), y(:), z(:)
163  real(kind=kreal) :: alpha
164 
165  integer(kind=kint) :: i
166 
167 
168  !$OMP PARALLEL
169  !$OMP DO
170  do i = 1, hecmesh%nn_internal * ndof
171  z(i) = alpha * x(i) + y(i)
172  end do
173  !$OMP END DO
174  !$OMP END PARALLEL
175  end subroutine hecmw_axpyz_r
176 
177  !C
178  !C***
179  !C*** hecmw_scale_R
180  !C***
181  !C
182  subroutine hecmw_scale_r (hecMESH, ndof, alpha, X)
183 
184  use hecmw_util
185 
186  implicit none
187  type (hecmwST_local_mesh) :: hecMESH
188  integer(kind=kint) :: ndof
189  real(kind=kreal) :: x(:)
190  real(kind=kreal) :: alpha
191 
192  integer(kind=kint) :: i
193 
194  !$OMP PARALLEL
195  !$OMP DO
196  do i = 1, hecmesh%nn_internal * ndof
197  x(i) = alpha * x(i)
198  end do
199  !$OMP END DO
200  !$OMP END PARALLEL
201 
202  end subroutine hecmw_scale_r
203 
204  !C
205  !C***
206  !C*** hecmw_copy_R
207  !C***
208  !C
209  subroutine hecmw_copy_r (hecMESH, ndof, X, Y)
210 
211  use hecmw_util
212 
213  implicit none
214  type (hecmwST_local_mesh) :: hecMESH
215  integer(kind=kint) :: ndof
216  real(kind=kreal) :: x(:), y(:)
217 
218  integer(kind=kint) :: i
219 
220 
221  !$OMP PARALLEL
222  !$OMP DO
223  do i = 1, hecmesh%nn_internal * ndof
224  y(i) = x(i)
225  end do
226  !$OMP END DO
227  !$OMP END PARALLEL
228 
229  end subroutine hecmw_copy_r
230 
231  !C
232  !C***
233  !C*** hecmw_time_statistics
234  !C***
235  !C
236  subroutine hecmw_time_statistics (hecMESH, time, &
237  t_max, t_min, t_avg, t_sd)
238  use hecmw_util
239  use m_hecmw_comm_f
240 
241  implicit none
242  type (hecmwST_local_mesh), intent(in) :: hecMESH
243  real(kind=kreal), intent(in) :: time
244  real(kind=kreal), intent(out) :: t_max
245  real(kind=kreal), intent(out), optional :: t_min, t_avg, t_sd
246  real(kind=kreal) :: t2_avg
247  integer(kind=kint) :: nprocs
248 
249  nprocs = hecmw_comm_get_size()
250 
251  t_max = time
252  call hecmw_allreduce_r1(hecmesh, t_max, hecmw_max)
253 
254  if (.not. present(t_min)) return
255  t_min = time
256  call hecmw_allreduce_r1(hecmesh, t_min, hecmw_min)
257 
258  if (.not. present(t_avg)) return
259  t_avg = time
260  call hecmw_allreduce_r1(hecmesh, t_avg, hecmw_sum)
261  t_avg = t_avg / nprocs
262 
263  if (.not. present(t_sd)) return
264  t2_avg = time*time
265  call hecmw_allreduce_r1(hecmesh, t2_avg, hecmw_sum)
266  t2_avg = t2_avg / nprocs
267 
268  t_sd = dsqrt(t2_avg - t_avg*t_avg)
269  end subroutine hecmw_time_statistics
270 
271 end module hecmw_solver_misc
subroutine hecmw_xpay_r(hecMESH, ndof, alpha, X, Y)
subroutine hecmw_innerproduct_r_nocomm(hecMESH, ndof, X, Y, sum)
subroutine hecmw_innerproduct_r(hecMESH, ndof, X, Y, sum, COMMtime)
subroutine hecmw_axpy_r(hecMESH, ndof, alpha, X, Y)
subroutine hecmw_scale_r(hecMESH, ndof, alpha, X)
subroutine hecmw_time_statistics(hecMESH, time, t_max, t_min, t_avg, t_sd)
subroutine hecmw_axpyz_r(hecMESH, ndof, alpha, X, Y, Z)
subroutine hecmw_innerproduct_i(hecMESH, ndof, X, Y, sum, COMMtime)
subroutine hecmw_copy_r(hecMESH, ndof, X, Y)
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), parameter hecmw_max
integer(kind=4), parameter kreal
integer(kind=kint), parameter hecmw_min
real(kind=kreal) function hecmw_wtime()
subroutine hecmw_allreduce_i1(hecMESH, s, ntag)
subroutine hecmw_allreduce_r1(hecMESH, s, ntag)