SCALE-RM
scale_prof.F90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
10 #include "scalelib.h"
11 module scale_prof
12  !-----------------------------------------------------------------------------
13  !
14  !++ used modules
15  !
16  use scale_precision
17  use scale_io
18  !-----------------------------------------------------------------------------
19  implicit none
20  private
21  !-----------------------------------------------------------------------------
22  !
23  !++ Public procedure
24  !
25  public :: prof_setup
26  public :: prof_finalize
27  public :: prof_setprefx
28  public :: prof_rapstart
29  public :: prof_rapend
30  public :: prof_rapreport
31 
32 #ifdef PAPI
33  public :: prof_papi_rapstart
34  public :: prof_papi_rapstop
35  public :: prof_papi_rapreport
36 #endif
37 
38  public :: prof_valcheck
39 
40  interface prof_valcheck
41  module procedure prof_valcheck_sp_1d
42  module procedure prof_valcheck_sp_2d
43  module procedure prof_valcheck_sp_3d
44  module procedure prof_valcheck_sp_4d
45  module procedure prof_valcheck_sp_5d
46  module procedure prof_valcheck_sp_6d
47  module procedure prof_valcheck_dp_1d
48  module procedure prof_valcheck_dp_2d
49  module procedure prof_valcheck_dp_3d
50  module procedure prof_valcheck_dp_4d
51  module procedure prof_valcheck_dp_5d
52  module procedure prof_valcheck_dp_6d
53  end interface prof_valcheck
54 
55  !-----------------------------------------------------------------------------
56  !
57  !++ Public parameters & variables
58  !
59  !-----------------------------------------------------------------------------
60  !
61  !++ Private procedure
62  !
63  private :: get_rapid
64 
65  !-----------------------------------------------------------------------------
66  !
67  !++ Private parameters & variables
68  !
69  integer, private, parameter :: PROF_rapnlimit = 300
70  character(len=H_SHORT), private :: PROF_prefix = ''
71  integer, private :: PROF_rapnmax = 0
72  character(len=H_SHORT), private :: PROF_rapname(PROF_rapnlimit)
73  integer, private :: PROF_grpnmax = 0
74  character(len=H_SHORT), private :: PROF_grpname(PROF_rapnlimit)
75  integer, private :: PROF_grpid (PROF_rapnlimit)
76  real(DP), private :: PROF_raptstr(PROF_rapnlimit)
77  real(DP), private :: PROF_rapttot(PROF_rapnlimit)
78  integer, private :: PROF_rapnstr(PROF_rapnlimit)
79  integer, private :: PROF_rapcnt (PROF_rapnlimit)
80  integer, private :: PROF_rapnend(PROF_rapnlimit)
81  integer, private :: PROF_raplevel(PROF_rapnlimit)
82 
83  integer, private, parameter :: PROF_default_rap_level = 2
84  integer, private :: PROF_rap_level = 2
85  logical, private :: PROF_mpi_barrier = .false.
86 
87 #ifdef PAPI
88  integer(DP), private :: PROF_PAPI_flops = 0
89  real(SP), private :: PROF_PAPI_real_time = 0.0
90  real(SP), private :: PROF_PAPI_proc_time = 0.0
91  real(SP), private :: PROF_PAPI_mflops = 0.0
92  integer, private :: PROF_PAPI_check
93 #endif
94 
95  character(len=7), private :: PROF_header
96  character(len=16), private :: PROF_item
97  real(DP), private :: PROF_max
98  real(DP), private :: PROF_min
99  real(DP), private :: PROF_sum
100 
101  logical, private :: PROF_barrier_flag
102 
103  !-----------------------------------------------------------------------------
104 contains
105  !-----------------------------------------------------------------------------
106  subroutine prof_setup
107  use scale_prc, only: &
108  prc_abort
109  implicit none
110 
111  namelist / param_prof / &
112  prof_rap_level, &
113  prof_mpi_barrier
114 
115  integer :: ierr
116 
117  log_newline
118  log_info("PROF_setup",*) 'Setup'
119 
120  !--- read namelist
121  rewind(io_fid_conf)
122  read(io_fid_conf,nml=param_prof,iostat=ierr)
123  if( ierr < 0 ) then !--- missing
124  log_info("PROF_setup",*) 'Not found namelist. Default used.'
125  elseif( ierr > 0 ) then !--- fatal error
126  log_error("PROF_setup",*) 'Not appropriate names in namelist PARAM_PROF. Check!'
127  call prc_abort
128  endif
129  log_nml(param_prof)
130 
131  log_newline
132  log_info("PROF_setup",*) 'Rap output level = ', prof_rap_level
133  log_info("PROF_setup",*) 'Add MPI_barrier in every rap? = ', prof_mpi_barrier
134 
135  prof_prefix = ''
136 
137  return
138  end subroutine prof_setup
139 
140  subroutine prof_finalize
141 
142  prof_rap_level = 2
143  prof_mpi_barrier = .false.
144 #ifdef PAPI
145  prof_papi_flops = 0
146  prof_papi_real_time = 0.0
147  prof_papi_proc_time = 0.0
148  prof_papi_mflops = 0.0
149 #endif
150 
151  prof_rapnmax = 0
152  prof_grpnmax = 0
153 
154  return
155  end subroutine prof_finalize
156 
157  !-----------------------------------------------------------------------------
158  subroutine prof_setprefx( &
159  prefxname )
160  implicit none
161 
162  character(len=*), intent(in) :: prefxname
163 
164  !---------------------------------------------------------------------------
165 
166  prof_prefix = prefxname
167 
168  return
169  end subroutine prof_setprefx
170 
171  !-----------------------------------------------------------------------------
173  subroutine prof_rapstart( rapname_base, level, disable_barrier )
174  use scale_prc, only: &
175  prc_mpibarrier, &
177  implicit none
178 
179  character(len=*), intent(in) :: rapname_base
180  integer, intent(in), optional :: level
181  logical, intent(in), optional :: disable_barrier
182 
183  character(len=H_SHORT) :: rapname
184 
185  integer :: id
186  integer :: level_
187  integer :: tn
188  integer :: i
189  logical :: disable_barrier_
190  !$ integer :: omp_get_thread_num
191  !---------------------------------------------------------------------------
192 
193  tn = 0
194  !$ tn = omp_get_thread_num()
195  if ( tn > 0 ) return
196 
197  if ( present(level) ) then
198  level_ = level
199  else
200  level_ = prof_default_rap_level
201  endif
202 
203  if ( present(disable_barrier) ) then
204  disable_barrier_ = disable_barrier
205  else
206  disable_barrier_ = .false.
207  endif
208 
209  if( level_ > prof_rap_level ) return
210 
211  if ( len_trim(prof_prefix) > 0 ) then
212  rapname = trim(prof_prefix)//" "//trim(rapname_base)
213  else
214  rapname = rapname_base
215  end if
216 
217  id = get_rapid( rapname, level_ )
218 
219  prof_rapcnt(id) = prof_rapcnt(id) + 1
220 
221  if ( prof_rapcnt(id) > 1 ) return
222 
223  if ( ( .not. disable_barrier_ ) .and. prof_mpi_barrier ) call prc_mpibarrier
224 
225  prof_raptstr(id) = prc_mpitime()
226  prof_rapnstr(id) = prof_rapnstr(id) + 1
227 
228  !LOG_INFO("PROF_rapstart",'(1x,A,I8)') rapname, PROF_rapnstr(id)
229  !call flush(IO_FID_LOG)
230 
231 #ifdef FAPP
232  i = index(rapname," ")
233  if ( i == 0 .or. i > len_trim(rapname)) then
234  call fapp_start( rapname, id, level_ )
235  else
236  call fapp_start( rapname(1:i-1)//"_"//trim(rapname(i+1:)), id, level_ )
237  end if
238 #endif
239 
240  return
241  end subroutine prof_rapstart
242 
243  !-----------------------------------------------------------------------------
245  subroutine prof_rapend( rapname_base, level, disable_barrier )
246  use scale_prc, only: &
247  prc_mpibarrier, &
249  implicit none
250 
251  character(len=*), intent(in) :: rapname_base
252  integer, intent(in), optional :: level
253  logical, intent(in), optional :: disable_barrier
254 
255  character(len=H_SHORT) :: rapname
256 
257  integer :: id
258  integer :: level_
259  integer :: tn
260  integer :: i
261  logical :: disable_barrier_
262  !$ integer :: omp_get_thread_num
263  !---------------------------------------------------------------------------
264 
265  tn = 0
266  !$ tn = omp_get_thread_num()
267  if ( tn > 0 ) return
268 
269  if ( present(level) ) then
270  if( level > prof_rap_level ) return
271  endif
272 
273  if ( len_trim(prof_prefix) > 0 ) then
274  rapname = trim(prof_prefix)//" "//trim(rapname_base)
275  else
276  rapname = rapname_base
277  end if
278 
279  if ( present(disable_barrier) ) then
280  disable_barrier_ = disable_barrier
281  else
282  disable_barrier_ = .false.
283  endif
284 
285 
286  level_ = -1
287  id = get_rapid( rapname, level_ )
288 
289  if( level_ > prof_rap_level ) return
290 
291  prof_rapcnt(id) = prof_rapcnt(id) - 1
292 
293  if ( prof_rapcnt(id) > 0 ) return
294 
295 #ifdef FAPP
296  i = index(rapname," ")
297  if ( i == 0 .or. i > len_trim(rapname)) then
298  call fapp_stop( rapname, id, level_ )
299  else
300  call fapp_stop( rapname(1:i-1)//"_"//trim(rapname(i+1:)), id, level_ )
301  end if
302 #endif
303 
304  prof_rapttot(id) = prof_rapttot(id) + ( prc_mpitime()-prof_raptstr(id) )
305  prof_rapnend(id) = prof_rapnend(id) + 1
306 
307  if ( ( .not. disable_barrier_ ) .and. prof_mpi_barrier ) call prc_mpibarrier
308 
309  return
310  end subroutine prof_rapend
311 
312  !-----------------------------------------------------------------------------
314  subroutine prof_rapreport
315  use scale_prc, only: &
316  prc_mpitimestat, &
318  implicit none
319 
320  real(dp) :: avgvar(prof_rapnlimit)
321  real(dp) :: maxvar(prof_rapnlimit)
322  real(dp) :: minvar(prof_rapnlimit)
323  integer :: maxidx(prof_rapnlimit)
324  integer :: minidx(prof_rapnlimit)
325 
326  integer :: idx(prof_rapnlimit)
327  integer :: cnt
328 
329  integer :: id, gid
330  integer :: fid
331  integer :: i, j
332  !---------------------------------------------------------------------------
333 
334  do id = 1, prof_rapnmax
335  if ( prof_rapnstr(id) /= prof_rapnend(id) ) then
336  log_warn("PROF_rapreport",*) 'Mismatch Report',id,prof_rapname(id),prof_rapnstr(id),prof_rapnend(id)
337  endif
338  enddo
339 
340  log_newline
341  log_info("PROF_rapreport",'(1x,A,I2,A)') 'Computational Time Report (Rap level = ', prof_rap_level, ')'
342 
343  if ( io_log_allnode ) then ! report for each node
344 
345  do gid = 1, prof_rapnmax
346  cnt = 0
347  do id = 1, prof_rapnmax
348  if ( prof_raplevel(id) <= prof_rap_level &
349  .AND. prof_grpid(id) == gid ) then
350  cnt = cnt + 1
351  idx(cnt) = id
352  end if
353  end do
354  do j = 1, cnt-1
355  do i = j+1, cnt
356  if ( prof_raplevel(idx(i)) < prof_raplevel(idx(j)) .or. &
357  ( prof_raplevel(idx(i)) == prof_raplevel(idx(j)) &
358  .and. prof_rapname(idx(i)) < prof_rapname(idx(j)) ) ) then
359  id = idx(i)
360  idx(i) = idx(j)
361  idx(j) = id
362  end if
363  end do
364  end do
365  do i = 1, cnt
366  id = idx(i)
367  log_info_cont('(1x,2A,I2,A,F10.3,A,I9)') &
368  prof_rapname(id), ' lev=', prof_raplevel(id), &
369  ': T=',prof_rapttot(id),' N=',prof_rapnstr(id)
370  enddo
371  enddo
372 
373  else
374 
375  call prc_mpitimestat( avgvar(1:prof_rapnmax), &
376  maxvar(1:prof_rapnmax), &
377  minvar(1:prof_rapnmax), &
378  maxidx(1:prof_rapnmax), &
379  minidx(1:prof_rapnmax), &
380  prof_rapttot(1:prof_rapnmax) )
381 
382  fid = -1
383  if ( io_log_suppress ) then ! report to STDOUT
384  if ( prc_ismaster ) then
385  write(*,*) 'INFO [PROF_rapreport] Computational Time Report'
386  fid = io_fid_stdout ! master node
387  endif
388  else
389  if ( io_l ) fid = io_fid_log
390  endif
391 
392  if ( fid > 0 ) then
393  do gid = 1, prof_rapnmax
394  cnt = 0
395  do id = 1, prof_rapnmax
396  if ( prof_raplevel(id) <= prof_rap_level &
397  .AND. prof_grpid(id) == gid ) then
398  cnt = cnt + 1
399  idx(cnt) = id
400  end if
401  end do
402  do j = 1, cnt-1
403  do i = j+1, cnt
404  if ( prof_raplevel(idx(i)) < prof_raplevel(idx(j)) .or. &
405  ( prof_raplevel(idx(i)) == prof_raplevel(idx(j)) &
406  .and. prof_rapname(idx(i)) < prof_rapname(idx(j)) ) ) then
407  id = idx(i)
408  idx(i) = idx(j)
409  idx(j) = id
410  end if
411  end do
412  end do
413  do i = 1, cnt
414  id = idx(i)
415  write(fid,'(1x,2A,I2,A,F10.3,2(A,F10.3,A,I6,A),A,I9)') &
416  prof_rapname(id), ' lev=', prof_raplevel(id), &
417  ': T(avg)=',avgvar(id), &
418  ', T(max)=',maxvar(id),'[',maxidx(id),']', &
419  ', T(min)=',minvar(id),'[',minidx(id),']', &
420  ', N=',prof_rapnstr(id)
421  end do
422  enddo
423  end if
424 
425  endif
426 
427  return
428  end subroutine prof_rapreport
429 
430 #ifdef PAPI
431  !-----------------------------------------------------------------------------
433  subroutine prof_papi_rapstart
434  implicit none
435  !---------------------------------------------------------------------------
436 
437  call papif_flops( prof_papi_real_time, prof_papi_proc_time, prof_papi_flops, prof_papi_mflops, prof_papi_check )
438 
439  return
440  end subroutine prof_papi_rapstart
441 
442  !-----------------------------------------------------------------------------
444  subroutine prof_papi_rapstop
445  implicit none
446  !---------------------------------------------------------------------------
447 
448  call papif_flops( prof_papi_real_time, prof_papi_proc_time, prof_papi_flops, prof_papi_mflops, prof_papi_check )
449 
450  return
451  end subroutine prof_papi_rapstop
452 
453  !-----------------------------------------------------------------------------
455  subroutine prof_papi_rapreport
456  use scale_prc, only: &
457  prc_mpitimestat, &
458  prc_nprocs, &
460  implicit none
461 
462  real(dp) :: avgvar(3)
463  real(dp) :: maxvar(3)
464  real(dp) :: minvar(3)
465  integer :: maxidx(3)
466  integer :: minidx(3)
467 
468  real(dp) :: zerosw
469  real(dp) :: prof_papi_gflop
470  real(dp) :: statistics(3)
471  !---------------------------------------------------------------------------
472 
473  prof_papi_gflop = real(prof_papi_flops,kind=8) / 1024.0_dp**3
474 
475  if ( io_log_allnode ) then ! report for each node
476 
477  log_newline
478  log_info("PROF_PAPI_rapreport",*) 'PAPI Report [Local PE information]'
479  log_info("PROF_PAPI_rapreport",'(1x,A,F15.3)') 'Real time [sec] : ', prof_papi_real_time
480  log_info("PROF_PAPI_rapreport",'(1x,A,F15.3)') 'CPU time [sec] : ', prof_papi_proc_time
481  log_info("PROF_PAPI_rapreport",'(1x,A,F15.3)') 'FLOP [GFLOP] : ', prof_papi_gflop
482  log_info("PROF_PAPI_rapreport",'(1x,A,F15.3)') 'FLOPS by PAPI [GFLOPS] : ', prof_papi_mflops/1024.0_dp
483  log_info("PROF_PAPI_rapreport",'(1x,A,F15.3)') 'FLOP / CPU Time [GFLOPS] : ', prof_papi_gflop/prof_papi_proc_time
484 
485  else
486  statistics(1) = real(prof_papi_real_time,kind=8)
487  statistics(2) = real(prof_papi_proc_time,kind=8)
488  statistics(3) = prof_papi_gflop
489 
490  call prc_mpitimestat( avgvar(1:3), &
491  maxvar(1:3), &
492  minvar(1:3), &
493  maxidx(1:3), &
494  minidx(1:3), &
495  statistics(1:3) )
496 
497  zerosw = 0.5_dp - sign(0.5_dp,maxvar(2)-1.d-12) ! if maxvar(2) = 0 then zerosw = 1
498 
499  log_newline
500  log_info("PROF_PAPI_rapreport",*) 'PAPI Report'
501  log_info("PROF_PAPI_rapreport",'(1x,A,A,F10.3,A,F10.3,A,I6,A,A,F10.3,A,I6,A)') &
502  'Real time [sec]',' T(avg)=',avgvar(1), &
503  ', T(max)=',maxvar(1),'[',maxidx(1),']',', T(min)=',minvar(1),'[',minidx(1),']'
504  log_info("PROF_PAPI_rapreport",'(1x,A,A,F10.3,A,F10.3,A,I6,A,A,F10.3,A,I6,A)') &
505  'CPU time [sec]',' T(avg)=',avgvar(2), &
506  ', T(max)=',maxvar(2),'[',maxidx(2),']',', T(min)=',minvar(2),'[',minidx(2),']'
507  log_info("PROF_PAPI_rapreport",'(1x,A,A,F10.3,A,F10.3,A,I6,A,A,F10.3,A,I6,A)') &
508  'FLOP [GFLOP]',' N(avg)=',avgvar(3), &
509  ', N(max)=',maxvar(3),'[',maxidx(3),']',', N(min)=',minvar(3),'[',minidx(3),']'
510  log_newline
511  log_info("PROF_PAPI_rapreport",'(1x,A,F15.3,A,I6,A)') &
512  'TOTAL FLOP [GFLOP] : ', avgvar(3)*prc_nprocs, '(',prc_nprocs,' PEs)'
513  log_info("PROF_PAPI_rapreport",'(1x,A,F15.3)') &
514  'FLOPS [GFLOPS] : ', avgvar(3)*prc_nprocs * ( 1.0_dp-zerosw ) / ( maxvar(2)+zerosw )
515  log_info("PROF_PAPI_rapreport",'(1x,A,F15.3)') &
516  'FLOPS per PE [GFLOPS] : ', avgvar(3) * ( 1.0_dp-zerosw ) / ( maxvar(2)+zerosw )
517  log_newline
518 
519  if ( io_log_suppress ) then ! report to STDOUT
520  if ( prc_ismaster ) then ! master node
521  write(*,*)
522  write(*,*) '*** PAPI Report'
523  write(*,'(1x,A,A,F10.3,A,F10.3,A,I6,A,A,F10.3,A,I6,A)') &
524  '*** Real time [sec]',' T(avg)=',avgvar(1), &
525  ', T(max)=',maxvar(1),'[',maxidx(1),']',', T(min)=',minvar(1),'[',minidx(1),']'
526  write(*,'(1x,A,A,F10.3,A,F10.3,A,I6,A,A,F10.3,A,I6,A)') &
527  '*** CPU time [sec]',' T(avg)=',avgvar(2), &
528  ', T(max)=',maxvar(2),'[',maxidx(2),']',', T(min)=',minvar(2),'[',minidx(2),']'
529  write(*,'(1x,A,A,F10.3,A,F10.3,A,I6,A,A,F10.3,A,I6,A)') &
530  '*** FLOP [GFLOP]',' N(avg)=',avgvar(3), &
531  ', N(max)=',maxvar(3),'[',maxidx(3),']',', N(min)=',minvar(3),'[',minidx(3),']'
532  write(*,*)
533  write(*,'(1x,A,F15.3,A,I6,A)') &
534  '*** TOTAL FLOP [GFLOP] : ', avgvar(3)*prc_nprocs, '(',prc_nprocs,' PEs)'
535  write(*,'(1x,A,F15.3)') &
536  '*** FLOPS [GFLOPS] : ', avgvar(3)*prc_nprocs * ( 1.0_dp-zerosw ) / ( maxvar(2)+zerosw )
537  write(*,'(1x,A,F15.3)') &
538  '*** FLOPS per PE [GFLOPS] : ', avgvar(3) * ( 1.0_dp-zerosw ) / ( maxvar(2)+zerosw )
539  endif
540  endif
541  endif
542 
543  return
544  end subroutine prof_papi_rapreport
545 #endif
546 
547  !-----------------------------------------------------------------------------
549  function get_rapid( rapname, level ) result(id)
550  use scale_prc, only: &
551  prc_abort
552  implicit none
553 
554  character(len=*), intent(in) :: rapname
555  integer, intent(inout) :: level
556 
557  character (len=H_SHORT) :: trapname
558 
559  integer :: lev
560  integer :: id
561  !---------------------------------------------------------------------------
562 
563  trapname = trim(rapname)
564 
565  do id = 1, prof_rapnmax
566  if ( trapname == prof_rapname(id) ) then
567  lev = prof_raplevel(id)
568 #ifdef QUICKDEBUG
569  if ( level > 0 .and. lev .ne. level ) then
570  log_error("PROF_get_rapid",*) 'level is different ', trim(rapname), lev, level
571  call prc_abort
572  end if
573 #endif
574  level = lev
575  return
576  endif
577  enddo
578 
579  prof_rapnmax = prof_rapnmax + 1
580  id = prof_rapnmax
581  prof_rapname(id) = trapname
582 
583  prof_rapnstr(id) = 0
584  prof_rapnend(id) = 0
585  prof_rapttot(id) = 0.0_dp
586  prof_rapcnt(id) = 0
587 
588  prof_grpid(id) = get_grpid(trapname)
589  prof_raplevel(id) = level
590 
591  return
592  end function get_rapid
593 
594  !-----------------------------------------------------------------------------
596  function get_grpid( rapname ) result(gid)
597  implicit none
598 
599  character(len=*), intent(in) :: rapname
600 
601  character(len=H_SHORT) :: grpname
602 
603  integer :: gid
604  integer :: idx
605  !---------------------------------------------------------------------------
606 
607  idx = index(rapname," ")
608  if ( idx > 1 ) then
609  grpname = rapname(1:idx-1)
610  else
611  grpname = rapname
612  endif
613 
614  do gid = 1, prof_grpnmax
615  if( grpname == prof_grpname(gid) ) return
616  enddo
617 
618  prof_grpnmax = prof_grpnmax + 1
619  gid = prof_grpnmax
620  prof_grpname(gid) = grpname
621 
622  return
623  end function get_grpid
624 
625  !-----------------------------------------------------------------------------
626  subroutine prof_valcheck_sp_1d( &
627  header, &
628  varname, &
629  var )
630  implicit none
631 
632  character(len=*), intent(in) :: header
633  character(len=*), intent(in) :: varname
634  real(sp), intent(in) :: var(:)
635  !---------------------------------------------------------------------------
636 
637  prof_header = trim(header)
638  prof_item = trim(varname)
639  prof_max = real(maxval(var),kind=dp)
640  prof_min = real(minval(var),kind=dp)
641  prof_sum = real(sum(var),kind=dp)
642  log_info("PROF_valcheck_SP_1D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
643  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
644 
645  return
646  end subroutine prof_valcheck_sp_1d
647 
648  !-----------------------------------------------------------------------------
649  subroutine prof_valcheck_sp_2d( &
650  header, &
651  varname, &
652  var )
653  implicit none
654 
655  character(len=*), intent(in) :: header
656  character(len=*), intent(in) :: varname
657  real(sp), intent(in) :: var(:,:)
658  !---------------------------------------------------------------------------
659 
660  prof_header = trim(header)
661  prof_item = trim(varname)
662  prof_max = real(maxval(var),kind=dp)
663  prof_min = real(minval(var),kind=dp)
664  prof_sum = real(sum(var),kind=dp)
665  log_info("PROF_valcheck_SP_2D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
666  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
667 
668  return
669  end subroutine prof_valcheck_sp_2d
670 
671  !-----------------------------------------------------------------------------
672  subroutine prof_valcheck_sp_3d( &
673  header, &
674  varname, &
675  var )
676  implicit none
677 
678  character(len=*), intent(in) :: header
679  character(len=*), intent(in) :: varname
680  real(sp), intent(in) :: var(:,:,:)
681  !---------------------------------------------------------------------------
682 
683  prof_header = trim(header)
684  prof_item = trim(varname)
685  prof_max = real(maxval(var),kind=dp)
686  prof_min = real(minval(var),kind=dp)
687  prof_sum = real(sum(var),kind=dp)
688  log_info("PROF_valcheck_SP_3D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
689  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
690 
691  return
692  end subroutine prof_valcheck_sp_3d
693 
694  !-----------------------------------------------------------------------------
695  subroutine prof_valcheck_sp_4d( &
696  header, &
697  varname, &
698  var )
699  implicit none
700 
701  character(len=*), intent(in) :: header
702  character(len=*), intent(in) :: varname
703  real(sp), intent(in) :: var(:,:,:,:)
704  !---------------------------------------------------------------------------
705 
706  prof_header = trim(header)
707  prof_item = trim(varname)
708  prof_max = real(maxval(var),kind=dp)
709  prof_min = real(minval(var),kind=dp)
710  prof_sum = real(sum(var),kind=dp)
711  log_info("PROF_valcheck_SP_4D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
712  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
713 
714  return
715  end subroutine prof_valcheck_sp_4d
716 
717  !-----------------------------------------------------------------------------
718  subroutine prof_valcheck_sp_5d( &
719  header, &
720  varname, &
721  var )
722  implicit none
723 
724  character(len=*), intent(in) :: header
725  character(len=*), intent(in) :: varname
726  real(sp), intent(in) :: var(:,:,:,:,:)
727  !---------------------------------------------------------------------------
728 
729  prof_header = trim(header)
730  prof_item = trim(varname)
731  prof_max = real(maxval(var),kind=dp)
732  prof_min = real(minval(var),kind=dp)
733  prof_sum = real(sum(var),kind=dp)
734  log_info("PROF_valcheck_SP_5D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
735  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
736 
737  return
738  end subroutine prof_valcheck_sp_5d
739 
740  !-----------------------------------------------------------------------------
741  subroutine prof_valcheck_sp_6d( &
742  header, &
743  varname, &
744  var )
745  implicit none
746 
747  character(len=*), intent(in) :: header
748  character(len=*), intent(in) :: varname
749  real(sp), intent(in) :: var(:,:,:,:,:,:)
750  !---------------------------------------------------------------------------
751 
752  prof_header = trim(header)
753  prof_item = trim(varname)
754  prof_max = real(maxval(var),kind=dp)
755  prof_min = real(minval(var),kind=dp)
756  prof_sum = real(sum(var),kind=dp)
757  log_info("PROF_valcheck_SP_6D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
758  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
759 
760  return
761  end subroutine prof_valcheck_sp_6d
762 
763  !-----------------------------------------------------------------------------
764  subroutine prof_valcheck_dp_1d( &
765  header, &
766  varname, &
767  var )
768  implicit none
769 
770  character(len=*), intent(in) :: header
771  character(len=*), intent(in) :: varname
772  real(dp), intent(in) :: var(:)
773  !---------------------------------------------------------------------------
774 
775  prof_header = trim(header)
776  prof_item = trim(varname)
777  prof_max = real(maxval(var),kind=dp)
778  prof_min = real(minval(var),kind=dp)
779  prof_sum = real(sum(var),kind=dp)
780  log_info("PROF_valcheck_DP_1D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
781  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
782 
783  return
784  end subroutine prof_valcheck_dp_1d
785 
786  !-----------------------------------------------------------------------------
787  subroutine prof_valcheck_dp_2d( &
788  header, &
789  varname, &
790  var )
791  implicit none
792 
793  character(len=*), intent(in) :: header
794  character(len=*), intent(in) :: varname
795  real(dp), intent(in) :: var(:,:)
796  !---------------------------------------------------------------------------
797 
798  prof_header = trim(header)
799  prof_item = trim(varname)
800  prof_max = real(maxval(var),kind=dp)
801  prof_min = real(minval(var),kind=dp)
802  prof_sum = real(sum(var),kind=dp)
803  log_info("PROF_valcheck_DP_2D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
804  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
805 
806  return
807  end subroutine prof_valcheck_dp_2d
808 
809  !-----------------------------------------------------------------------------
810  subroutine prof_valcheck_dp_3d( &
811  header, &
812  varname, &
813  var )
814  implicit none
815 
816  character(len=*), intent(in) :: header
817  character(len=*), intent(in) :: varname
818  real(dp), intent(in) :: var(:,:,:)
819  !---------------------------------------------------------------------------
820 
821  prof_header = trim(header)
822  prof_item = trim(varname)
823  prof_max = real(maxval(var),kind=dp)
824  prof_min = real(minval(var),kind=dp)
825  prof_sum = real(sum(var),kind=dp)
826  log_info("PROF_valcheck_DP_3D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
827  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
828 
829  return
830  end subroutine prof_valcheck_dp_3d
831 
832  !-----------------------------------------------------------------------------
833  subroutine prof_valcheck_dp_4d( &
834  header, &
835  varname, &
836  var )
837  implicit none
838 
839  character(len=*), intent(in) :: header
840  character(len=*), intent(in) :: varname
841  real(dp), intent(in) :: var(:,:,:,:)
842  !---------------------------------------------------------------------------
843 
844  prof_header = trim(header)
845  prof_item = trim(varname)
846  prof_max = real(maxval(var),kind=dp)
847  prof_min = real(minval(var),kind=dp)
848  prof_sum = real(sum(var),kind=dp)
849  log_info("PROF_valcheck_DP_4D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
850  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
851 
852  return
853  end subroutine prof_valcheck_dp_4d
854 
855  !-----------------------------------------------------------------------------
856  subroutine prof_valcheck_dp_5d( &
857  header, &
858  varname, &
859  var )
860  implicit none
861 
862  character(len=*), intent(in) :: header
863  character(len=*), intent(in) :: varname
864  real(dp), intent(in) :: var(:,:,:,:,:)
865  !---------------------------------------------------------------------------
866 
867  prof_header = trim(header)
868  prof_item = trim(varname)
869  prof_max = real(maxval(var),kind=dp)
870  prof_min = real(minval(var),kind=dp)
871  prof_sum = real(sum(var),kind=dp)
872  log_info("PROF_valcheck_DP_5D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
873  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
874 
875  return
876  end subroutine prof_valcheck_dp_5d
877 
878  !-----------------------------------------------------------------------------
879  subroutine prof_valcheck_dp_6d( &
880  header, &
881  varname, &
882  var )
883  implicit none
884 
885  character(len=*), intent(in) :: header
886  character(len=*), intent(in) :: varname
887  real(dp), intent(in) :: var(:,:,:,:,:,:)
888  !---------------------------------------------------------------------------
889 
890  prof_header = trim(header)
891  prof_item = trim(varname)
892  prof_max = real(maxval(var),kind=dp)
893  prof_min = real(minval(var),kind=dp)
894  prof_sum = real(sum(var),kind=dp)
895  log_info("PROF_valcheck_DP_6D",'(1x,A,A7,A,A16,3(A,ES24.16))') &
896  '+',prof_header,'[',prof_item,'] max=',prof_max,',min=',prof_min,',sum=',prof_sum
897 
898  return
899  end subroutine prof_valcheck_dp_6d
900 
901 end module scale_prof
scale_prof::prof_finalize
subroutine, public prof_finalize
Definition: scale_prof.F90:141
scale_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:349
scale_precision
module PRECISION
Definition: scale_precision.F90:14
scale_prof::prof_rapstart
subroutine, public prof_rapstart(rapname_base, level, disable_barrier)
Start raptime.
Definition: scale_prof.F90:174
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_prof::get_grpid
integer function get_grpid(rapname)
Get group ID.
Definition: scale_prof.F90:597
scale_prc::prc_mpibarrier
subroutine, public prc_mpibarrier
Barrier MPI.
Definition: scale_prc.F90:827
scale_io
module STDIO
Definition: scale_io.F90:10
scale_prof::prof_setprefx
subroutine, public prof_setprefx(prefxname)
Definition: scale_prof.F90:160
scale_io::io_fid_log
integer, public io_fid_log
Log file ID.
Definition: scale_io.F90:58
scale_prof
module profiler
Definition: scale_prof.F90:11
scale_prc::prc_nprocs
integer, public prc_nprocs
myrank in local communicator
Definition: scale_prc.F90:89
scale_io::io_log_suppress
logical, public io_log_suppress
suppress all of log output?
Definition: scale_io.F90:65
scale_prof::prof_setup
subroutine, public prof_setup
Definition: scale_prof.F90:107
scale_io::io_l
logical, public io_l
output log or not? (this process)
Definition: scale_io.F90:63
scale_prof::prof_rapend
subroutine, public prof_rapend(rapname_base, level, disable_barrier)
Save raptime.
Definition: scale_prof.F90:246
scale_prc::prc_mpitimestat
subroutine, public prc_mpitimestat(avgvar, maxvar, minvar, maxidx, minidx, var)
Calc global statistics for timer.
Definition: scale_prc.F90:864
scale_io::io_fid_stdout
integer, parameter, public io_fid_stdout
Definition: scale_io.F90:56
scale_io::io_log_allnode
logical, public io_log_allnode
output log for each node?
Definition: scale_io.F90:67
scale_prc::prc_mpitime
real(dp) function, public prc_mpitime()
Get MPI time.
Definition: scale_prc.F90:842
scale_io::io_fid_conf
integer, public io_fid_conf
Config file ID.
Definition: scale_io.F90:57
scale_prc::prc_ismaster
logical, public prc_ismaster
master process in local communicator?
Definition: scale_prc.F90:91
scale_prof::prof_rapreport
subroutine, public prof_rapreport
Report raptime.
Definition: scale_prof.F90:315