d725cd968a3307f7387c29419b1a7adb5a22292e
[alexxy/gromacs.git] / src / gromacs / ewald / pme_gpu_timings.cpp
1 /*
2  * This file is part of the GROMACS molecular simulation package.
3  *
4  * Copyright (c) 2016,2017,2018,2019, by the GROMACS development team, led by
5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6  * and including many others, as listed in the AUTHORS file in the
7  * top-level source directory and at http://www.gromacs.org.
8  *
9  * GROMACS is free software; you can redistribute it and/or
10  * modify it under the terms of the GNU Lesser General Public License
11  * as published by the Free Software Foundation; either version 2.1
12  * of the License, or (at your option) any later version.
13  *
14  * GROMACS is distributed in the hope that it will be useful,
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
17  * Lesser General Public License for more details.
18  *
19  * You should have received a copy of the GNU Lesser General Public
20  * License along with GROMACS; if not, see
21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
23  *
24  * If you want to redistribute modifications to GROMACS, please
25  * consider that scientific software is very special. Version
26  * control is crucial - bugs must be traceable. We will be happy to
27  * consider code for inclusion in the official distribution, but
28  * derived work must not be called official GROMACS. Details are found
29  * in the README & COPYING files - if they are missing, get the
30  * official version at http://www.gromacs.org.
31  *
32  * To help us fund GROMACS development, we humbly ask that you cite
33  * the research papers on the package. Check out http://www.gromacs.org.
34  */
35
36 /*! \internal \file
37  *  \brief Implements PME GPU timing events wrappers.
38  *
39  *  \author Aleksei Iupinov <a.yupinov@gmail.com>
40  * \ingroup module_ewald
41  */
42
43 #include "gmxpre.h"
44
45 #include "pme_gpu_timings.h"
46
47 #include "gromacs/utility/gmxassert.h"
48
49 #include "pme_gpu_internal.h"
50 #include "pme_gpu_types_host.h"
51 #include "pme_gpu_types_host_impl.h"
52
53 /*! \brief
54  * Tells if CUDA-based performance tracking is enabled for PME.
55  *
56  * \param[in] pmeGpu         The PME GPU data structure.
57  * \returns                  True if timings are enabled, false otherwise.
58  */
59 inline bool pme_gpu_timings_enabled(const PmeGpu* pmeGpu)
60 {
61     return pmeGpu->archSpecific->useTiming;
62 }
63
64 void pme_gpu_start_timing(const PmeGpu* pmeGpu, size_t PMEStageId)
65 {
66     if (pme_gpu_timings_enabled(pmeGpu))
67     {
68         GMX_ASSERT(PMEStageId < pmeGpu->archSpecific->timingEvents.size(),
69                    "Wrong PME GPU timing event index");
70         pmeGpu->archSpecific->timingEvents[PMEStageId].openTimingRegion(pmeGpu->archSpecific->pmeStream);
71     }
72 }
73
74 CommandEvent* pme_gpu_fetch_timing_event(const PmeGpu* pmeGpu, size_t PMEStageId)
75 {
76     CommandEvent* timingEvent = nullptr;
77     if (pme_gpu_timings_enabled(pmeGpu))
78     {
79         GMX_ASSERT(PMEStageId < pmeGpu->archSpecific->timingEvents.size(),
80                    "Wrong PME GPU timing event index");
81         timingEvent = pmeGpu->archSpecific->timingEvents[PMEStageId].fetchNextEvent();
82     }
83     return timingEvent;
84 }
85
86 void pme_gpu_stop_timing(const PmeGpu* pmeGpu, size_t PMEStageId)
87 {
88     if (pme_gpu_timings_enabled(pmeGpu))
89     {
90         GMX_ASSERT(PMEStageId < pmeGpu->archSpecific->timingEvents.size(),
91                    "Wrong PME GPU timing event index");
92         pmeGpu->archSpecific->timingEvents[PMEStageId].closeTimingRegion(pmeGpu->archSpecific->pmeStream);
93     }
94 }
95
96 void pme_gpu_get_timings(const PmeGpu* pmeGpu, gmx_wallclock_gpu_pme_t* timings)
97 {
98     if (pme_gpu_timings_enabled(pmeGpu))
99     {
100         GMX_RELEASE_ASSERT(timings, "Null GPU timing pointer");
101         for (size_t i = 0; i < pmeGpu->archSpecific->timingEvents.size(); i++)
102         {
103             timings->timing[i].t = pmeGpu->archSpecific->timingEvents[i].getTotalTime();
104             timings->timing[i].c = pmeGpu->archSpecific->timingEvents[i].getCallCount();
105         }
106     }
107 }
108
109 void pme_gpu_update_timings(const PmeGpu* pmeGpu)
110 {
111     if (pme_gpu_timings_enabled(pmeGpu))
112     {
113         for (const size_t& activeTimer : pmeGpu->archSpecific->activeTimers)
114         {
115             pmeGpu->archSpecific->timingEvents[activeTimer].getLastRangeTime();
116         }
117     }
118 }
119
120 void pme_gpu_reinit_timings(const PmeGpu* pmeGpu)
121 {
122     if (pme_gpu_timings_enabled(pmeGpu))
123     {
124         pmeGpu->archSpecific->activeTimers.clear();
125         pmeGpu->archSpecific->activeTimers.insert(gtPME_SPLINEANDSPREAD);
126         // TODO: no separate gtPME_SPLINE and gtPME_SPREAD as they are not used currently
127         if (pme_gpu_performs_FFT(pmeGpu))
128         {
129             pmeGpu->archSpecific->activeTimers.insert(gtPME_FFT_C2R);
130             pmeGpu->archSpecific->activeTimers.insert(gtPME_FFT_R2C);
131         }
132         if (pme_gpu_performs_solve(pmeGpu))
133         {
134             pmeGpu->archSpecific->activeTimers.insert(gtPME_SOLVE);
135         }
136         if (pme_gpu_performs_gather(pmeGpu))
137         {
138             pmeGpu->archSpecific->activeTimers.insert(gtPME_GATHER);
139         }
140     }
141 }
142
143 void pme_gpu_reset_timings(const PmeGpu* pmeGpu)
144 {
145     if (pme_gpu_timings_enabled(pmeGpu))
146     {
147         for (size_t i = 0; i < pmeGpu->archSpecific->timingEvents.size(); i++)
148         {
149             pmeGpu->archSpecific->timingEvents[i].reset();
150         }
151     }
152 }