Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
statistics.f90
Go to the documentation of this file.
1! Copyright (c) 2021, The Neko Authors
2! All rights reserved.
3!
4! Redistribution and use in source and binary forms, with or without
5! modification, are permitted provided that the following conditions
6! are met:
7!
8! * Redistributions of source code must retain the above copyright
9! notice, this list of conditions and the following disclaimer.
10!
11! * Redistributions in binary form must reproduce the above
12! copyright notice, this list of conditions and the following
13! disclaimer in the documentation and/or other materials provided
14! with the distribution.
15!
16! * Neither the name of the authors nor the names of its
17! contributors may be used to endorse or promote products derived
18! from this software without specific prior written permission.
19!
20! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31! POSSIBILITY OF SUCH DAMAGE.
32!
34module stats
35 use num_types, only : rp, dp
36 use stats_quant, only : stats_quant_t
37 use logger, only : log_size, neko_log
38 use mpi_f08, only : mpi_barrier, mpi_wtime
39 use comm, only : neko_comm
40 implicit none
41 private
42
44 type, private :: quantp_t
45 class(stats_quant_t), pointer :: quantp
46 end type quantp_t
47
49 type, public :: stats_t
50 type(quantp_t), allocatable :: quant_list(:)
51 integer :: n
52 integer :: size
53 real(kind=rp) :: t_begin
54 real(kind=rp) :: t_diff
55 integer :: samp_interval
56 contains
57 procedure, pass(this) :: init => stats_init
58 procedure, pass(this) :: free => stats_free
59 procedure, pass(this) :: add => stats_add
60 procedure, pass(this) :: eval => stats_eval
61 end type stats_t
62
63contains
64
66 subroutine stats_init(this, T_begin, samp_interval, size)
67 class(stats_t), intent(inout) :: this
68 real(kind=rp), intent(in) :: t_begin
69 integer, intent(in) :: samp_interval
70 integer, intent(inout), optional ::size
71 integer :: n, i
72
73 call this%free()
74
75 if (present(size)) then
76 n = size
77 else
78 n = 1
79 end if
80
81 allocate(this%quant_list(n))
82
83 do i = 1, n
84 this%quant_list(i)%quantp => null()
85 end do
86
87 this%n = 0
88 this%size = n
89 this%T_begin = t_begin
90 this%samp_interval = samp_interval
91 this%t_diff = 0.0
92
93 end subroutine stats_init
94
96 subroutine stats_free(this)
97 class(stats_t), intent(inout) :: this
98
99 if (allocated(this%quant_list)) then
100 deallocate(this%quant_list)
101 end if
102
103 this%n = 0
104 this%size = 0
105 end subroutine stats_free
106
108 subroutine stats_add(this, quant)
109 class(stats_t), intent(inout) :: this
110 class(stats_quant_t), intent(inout), target :: quant
111 type(quantp_t), allocatable :: tmp(:)
112
113 if (this%n .ge. this%size) then
114 allocate(tmp(this%size * 2))
115 tmp(1:this%size) = this%quant_list
116 call move_alloc(tmp, this%quant_list)
117 this%size = this%size * 2
118 end if
119
120 this%n = this%n + 1
121 this%quant_list(this%n)%quantp => quant
122 end subroutine stats_add
123
125 subroutine stats_eval(this, t, dt, tstep)
126 class(stats_t), intent(inout) :: this
127 real(kind=rp), intent(in) :: t
128 real(kind=rp), intent(in) :: dt
129 integer, intent(in) :: tstep
130 integer :: i, ierr
131 character(len=LOG_SIZE) :: log_buf
132 real(kind=rp) :: sample_start_time, sample_end_time
133 real(kind=dp) :: sample_time
134
135 if (t .ge. this%T_begin .and. this%n .gt. 0) then
136 this%t_diff = this%t_diff + dt
137 ! There is technically an issue here for the last sample if we
138 ! reset the stats If the reset is not on a multiple of
139 ! samp_interval the weight of the last sample is wrong.
140 if (mod(tstep,this%samp_interval) .eq. 0) then
141 call neko_log%section('Statistics')
142 call mpi_barrier(neko_comm, ierr)
143 sample_start_time = mpi_wtime()
144 do i = 1, this%n
145 call this%quant_list(i)%quantp%update(this%t_diff)
146 end do
147 this%t_diff = 0.0_rp
148 call mpi_barrier(neko_comm, ierr)
149 sample_end_time = mpi_wtime()
150 sample_time = sample_end_time - sample_start_time
151 write(log_buf,'(A17,1x,F10.6,A,F9.6)') 'Sampling at time:', t, &
152 ' Sampling time (s): ', sample_time
153 call neko_log%message(log_buf)
154 call neko_log%end_section()
155 end if
156 end if
157 end subroutine stats_eval
158
159end module stats
Definition comm.F90:1
type(mpi_comm), public neko_comm
MPI communicator.
Definition comm.F90:42
Logging routines.
Definition log.f90:34
type(log_t), public neko_log
Global log stream.
Definition log.f90:70
integer, parameter, public log_size
Definition log.f90:46
integer, parameter, public dp
Definition num_types.f90:9
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Defines a statistical quantity.
Defines a container for all statistics.
subroutine stats_eval(this, t, dt, tstep)
Evaluated all statistical quantities.
subroutine stats_free(this)
Deallocate.
subroutine stats_init(this, t_begin, samp_interval, size)
Initialize statistics, computed after T_begin.
subroutine stats_add(this, quant)
Add a statistic quantitiy quant to the backend.
Pointer to an arbitrary quantitiy.
Statistics backend.
Abstract type defining a statistical quantity.