1!! Copyright (C) 2005-2006 Heiko Appel, Florian Lorenzen
2!!
3!! This program is free software; you can redistribute it and/or modify
4!! it under the terms of the GNU General Public License as published by
5!! the Free Software Foundation; either version 2, or (at your option)
6!! any later version.
7!!
8!! This program is distributed in the hope that it will be useful,
9!! but WITHOUT ANY WARRANTY; without even the implied warranty of
10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
11!! GNU General Public License for more details.
12!!
13!! You should have received a copy of the GNU General Public License
14!! along with this program; if not, write to the Free Software
15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
16!! 02110-1301, USA.
17!!
18
19! Routines to support MPI debugging.
20
21#include "global.h"
22
23module mpi_debug_oct_m
24  use global_oct_m
25  use messages_oct_m
26  use mpi_oct_m
27
28  implicit none
29
30  private
31
32  public ::               &
33    mpi_debug_statistics, &
34    mpi_debug_in,         &
35    mpi_debug_out
36#if defined(HAVE_MPI)
37  public ::                      &
38    TSD_MPI_Barrier,             &
39    TSZ_MPI_Barrier,             &
40    TSI_MPI_Barrier,             &
41    TSD_MPI_Scatterv,            &
42    TSZ_MPI_Scatterv,            &
43    TSI_MPI_Scatterv,            &
44    TSD_MPI_Gatherv,             &
45    TSZ_MPI_Gatherv,             &
46    TSI_MPI_Gatherv,             &
47    TSD_MPI_Alltoallv,           &
48    TSZ_MPI_Alltoallv,           &
49    TSI_MPI_Alltoallv,           &
50    TSD_MPI_Allgatherv,          &
51    TSZ_MPI_Allgatherv,          &
52    TSI_MPI_Allgatherv,          &
53    TSD_MPI_Bcast,               &
54    TSZ_MPI_Bcast,               &
55    TSI_MPI_Bcast,               &
56    TSD_MPI_Allreduce,           &
57    TSZ_MPI_Allreduce,           &
58    TSI_MPI_Allreduce,           &
59    TSD_MPI_Alltoall,            &
60    TSZ_MPI_Alltoall,            &
61    TSI_MPI_Alltoall,            &
62    TSD_MPI_Allgather,           &
63    TSZ_MPI_Allgather,           &
64    TSI_MPI_Allgather
65#endif
66  integer, parameter :: C_NUM_MPI_ROUTINES = 11
67
68  integer, public, parameter ::  &
69    C_MPI_BARRIER    = 1,        &
70    C_MPI_SCATTERV   = 2,        &
71    C_MPI_GATHERV    = 3,        &
72    C_MPI_ALLTOALLV  = 4,        &
73    C_MPI_ALLGATHERV = 5,        &
74    C_MPI_BCAST      = 6,        &
75    C_MPI_ALLREDUCE  = 7,        &
76    C_MPI_ALLTOALL   = 8,        &
77    C_MPI_ALLGATHER  = 9,        &
78    C_MPI_FILE_READ  = 10,       &
79    C_MPI_FILE_WRITE = 11
80
81  character(len=14), dimension(C_NUM_MPI_ROUTINES), public :: mpi_rlabel = &
82    (/                           &
83    'MPI_BARRIER   ',            &
84    'MPI_SCATTERV  ',            &
85    'MPI_GATHERV   ',            &
86    'MPI_ALLTOALLV ',            &
87    'MPI_ALLGATHERV',            &
88    'MPI_BCAST     ',            &
89    'MPI_ALLREDUCE ',            &
90    'MPI_ALLTOALL  ',            &
91    'MPI_ALLGATHER ',            &
92    'MPI_FILE_READ ',            &
93    'MPI_FILE_WRITE'             &
94    /)
95
96  integer, public :: call_counter(C_NUM_MPI_ROUTINES) = 0
97  real(8), public :: sec_accum(C_NUM_MPI_ROUTINES)    = M_ZERO
98
99  real(8) :: sec_in
100
101contains
102
103  ! ---------------------------------------------------------
104  subroutine mpi_debug_statistics()
105    integer :: j
106    real(8) :: usec_call(C_NUM_MPI_ROUTINES)
107
108    if(.not.debug%info) return
109
110    message(1) = ''
111    message(2) = hyphens
112    message(3) = ''
113    write(message(4), '(23x,a,4x,a,8x,a)') 'total time', 'calls', 'usec/call'
114    do j = 1, C_NUM_MPI_ROUTINES
115      if (call_counter(j) <= 0) then
116        usec_call(j) = 0
117      else
118        usec_call(j) = (sec_accum(j)*1000000)/call_counter(j)
119      end if
120
121      write(message(j+4),'(a,f13.6,6x,i4,6x,f13.0)') &
122        mpi_rlabel(j)//' : ', sec_accum(j),          &
123        call_counter(j), usec_call(j)
124    end do
125    message(C_NUM_MPI_ROUTINES+5) = ''
126    message(C_NUM_MPI_ROUTINES+6) = hyphens
127    call messages_debug(C_NUM_MPI_ROUTINES+6)
128  end subroutine mpi_debug_statistics
129
130
131  ! ---------------------------------------------------------
132  subroutine mpi_debug_in(comm, index)
133    integer, intent(in) :: comm, index
134
135    if(.not.debug%info) return
136
137    call_counter(index) = call_counter(index) + 1
138#if defined(HAVE_MPI)
139    sec_in              = MPI_Wtime()
140#endif
141    write(message(1),'(a,f18.6,a,z8.8,a,i6.6,a,f13.6)') '* MPI_I ', &
142      sec_in, ' '//mpi_rlabel(index)//' : 0x', comm, ' | ',  &
143      call_counter(index), ' - ', sec_accum(index)
144    call messages_debug(1)
145  end subroutine mpi_debug_in
146
147
148  ! ---------------------------------------------------------
149  subroutine mpi_debug_out(comm, index)
150    integer, intent(in) :: comm, index
151
152    real(8) :: sec_out, sec_diff
153
154    if(.not.debug%info) return
155
156#if defined(HAVE_MPI)
157    sec_out = MPI_Wtime()
158#endif
159    call mpi_time_accum(index, sec_out, sec_diff)
160    write(message(1),'(a,f18.6,a,z8.8,a,i6.6,a,f13.6,a,f13.6)')         &
161      '* MPI_O ', sec_out, ' '//mpi_rlabel(index)//' : 0x', comm, ' | ', &
162      call_counter(index), ' - ', sec_accum(index), ' - ', sec_diff
163    call messages_debug(1)
164  end subroutine mpi_debug_out
165
166
167  ! ---------------------------------------------------------
168  subroutine mpi_time_accum(index, sec, sec_diff)
169    integer, intent(in)  :: index
170    real(8), intent(in)  :: sec
171    real(8), intent(out) :: sec_diff
172
173    sec_diff         = sec - sec_in
174    sec_accum(index) = sec_accum(index) + sec_diff
175  end subroutine mpi_time_accum
176
177#if defined(HAVE_MPI)
178#include "undef.F90"
179#include "real.F90"
180#include "mpi_debug_inc.F90"
181
182#include "undef.F90"
183#include "complex.F90"
184#include "mpi_debug_inc.F90"
185
186#include "undef.F90"
187#include "integer.F90"
188#include "mpi_debug_inc.F90"
189#endif
190end module mpi_debug_oct_m
191
192
193!! Local Variables:
194!! mode: f90
195!! coding: utf-8
196!! End:
197