SCALE-RM
Data Types | Functions/Subroutines | Variables
scale_statistics Module Reference

module Statistics More...

Functions/Subroutines

subroutine, public statistics_setup
 Setup. More...
 
subroutine statistics_total_2d (IA, IS, IE, JA, JS, JE, var, varname, area, total, log_suppress, global, mean, sum)
 Calc domain sum and area-weighted mean. More...
 
subroutine statistics_total_3d (KA, KS, KE, IA, IS, IE, JA, JS, JE, var, varname, vol, total, log_suppress, global, mean, sum)
 Calc domain sum and volume-weighted mean. More...
 
subroutine statistics_horizontal_mean_2d (IA, IS, IE, JA, JS, JE, var, area, varmean)
 Calc horizontal mean value. More...
 
subroutine statistics_horizontal_mean_3d (KA, KS, KE, IA, IS, IE, JA, JS, JE, var, area, varmean)
 
subroutine statistics_horizontal_min_2d (IA, IS, IE, JA, JS, JE, var, varmin)
 Calc horizontal minimum value. More...
 
subroutine statistics_horizontal_min_3d (KA, KS, KE, IA, IS, IE, JA, JS, JE, var, varmin)
 
subroutine statistics_horizontal_max_2d (IA, IS, IE, JA, JS, JE, var, varmax)
 Calc horizontal maximum value. More...
 
subroutine statistics_horizontal_max_3d (KA, KS, KE, IA, IS, IE, JA, JS, JE, var, varmax)
 
subroutine statistics_detail_3d (KA, KS, KE, IA, IS, IE, JA, JS, JE, VA, varname, var, local)
 Search global maximum & minimum value. More...
 
subroutine statistics_detail_2d (IA, IS, IE, JA, JS, JE, VA, varname, var, local)
 
real(rp) function statistics_summation_1d (IA, ARRAY, UNDEF)
 
real(rp) function, public statistics_covariance (IA1, IA2, ARRAY1, ARRAY2, UNDEF)
 
real(rp) function, public statistics_correlation (IA1, IA2, ARRAY1, ARRAY2, UNDEF)
 
real function, public statistics_regression (IA1, IA2, ARRAY1, ARRAY2, UNDEF)
 
real(rp) function, public statistics_lag_correlation (IA1, IA2, ARRAY1, ARRAY2, LAG, UNDEF)
 
real(rp) function, public statistics_partial_correlation (IA1, IA2, IA3, ARRAY1, ARRAY2, ARRAY3, UNDEF)
 
real(rp) function, dimension(ia), public statistics_undef_replace (IA, ARRAY, UNDEF1, UNDEF2)
 
real(rp) function, dimension(ia1), public statistics_undef_embed (IA1, IA2, ARRAY1, ARRAY2, UNDEF)
 
integer function, public statistics_undef_arraysize (IA1, IA2, ARRAY1, ARRAY2, UNDEF)
 

Variables

logical, public statistics_checktotal = .false.
 calc&report variable totals to logfile? More...
 

Detailed Description

module Statistics

Description
global statistics module
Author
Team SCALE
NAMELIST
  • PARAM_STATISTICS
    nametypedefault valuecomment
    STATISTICS_CHECKTOTAL logical .false. calc&report variable totals to logfile?
    STATISTICS_USE_GLOBALCOMM logical .false. calculate total with global communication?

History Output
No history output

Function/Subroutine Documentation

◆ statistics_setup()

subroutine, public scale_statistics::statistics_setup

Setup.

Definition at line 126 of file scale_statistics.F90.

126  use scale_prc, only: &
127  prc_abort
128  use scale_comm_cartesc, only: &
129  comm_setup
130  implicit none
131 
132  namelist / param_statistics / &
133  statistics_checktotal, &
134  statistics_use_globalcomm
135 
136  integer :: ierr
137  !---------------------------------------------------------------------------
138 
139  call comm_setup
140 
141  log_newline
142  log_info("STATISTICS_setup",*) 'Setup'
143 
144  !--- read namelist
145  rewind(io_fid_conf)
146  read(io_fid_conf,nml=param_statistics,iostat=ierr)
147  if( ierr < 0 ) then !--- missing
148  log_info("STATISTICS_setup",*) 'Not found namelist. Default used.'
149  elseif( ierr > 0 ) then !--- fatal error
150  log_error("STATISTICS_setup",*) 'Not appropriate names in namelist PARAM_STATISTICS. Check!'
151  call prc_abort
152  endif
153  log_nml(param_statistics)
154 
155  log_newline
156  log_info("STATISTICS_setup",*) 'Caluculate total statistics for monitoring? : ', statistics_checktotal
157  if ( statistics_use_globalcomm ) then
158  log_info_cont(*) '=> The total is calculated for the global domain.'
159  else
160  log_info_cont(*) '=> The total is calculated for the local domain.'
161  endif
162 
163  return

References scale_comm_cartesc::comm_setup(), scale_io::io_fid_conf, scale_prc::prc_abort(), and statistics_checktotal.

Referenced by mod_rm_driver::rm_driver(), and mod_rm_prep::rm_prep().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ statistics_total_2d()

subroutine scale_statistics::statistics_total_2d ( integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
real(rp), dimension(ia,ja), intent(in)  var,
character(len=*), intent(in)  varname,
real(rp), dimension(ia,ja), intent(in)  area,
real(rp), intent(in)  total,
logical, intent(in), optional  log_suppress,
logical, intent(in), optional  global,
real(rp), intent(out), optional  mean,
real(dp), intent(out), optional  sum 
)

Calc domain sum and area-weighted mean.

Parameters
[in]var3D value
[in]varnamename of item
[in]areaarea of the grid cell
[in]totaltotal area
[in]log_suppresssuppress log output
[in]globalglobal or local sum
[out]meanarea-weighted mean
[out]sumdomain sum

Definition at line 175 of file scale_statistics.F90.

175  use scale_prc, only: &
176  prc_myrank, &
177  prc_abort
178  use scale_const, only: &
179  eps => const_eps, &
180  undef => const_undef
181  implicit none
182 
183  integer, intent(in) :: IA, IS, IE
184  integer, intent(in) :: JA, JS, JE
185 
186  real(RP), intent(in) :: var(IA,JA)
187  character(len=*), intent(in) :: varname
188  real(RP), intent(in) :: area(IA,JA)
189  real(RP), intent(in) :: total
190 
191  logical, intent(in), optional :: log_suppress
192  logical, intent(in), optional :: global
193  real(RP), intent(out), optional :: mean
194  real(DP), intent(out), optional :: sum
195 
196  real(DP) :: statval
197  real(DP) :: sendbuf(2), recvbuf(2)
198  real(DP) :: sum_, mean_
199 
200  logical :: suppress_, global_
201  integer :: ierr
202  integer :: i, j
203  !---------------------------------------------------------------------------
204 
205  statval = 0.0_dp
206  !$acc update host(var(IS,JS)) if(acc_is_present(var))
207  if ( var(is,js) /= undef ) then
208  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2) reduction(+:statval)
209  !$acc kernels copyin(var, area) if(acc_is_present(var))
210  !$acc loop reduction(+:statval)
211  do j = js, je
212  !$acc loop reduction(+:statval)
213  do i = is, ie
214  statval = statval + var(i,j) * area(i,j)
215  end do
216  end do
217  !$acc end kernels
218  end if
219 
220  if ( .NOT. ( statval > -1.0_dp .OR. statval < 1.0_dp ) ) then ! must be NaN
221  log_error("STATISTICS_total_2D",*) 'NaN is detected for ', trim(varname), ' in rank ', prc_myrank
222  call prc_abort
223  endif
224 
225  if ( present(log_suppress) ) then
226  suppress_ = log_suppress
227  else
228  suppress_ = .false.
229  end if
230 
231  if ( present(global) ) then
232  global_ = global
233  else
234  global_ = statistics_use_globalcomm
235  end if
236 
237  if ( global_ ) then
238  call prof_rapstart('COMM_Allreduce', 2)
239  sendbuf(1) = statval
240  sendbuf(2) = total
241  ! All reduce
242  call mpi_allreduce( sendbuf(:), recvbuf(:), &
243  2, &
244  mpi_double_precision, &
245  mpi_sum, &
246  prc_local_comm_world, &
247  ierr )
248  call prof_rapend ('COMM_Allreduce', 2)
249 
250  if ( recvbuf(2) < eps ) then
251  sum_ = undef
252  mean_ = undef
253  else
254  sum_ = recvbuf(1)
255  mean_ = recvbuf(1) / recvbuf(2)
256  end if
257  ! statistics over the all node
258  if ( .not. suppress_ ) then ! if varname is empty, suppress output
259  log_info("STATISTICS_total_2D",'(1x,A,A24,A,ES24.17)') &
260  '[', trim(varname), '] MEAN(global) = ', mean_
261  endif
262  else
263  if ( total < eps ) then
264  sum_ = undef
265  mean_ = undef
266  else
267  sum_ = statval
268  mean_ = statval / total
269  end if
270 
271  ! statistics on each node
272  if ( .not. suppress_ ) then ! if varname is empty, suppress output
273  log_info("STATISTICS_total_2D",'(1x,A,A24,A,ES24.17)') &
274  '[', trim(varname), '] MEAN(local) = ', mean_
275  endif
276  endif
277 
278  if ( present(mean) ) mean = mean_
279  if ( present(sum ) ) sum = sum_
280 
281  return

References scale_const::const_eps, scale_const::const_undef, scale_prc::prc_abort(), scale_prc::prc_myrank, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_total_3d()

subroutine scale_statistics::statistics_total_3d ( integer, intent(in)  KA,
integer, intent(in)  KS,
integer, intent(in)  KE,
integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
real(rp), dimension(ka,ia,ja), intent(in)  var,
character(len=*), intent(in)  varname,
real(rp), dimension(ka,ia,ja), intent(in)  vol,
real(rp), intent(in)  total,
logical, intent(in), optional  log_suppress,
logical, intent(in), optional  global,
real(rp), intent(out), optional  mean,
real(dp), intent(out), optional  sum 
)

Calc domain sum and volume-weighted mean.

Parameters
[in]var3D value
[in]varnamename of item
[in]volvolume of the grid cell
[in]totaltotal volume
[in]log_suppresssuppress log output
[in]globalglobal or local sum
[out]meanvolume/area-weighted total
[out]sumdomain sum

Definition at line 293 of file scale_statistics.F90.

293  use scale_prc, only: &
294  prc_myrank, &
295  prc_abort
296  use scale_const, only: &
297  eps => const_eps , &
298  undef => const_undef
299  implicit none
300 
301  integer, intent(in) :: KA, KS, KE
302  integer, intent(in) :: IA, IS, IE
303  integer, intent(in) :: JA, JS, JE
304 
305  real(RP), intent(in) :: var(KA,IA,JA)
306  character(len=*), intent(in) :: varname
307  real(RP), intent(in) :: vol(KA,IA,JA)
308  real(RP), intent(in) :: total
309 
310  logical, intent(in), optional :: log_suppress
311  logical, intent(in), optional :: global
312  real(RP), intent(out), optional :: mean
313  real(DP), intent(out), optional :: sum
314 
315  real(DP) :: statval
316  real(DP) :: sendbuf(2), recvbuf(2)
317  real(DP) :: mean_, sum_
318 
319  logical :: suppress_, global_
320 
321  real(DP) :: work
322  integer :: ierr
323  integer :: k, i, j
324  !---------------------------------------------------------------------------
325 
326  statval = 0.0_dp
327  !$acc update host(var(KS,IS,JS)) if(acc_is_present(var))
328  if ( var(ks,is,js) /= undef ) then
329  !$omp parallel do OMP_SCHEDULE_ reduction(+:statval) &
330  !$omp private(work)
331  !$acc kernels copyin(var, vol) if(acc_is_present(var))
332  !$acc loop reduction(statval)
333  do j = je, js, -1
334  !$acc loop reduction(statval)
335  do i = ie, is, -1
336 #ifdef _OPENACC
337  !$acc loop reduction(statval)
338  do k = ke, ks, -1
339  statval = statval + var(k,i,j) * vol(k,i,j)
340  enddo
341 #else
342  work = 0.0_rp
343  do k = ke, ks, -1
344  work = work + var(k,i,j) * vol(k,i,j)
345  enddo
346  statval = statval + work
347 #endif
348  enddo
349  enddo
350  !$acc end kernels
351  end if
352 
353  if ( .NOT. ( statval > -1.0_dp .OR. statval < 1.0_dp ) ) then ! must be NaN
354  log_error("STATISTICS_total_3D",*) 'NaN is detected for ', trim(varname), ' in rank ', prc_myrank
355  call prc_abort
356  endif
357 
358  if ( present(log_suppress) ) then
359  suppress_ = log_suppress
360  else
361  suppress_ = .false.
362  end if
363 
364  if ( present(global) ) then
365  global_ = global
366  else
367  global_ = statistics_use_globalcomm
368  end if
369 
370  if ( global_ ) then
371  call prof_rapstart('COMM_Allreduce', 2)
372  sendbuf(1) = statval
373  sendbuf(2) = total
374  ! All reduce
375  call mpi_allreduce( sendbuf(:), recvbuf(:), &
376  2, &
377  mpi_double_precision, &
378  mpi_sum, &
379  prc_local_comm_world, &
380  ierr )
381  call prof_rapend ('COMM_Allreduce', 2)
382 
383  if ( recvbuf(2) < eps ) then
384  sum_ = undef
385  mean_ = undef
386  else
387  sum_ = recvbuf(1)
388  mean_ = recvbuf(1) / recvbuf(2)
389  end if
390  ! statistics over the all node
391  if ( .not. suppress_ ) then ! if varname is empty, suppress output
392  log_info("STATISTICS_total_3D",'(1x,A,A24,A,ES24.17)') &
393  '[', trim(varname), '] MEAN(global) = ', mean_
394  endif
395  else
396  if ( total < eps ) then
397  sum_ = undef
398  mean_ = undef
399  else
400  sum_ = statval
401  mean_ = statval / total
402  end if
403 
404  ! statistics on each node
405  if ( .not. suppress_ ) then ! if varname is empty, suppress output
406  log_info("STATISTICS_total_3D",'(1x,A,A24,A,ES24.17)') &
407  '[', trim(varname), '] MEAN(local) = ', mean_
408  endif
409  endif
410 
411  if ( present(mean) ) mean = mean_
412  if ( present(sum ) ) sum = sum_
413 
414  return

References scale_const::const_eps, scale_const::const_undef, scale_prc::prc_abort(), scale_prc::prc_myrank, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_horizontal_mean_2d()

subroutine scale_statistics::statistics_horizontal_mean_2d ( integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
real(rp), dimension (ia,ja), intent(in)  var,
real(rp), dimension(ia,ja), intent(in)  area,
real(rp), intent(out)  varmean 
)

Calc horizontal mean value.

Definition at line 423 of file scale_statistics.F90.

423  use scale_const, only: &
424  undef => const_undef
425  integer, intent(in) :: IA, IS, IE
426  integer, intent(in) :: JA, JS, JE
427  real(RP), intent(in) :: var (IA,JA)
428  real(RP), intent(in) :: area(IA,JA)
429  real(RP), intent(out) :: varmean
430 
431  real(DP) :: statval (2)
432  real(DP) :: allstatval(2)
433  real(DP) :: s1, s2
434 
435  integer :: ierr
436  integer :: i, j
437  !---------------------------------------------------------------------------
438 
439  s1 = 0.0_dp
440  s2 = 0.0_dp
441  !$omp parallel do reduction(+:s1,s2)
442  !$acc kernels copyin(area, var) if(acc_is_present(var))
443  !$acc loop reduction(+:s1,s2)
444  do j = js, je
445  !$acc loop reduction(+:s1,s2)
446  do i = is, ie
447  if ( var(i,j) /= undef ) then
448  s1 = s1 + area(i,j) * var(i,j)
449  s2 = s2 + area(i,j)
450  endif
451  enddo
452  enddo
453  !$acc end kernels
454  statval(1) = s1
455  statval(2) = s2
456 
457  call prof_rapstart('COMM_Allreduce', 2)
458  ! All reduce
459  call mpi_allreduce( statval(:), &
460  allstatval(:), &
461  2, &
462  mpi_double_precision, &
463  mpi_sum, &
464  prc_local_comm_world, &
465  ierr )
466  call prof_rapend ('COMM_Allreduce', 2)
467 
468  if ( allstatval(2) > 0.0_dp ) then
469  varmean = allstatval(1) / allstatval(2)
470  else
471  varmean = undef
472  end if
473 
474  return

References scale_const::const_undef, scale_prc::prc_local_comm_world, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_horizontal_mean_3d()

subroutine scale_statistics::statistics_horizontal_mean_3d ( integer, intent(in)  KA,
integer, intent(in)  KS,
integer, intent(in)  KE,
integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
real(rp), dimension (ka,ia,ja), intent(in)  var,
real(rp), dimension( ia,ja), intent(in)  area,
real(rp), dimension(ka), intent(out)  varmean 
)

Definition at line 481 of file scale_statistics.F90.

481  use scale_const, only: &
482  undef => const_undef
483  integer, intent(in) :: KA, KS, KE
484  integer, intent(in) :: IA, IS, IE
485  integer, intent(in) :: JA, JS, JE
486 
487  real(RP), intent(in) :: var (KA,IA,JA)
488  real(RP), intent(in) :: area( IA,JA)
489  real(RP), intent(out) :: varmean(KA)
490 
491  real(DP) :: statval (KS:KE,2)
492  real(DP) :: allstatval(KS:KE,2)
493 
494  integer :: ierr
495  integer :: k, i, j
496  !---------------------------------------------------------------------------
497 
498  !$acc data copyin(var, area) copyout(varmean) create(statval, allstatval) if(acc_is_present(var))
499 
500  !$acc kernels if(acc_is_present(var))
501  statval(:,:) = 0.0_dp
502  !$acc end kernels
503 
504 ! !$omp parallel do reduction(+:statval)
505  !$acc kernels if(acc_is_present(var))
506  !$acc loop independent
507  do j = js, je
508  !$acc loop independent
509  do i = is, ie
510  do k = ks, ke
511  if ( var(k,i,j) /= undef ) then
512  !$acc atomic update
513  statval(k,1) = statval(k,1) + area(i,j) * var(k,i,j)
514  !$acc end atomic
515  !$acc atomic update
516  statval(k,2) = statval(k,2) + area(i,j)
517  !$acc end atomic
518  endif
519  enddo
520  enddo
521  enddo
522  !$acc end kernels
523 
524  call prof_rapstart('COMM_Allreduce', 2)
525  ! All reduce
526  !$acc host_data use_device(statval, allstatval) if(acc_is_present(var))
527  call mpi_allreduce( statval(:,:), &
528  allstatval(:,:), &
529  (ke-ks+1)*2, &
530  mpi_double_precision, &
531  mpi_sum, &
532  prc_local_comm_world, &
533  ierr )
534  !$acc end host_data
535  call prof_rapend ('COMM_Allreduce', 2)
536 
537  !$acc kernels if(acc_is_present(var))
538  do k = ks, ke
539  if ( allstatval(k,2) > 0.0_dp ) then
540  varmean(k) = allstatval(k,1) / allstatval(k,2)
541  else
542  varmean(k) = undef
543  end if
544  enddo
545  !$acc end kernels
546  !$acc kernels if(acc_is_present(var))
547  do k = 1, ks-1
548  varmean(k) = undef
549  end do
550  !$acc end kernels
551  !$acc kernels if(acc_is_present(var))
552  do k = ke+1, ka
553  varmean(k) = undef
554  end do
555  !$acc end kernels
556 
557  !$acc end data
558 
559  return

References scale_const::const_undef, scale_prc::prc_local_comm_world, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_horizontal_min_2d()

subroutine scale_statistics::statistics_horizontal_min_2d ( integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
real(rp), dimension(ia,ja), intent(in)  var,
real(rp), intent(out)  varmin 
)

Calc horizontal minimum value.

Definition at line 568 of file scale_statistics.F90.

568  use scale_const, only: &
569  undef => const_undef, &
570  huge => const_huge
571  use scale_comm_cartesc, only: &
573  integer, intent(in) :: IA, IS, IE
574  integer, intent(in) :: JA, JS, JE
575 
576  real(RP), intent(in) :: var(IA,JA)
577  real(RP), intent(out) :: varmin
578 
579  real(RP) :: statval
580  real(RP) :: allstatval
581 
582  integer :: ierr
583  integer :: i, j
584  !---------------------------------------------------------------------------
585 
586  statval = huge
587  !$omp parallel do reduction(min:statval)
588  !$acc kernels copyin(var) if(acc_is_present(var))
589  !$acc loop reduction(min:statval)
590  do j = js, je
591  !$acc loop reduction(min:statval)
592  do i = is, ie
593  if ( var(i,j) /= undef .and. var(i,j) < statval ) then
594  statval = var(i,j)
595  endif
596  enddo
597  enddo
598  !$acc end kernels
599 
600  call prof_rapstart('COMM_Allreduce', 2)
601  ! All reduce
602  call mpi_allreduce( statval, &
603  allstatval, &
604  1, &
605  comm_datatype, &
606  mpi_min, &
607  prc_local_comm_world, &
608  ierr )
609  call prof_rapend ('COMM_Allreduce', 2)
610 
611  if ( allstatval < huge ) then
612  varmin = allstatval
613  else
614  varmin = undef
615  end if
616 
617  return

References scale_comm_cartesc::comm_datatype, scale_const::const_huge, scale_const::const_undef, scale_prc::prc_local_comm_world, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_horizontal_min_3d()

subroutine scale_statistics::statistics_horizontal_min_3d ( integer, intent(in)  KA,
integer, intent(in)  KS,
integer, intent(in)  KE,
integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
real(rp), dimension(ka,ia,ja), intent(in)  var,
real(rp), dimension(ka), intent(out)  varmin 
)

Definition at line 624 of file scale_statistics.F90.

624  use scale_const, only: &
625  undef => const_undef, &
626  huge => const_huge
627  use scale_comm_cartesc, only: &
629  integer, intent(in) :: KA, KS, KE
630  integer, intent(in) :: IA, IS, IE
631  integer, intent(in) :: JA, JS, JE
632 
633  real(RP), intent(in) :: var(KA,IA,JA)
634  real(RP), intent(out) :: varmin(KA)
635 
636  real(RP) :: statval (KA)
637  real(RP) :: allstatval(KA)
638 
639  integer :: ierr
640  integer :: k, i, j
641  !---------------------------------------------------------------------------
642 
643  !$acc data copyin(var) copyout(varmin) create(statval, allstatval) if(acc_is_present(var))
644 
645  !$acc kernels if(acc_is_present(var))
646  statval(:) = huge
647  !$acc end kernels
648 
649 ! !$omp parallel do reduction(min:statval)
650  !$acc kernels if(acc_is_present(var))
651  !$acc loop independent
652  do j = js, je
653  !$acc loop independent
654  do i = is, ie
655  do k = ks, ke
656 #ifdef _OPENACC
657  if ( var(k,i,j) /= undef )then
658  !$acc atomic update
659  statval(k) = min( var(k,i,j), statval(k) )
660  !$acc end atomic
661  endif
662 #else
663  if ( var(k,i,j) /= undef .and. var(k,i,j) < statval(k) ) then
664  statval(k) = var(k,i,j)
665  endif
666 #endif
667  enddo
668  enddo
669  enddo
670  !$acc end kernels
671 
672  call prof_rapstart('COMM_Allreduce', 2)
673  ! All reduce
674  !$acc host_data use_device(statval, allstatval) if(acc_is_present(var))
675  call mpi_allreduce( statval(ks:ke), &
676  allstatval(ks:ke), &
677  ke-ks+1, &
678  comm_datatype, &
679  mpi_min, &
680  prc_local_comm_world, &
681  ierr )
682  !$acc end host_data
683  call prof_rapend ('COMM_Allreduce', 2)
684 
685  !$acc kernels if(acc_is_present(var))
686  do k = ks, ke
687  if ( allstatval(k) < huge ) then
688  varmin(k) = allstatval(k)
689  else
690  varmin(k) = undef
691  end if
692  enddo
693  !$acc end kernels
694  !$acc kernels if(acc_is_present(var))
695  do k = 1, ks-1
696  varmin(k) = undef
697  end do
698  !$acc end kernels
699  !$acc kernels if(acc_is_present(var))
700  do k = ke+1, ka
701  varmin(k) = undef
702  end do
703  !$acc end kernels
704 
705  !$acc end data
706 
707  return

References scale_comm_cartesc::comm_datatype, scale_const::const_huge, scale_const::const_undef, scale_prc::prc_local_comm_world, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_horizontal_max_2d()

subroutine scale_statistics::statistics_horizontal_max_2d ( integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
real(rp), dimension(ia,ja), intent(in)  var,
real(rp), intent(out)  varmax 
)

Calc horizontal maximum value.

Definition at line 716 of file scale_statistics.F90.

716  use scale_const, only: &
717  undef => const_undef, &
718  huge => const_huge
719  use scale_comm_cartesc, only: &
721  integer, intent(in) :: IA, IS, IE
722  integer, intent(in) :: JA, JS, JE
723 
724  real(RP), intent(in) :: var(IA,JA)
725  real(RP), intent(out) :: varmax
726 
727  real(RP) :: statval
728  real(RP) :: allstatval
729 
730  integer :: ierr
731  integer :: i, j
732  !---------------------------------------------------------------------------
733 
734  statval = - huge
735  !$omp parallel do reduction(max:statval)
736  !$acc kernels copyin(var) if(acc_is_present(var))
737  !$acc loop reduction(max:statval)
738  do j = js, je
739  !$acc loop reduction(max:statval)
740  do i = is, ie
741  if ( var(i,j) /= undef .and. var(i,j) > statval ) then
742  statval = var(i,j)
743  endif
744  enddo
745  enddo
746  !$acc end kernels
747 
748  call prof_rapstart('COMM_Allreduce', 2)
749  ! All reduce
750  call mpi_allreduce( statval, &
751  allstatval, &
752  1, &
753  comm_datatype, &
754  mpi_max, &
755  prc_local_comm_world, &
756  ierr )
757  call prof_rapend ('COMM_Allreduce', 2)
758 
759  if ( allstatval > - huge ) then
760  varmax = allstatval
761  else
762  varmax = undef
763  end if
764 
765  return

References scale_comm_cartesc::comm_datatype, scale_const::const_huge, scale_const::const_undef, scale_prc::prc_local_comm_world, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_horizontal_max_3d()

subroutine scale_statistics::statistics_horizontal_max_3d ( integer, intent(in)  KA,
integer, intent(in)  KS,
integer, intent(in)  KE,
integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
real(rp), dimension(ka,ia,ja), intent(in)  var,
real(rp), dimension(ka), intent(out)  varmax 
)

Definition at line 772 of file scale_statistics.F90.

772  use scale_const, only: &
773  undef => const_undef, &
774  huge => const_huge
775  use scale_comm_cartesc, only: &
777  integer, intent(in) :: KA, KS, KE
778  integer, intent(in) :: IA, IS, IE
779  integer, intent(in) :: JA, JS, JE
780 
781  real(RP), intent(in) :: var(KA,IA,JA)
782  real(RP), intent(out) :: varmax(KA)
783 
784  real(RP) :: statval (KA)
785  real(RP) :: allstatval(KA)
786 
787  integer :: ierr
788  integer :: k, i, j
789  !---------------------------------------------------------------------------
790 
791  !$acc data copyin(var) copyout(varmax) create(statval, allstatval) if(acc_is_present(var))
792 
793  !$acc kernels if(acc_is_present(var))
794  statval(:) = - huge
795  !$acc end kernels
796 
797 ! !$omp parallel do reduction(max:statval)
798  !$acc kernels if(acc_is_present(var))
799  !$acc loop independent
800  do j = js, je
801  !$acc loop independent
802  do i = is, ie
803  do k = ks, ke
804 #ifdef _OPENACC
805  if ( var(k,i,j) /= undef ) then
806  !$acc atomic update
807  statval(k) = max( var(k,i,j), statval(k) )
808  !$acc end atomic
809  endif
810 #else
811  if ( var(k,i,j) /= undef .and. var(k,i,j) > statval(k) ) then
812  statval(k) = var(k,i,j)
813  endif
814 #endif
815  enddo
816  enddo
817  enddo
818  !$acc end kernels
819 
820  call prof_rapstart('COMM_Allreduce', 2)
821  ! All reduce
822  !$acc host_data use_device(statval, allstatval) if(acc_is_present(var))
823  call mpi_allreduce( statval(ks:ke), &
824  allstatval(ks:ke), &
825  ke-ks+1, &
826  comm_datatype, &
827  mpi_max, &
828  prc_local_comm_world, &
829  ierr )
830  !$acc end host_data
831  call prof_rapend ('COMM_Allreduce', 2)
832 
833  !$acc kernels if(acc_is_present(var))
834  do k = ks, ke
835  if ( allstatval(k) > - huge ) then
836  varmax(k) = allstatval(k)
837  else
838  varmax(k) = undef
839  end if
840  enddo
841  !$acc end kernels
842  !$acc kernels if(acc_is_present(var))
843  do k = 1, ks-1
844  varmax(k) = undef
845  end do
846  !$acc end kernels
847  !$acc kernels if(acc_is_present(var))
848  do k = ke+1, ka
849  varmax(k) = undef
850  end do
851  !$acc end kernels
852 
853  !$acc end data
854 
855  return

References scale_comm_cartesc::comm_datatype, scale_const::const_huge, scale_const::const_undef, scale_prc::prc_local_comm_world, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_detail_3d()

subroutine scale_statistics::statistics_detail_3d ( integer, intent(in)  KA,
integer, intent(in)  KS,
integer, intent(in)  KE,
integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
integer, intent(in)  VA,
character(len=*), dimension(va), intent(in)  varname,
real(rp), dimension(ka,ia,ja,va), intent(in)  var,
logical, intent(in), optional  local 
)

Search global maximum & minimum value.

Parameters
[in]varnamename of item
[in]varvalues
[in]localcalc in local node

Definition at line 864 of file scale_statistics.F90.

864  use scale_prc, only: &
865  prc_nprocs
866  use scale_comm_cartesc, only: &
868  implicit none
869 
870  integer, intent(in) :: KA, KS, KE
871  integer, intent(in) :: IA, IS, IE
872  integer, intent(in) :: JA, JS, JE
873  integer, intent(in) :: VA
874 
875  character(len=*), intent(in) :: varname(VA)
876  real(RP), intent(in) :: var(KA,IA,JA,VA)
877 
878  logical, intent(in), optional :: local
879 
880  real(RP) :: statval_l ( VA,2)
881  integer :: statidx_l (3,VA,2)
882  real(RP) :: statval ( VA,2,0:PRC_nprocs-1)
883  integer :: statidx (3,VA,2,0:PRC_nprocs-1)
884  real(RP) :: allstatval(VA,2)
885  integer :: allstatidx(VA,2)
886  logical :: do_globalcomm
887 
888  integer :: k, i, j
889  integer :: ierr
890  integer :: v, p
891  !---------------------------------------------------------------------------
892 
893  do_globalcomm = statistics_use_globalcomm
894  if ( present(local) ) do_globalcomm = ( .not. local )
895 
896  !$acc update host(var) if ( acc_is_present(var) )
897 
898  log_newline
899  log_info("STATISTICS_detail_3D",*) 'Variable Statistics '
900  do v = 1, va
901  statval_l( v,:) = var(ks,is,js,v)
902  statidx_l(1,v,:) = ks
903  statidx_l(2,v,:) = is
904  statidx_l(3,v,:) = js
905  do j = js, je
906  do i = is, ie
907  do k = ks, ke
908  if ( var(k,i,j,v) > statval_l(v,1) ) then
909  statval_l( v,1) = var(k,i,j,v)
910  statidx_l(1,v,1) = k
911  statidx_l(2,v,1) = i
912  statidx_l(3,v,1) = j
913  end if
914  if ( var(k,i,j,v) < statval_l(v,2) ) then
915  statval_l( v,2) = var(k,i,j,v)
916  statidx_l(1,v,2) = k
917  statidx_l(2,v,2) = i
918  statidx_l(3,v,2) = j
919  end if
920  end do
921  end do
922  end do
923  enddo
924 
925  if ( do_globalcomm ) then
926  call prof_rapstart('COMM_Bcast', 2)
927 
928  call mpi_allgather( statval_l(:,:), &
929  va*2, &
930  comm_datatype, &
931  statval(:,:,:), &
932  va*2, &
933  comm_datatype, &
934  prc_local_comm_world, &
935  ierr )
936 
937  call mpi_allgather( statidx_l(:,:,:), &
938  3*va*2, &
939  mpi_integer, &
940  statidx(:,:,:,:), &
941  3*va*2, &
942  mpi_integer, &
943  prc_local_comm_world, &
944  ierr )
945 
946  call prof_rapend ('COMM_Bcast', 2)
947 
948  do v = 1, va
949  allstatval(v,1) = statval(v,1,0)
950  allstatval(v,2) = statval(v,2,0)
951  allstatidx(v,:) = 0
952  do p = 1, prc_nprocs-1
953  if ( statval(v,1,p) > allstatval(v,1) ) then
954  allstatval(v,1) = statval(v,1,p)
955  allstatidx(v,1) = p
956  end if
957  if ( statval(v,2,p) < allstatval(v,2) ) then
958  allstatval(v,2) = statval(v,2,p)
959  allstatidx(v,2) = p
960  end if
961  end do
962  log_info_cont(*) '[', trim(varname(v)), ']'
963  log_info_cont('(1x,A,ES17.10,A,4(I5,A))') ' MAX =', &
964  allstatval(v,1), ' (rank=', &
965  allstatidx(v,1), '; ', &
966  statidx(1,v,1,allstatidx(v,1)),',', &
967  statidx(2,v,1,allstatidx(v,1)),',', &
968  statidx(3,v,1,allstatidx(v,1)),')'
969  log_info_cont('(1x,A,ES17.10,A,4(I5,A))') ' MIN =', &
970  allstatval(v,2), ' (rank=', &
971  allstatidx(v,2), '; ', &
972  statidx(1,v,2,allstatidx(v,2)),',', &
973  statidx(2,v,2,allstatidx(v,2)),',', &
974  statidx(3,v,2,allstatidx(v,2)),')'
975  enddo
976  else
977  ! statistics on each node
978  do v = 1, va
979  log_info_cont(*) '[', trim(varname(v)), ']'
980  log_info_cont('(1x,A,ES17.10,A,3(I5,A))') 'MAX = ', &
981  statval_l( v,1),' (', &
982  statidx_l(1,v,1),',', &
983  statidx_l(2,v,1),',', &
984  statidx_l(3,v,1),')'
985  log_info_cont('(1x,A,ES17.10,A,3(I5,A))') 'MIN = ', &
986  statval_l( v,2),' (', &
987  statidx_l(1,v,2),',', &
988  statidx_l(2,v,2),',', &
989  statidx_l(3,v,2),')'
990  enddo
991  endif
992 
993  log_newline
994 
995  return

References scale_comm_cartesc::comm_datatype, scale_prc::prc_nprocs, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_detail_2d()

subroutine scale_statistics::statistics_detail_2d ( integer, intent(in)  IA,
integer, intent(in)  IS,
integer, intent(in)  IE,
integer, intent(in)  JA,
integer, intent(in)  JS,
integer, intent(in)  JE,
integer, intent(in)  VA,
character(len=*), dimension(va), intent(in)  varname,
real(rp), dimension(ia,ja,va), intent(in)  var,
logical, intent(in), optional  local 
)
Parameters
[in]varnamename of item
[in]varvalues

Definition at line 1002 of file scale_statistics.F90.

1002  use scale_prc, only: &
1003  prc_nprocs
1004  use scale_comm_cartesc, only: &
1006  implicit none
1007 
1008  integer, intent(in) :: IA, IS, IE
1009  integer, intent(in) :: JA, JS, JE
1010  integer, intent(in) :: VA
1011 
1012  character(len=*), intent(in) :: varname(VA)
1013  real(RP), intent(in) :: var(IA,JA,VA)
1014 
1015  logical, intent(in), optional :: local ! calc in local node
1016 
1017  real(RP) :: statval_l ( VA,2)
1018  integer :: statidx_l (2,VA,2)
1019  real(RP) :: statval ( VA,2,0:PRC_nprocs-1)
1020  integer :: statidx (2,VA,2,0:PRC_nprocs-1)
1021  real(RP) :: allstatval(VA,2)
1022  integer :: allstatidx(VA,2)
1023  logical :: do_globalcomm
1024 
1025  integer :: i, j
1026  integer :: ierr
1027  integer :: v, p
1028  !---------------------------------------------------------------------------
1029 
1030  do_globalcomm = statistics_use_globalcomm
1031  if ( present(local) ) do_globalcomm = ( .not. local )
1032 
1033  !$acc update host(var) if ( acc_is_present(var) )
1034 
1035  log_newline
1036  log_info("STATISTICS_detail_2D",*) 'Variable Statistics '
1037  do v = 1, va
1038  statval_l( v,:) = var(is,js,v)
1039  statidx_l(1,v,:) = is
1040  statidx_l(2,v,:) = js
1041  do j = js, je
1042  do i = is, ie
1043  if ( var(i,j,v) > statval_l(v,1) ) then
1044  statval_l( v,1) = var(i,j,v)
1045  statidx_l(1,v,1) = i
1046  statidx_l(2,v,1) = j
1047  end if
1048  if ( var(i,j,v) < statval_l(v,2) ) then
1049  statval_l( v,2) = var(i,j,v)
1050  statidx_l(1,v,2) = i
1051  statidx_l(2,v,2) = j
1052  end if
1053  end do
1054  end do
1055  enddo
1056 
1057  if ( do_globalcomm ) then
1058  call prof_rapstart('COMM_Bcast', 2)
1059 
1060  call mpi_allgather( statval_l(:,:), &
1061  va*2, &
1062  comm_datatype, &
1063  statval(:,:,:), &
1064  va*2, &
1065  comm_datatype, &
1066  prc_local_comm_world, &
1067  ierr )
1068 
1069  call mpi_allgather( statidx_l(:,:,:), &
1070  2*va*2, &
1071  mpi_integer, &
1072  statidx(:,:,:,:), &
1073  2*va*2, &
1074  mpi_integer, &
1075  prc_local_comm_world, &
1076  ierr )
1077 
1078  call prof_rapend ('COMM_Bcast', 2)
1079 
1080  do v = 1, va
1081  allstatval(v,1) = statval(v,1,0)
1082  allstatval(v,2) = statval(v,2,0)
1083  allstatidx(v,:) = 0
1084  do p = 1, prc_nprocs-1
1085  if ( statval(v,1,p) > allstatval(v,1) ) then
1086  allstatval(v,1) = statval(v,1,p)
1087  allstatidx(v,1) = p
1088  end if
1089  if ( statval(v,2,p) < allstatval(v,2) ) then
1090  allstatval(v,2) = statval(v,2,p)
1091  allstatidx(v,2) = p
1092  end if
1093  end do
1094  log_info_cont(*) '[', trim(varname(v)), ']'
1095  log_info_cont('(1x,A,ES17.10,A,3(I5,A))') ' MAX =', &
1096  allstatval(v,1), ' (rank=', &
1097  allstatidx(v,1), '; ', &
1098  statidx(1,v,1,allstatidx(v,1)),',', &
1099  statidx(2,v,1,allstatidx(v,1)),')'
1100  log_info_cont('(1x,A,ES17.10,A,3(I5,A))') ' MIN =', &
1101  allstatval(v,2), ' (rank=', &
1102  allstatidx(v,2), '; ', &
1103  statidx(1,v,2,allstatidx(v,2)),',', &
1104  statidx(2,v,2,allstatidx(v,2)),')'
1105  enddo
1106  else
1107  ! statistics on each node
1108  do v = 1, va
1109  log_info_cont(*) '[', trim(varname(v)), ']'
1110  log_info_cont('(1x,A,ES17.10,A,2(I5,A))') 'MAX = ', &
1111  statval_l( v,1),' (', &
1112  statidx_l(1,v,1),',', &
1113  statidx_l(2,v,1),')'
1114  log_info_cont('(1x,A,ES17.10,A,2(I5,A))') 'MIN = ', &
1115  statval_l( v,2),' (', &
1116  statidx_l(1,v,2),',', &
1117  statidx_l(2,v,2),')'
1118  enddo
1119  endif
1120 
1121  log_newline
1122 
1123  return

References scale_comm_cartesc::comm_datatype, scale_prc::prc_nprocs, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Here is the call graph for this function:

◆ statistics_summation_1d()

real(rp) function scale_statistics::statistics_summation_1d ( integer, intent(in)  IA,
real(rp), dimension(ia), intent(in)  ARRAY,
real(rp), intent(in), optional  UNDEF 
)

Definition at line 1136 of file scale_statistics.F90.

1136  implicit none
1137 
1138  ! arguments
1139  integer, intent(in) :: IA
1140  real(RP), intent(in) :: ARRAY(IA)
1141 
1142  real(RP), intent(in), optional :: UNDEF
1143 
1144  ! function result
1145  real(RP) :: STATISTICS_summation_1D
1146 
1147  ! work
1148  integer :: i
1149 
1150  real(RP) :: tmp_ARRAY(IA)
1151  real(RP) :: tmp
1152  real(RP) :: res
1153  !---------------------------------------------------------------------------
1154 
1155  if( present( undef ) ) then
1156  if( any( abs( array(:) - undef ) > eps ) ) then
1157  where( abs( array(:) - undef ) > eps )
1158  tmp_array(:) = array(:)
1159  else where
1160  tmp_array(:) = 0.0_rp
1161  end where
1162  else
1163  statistics_summation_1d = undef
1164  return ! end function
1165  end if
1166  else
1167  tmp_array(:) = array(:)
1168  end if
1169 
1170  statistics_summation_1d = 0.0_rp
1171 
1172  tmp = 0.0_rp
1173  res = 0.0_rp
1174 
1175  ! Kahan's Compensated Summation (Kahan 1965)
1176  do i = 1, ia
1177  tmp = statistics_summation_1d + ( tmp_array(i) + res )
1178  res = ( tmp_array(i) + res ) - ( tmp - statistics_summation_1d )
1179  statistics_summation_1d = tmp
1180  end do
1181 
1182  return

References scale_precision::rp.

◆ statistics_covariance()

real(rp) function, public scale_statistics::statistics_covariance ( integer, intent(in)  IA1,
integer, intent(in)  IA2,
real(rp), dimension(ia1), intent(in)  ARRAY1,
real(rp), dimension(ia2), intent(in)  ARRAY2,
real(rp), intent(in), optional  UNDEF 
)

Definition at line 1671 of file scale_statistics.F90.

1671  implicit none
1672 
1673  ! argument
1674  integer, intent(in) :: IA1, IA2
1675  real(RP), intent(in) :: ARRAY1(IA1)
1676  real(RP), intent(in) :: ARRAY2(IA2)
1677 
1678  real(RP), intent(in), optional :: UNDEF
1679 
1680  ! function result
1681  real(RP) :: STATISTICS_covariance
1682 
1683  ! works
1684  real(RP) :: tmp_ARRAY1(IA1)
1685  real(RP) :: tmp_ARRAY2(IA2)
1686 
1687  real(RP), allocatable :: fixed_ARRAY1(:)
1688  real(RP), allocatable :: fixed_ARRAY2(:)
1689 
1690  integer :: n
1691  !---------------------------------------------------------------------------
1692 
1693  if( present( undef ) ) then
1694  n = statistics_undef_arraysize( ia1, ia2, array1(:), array2(:), undef )
1695 
1696  if( n > 1 ) then
1697  ! embed undefined value each other
1698  tmp_array1 = statistics_undef_embed( ia1, ia2, array1(:), array2(:), undef )
1699  tmp_array2 = statistics_undef_embed( ia2, ia1, array2(:), array1(:), undef )
1700 
1701  allocate( fixed_array1(n) )
1702  allocate( fixed_array2(n) )
1703 
1704  ! generate arrays without undefined value
1705  fixed_array1(:) = pack( tmp_array1(:), abs( tmp_array1(:) - undef ) > eps )
1706  fixed_array2(:) = pack( tmp_array2(:), abs( tmp_array2(:) - undef ) > eps )
1707 
1708  statistics_covariance = statistics_summation( n, ( fixed_array1(:) - statistics_average( n, array1(:), undef ) ) &
1709  * ( fixed_array2(:) - statistics_average( n, array2(:), undef ) ) ) &
1710  / real( n-1, kind=rp )
1711 
1712  deallocate( fixed_array1 )
1713  deallocate( fixed_array2 )
1714  else
1715  statistics_covariance = undef
1716  end if
1717  else
1718  n = statistics_undef_arraysize( ia1, ia2, array1(:), array2(:) )
1719  ! simple covariance
1720  statistics_covariance = statistics_summation( n, ( array1(:) - statistics_average( n, array1(:) ) ) &
1721  * ( array2(:) - statistics_average( n, array2(:) ) ) ) &
1722  / real( n-1, kind=rp )
1723  end if
1724 
1725  return

References scale_precision::rp, statistics_undef_arraysize(), and statistics_undef_embed().

Referenced by statistics_correlation().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ statistics_correlation()

real(rp) function, public scale_statistics::statistics_correlation ( integer, intent(in)  IA1,
integer, intent(in)  IA2,
real(rp), dimension(ia1), intent(in)  ARRAY1,
real(rp), dimension(ia2), intent(in)  ARRAY2,
real(rp), intent(in), optional  UNDEF 
)

Definition at line 1739 of file scale_statistics.F90.

1739  implicit none
1740 
1741  ! arguments
1742  integer, intent(in) :: IA1, IA2
1743  real(RP), intent(in) :: ARRAY1(IA1)
1744  real(RP), intent(in) :: ARRAY2(IA2)
1745 
1746  real(RP), intent(in), optional :: UNDEF
1747 
1748  ! function result
1749  real(RP) :: STATISTICS_correlation
1750  !---------------------------------------------------------------------------
1751 
1752  if( present( undef ) ) then
1753  if( abs( statistics_covariance( ia1, ia2, array1(:), array2(:), undef ) - undef ) > eps .and. &
1754  abs( statistics_stddev( ia1, array1(:), undef ) - undef ) > eps .and. &
1755  abs( statistics_stddev( ia2, array2(:), undef ) - undef ) > eps .and. &
1756  abs( statistics_stddev( ia1, array1(:), undef ) ) > 0.0_rp .and. &
1757  abs( statistics_stddev( ia2, array2(:), undef ) ) > 0.0_rp ) then
1758  ! correlation coefficient without undefined value
1759  statistics_correlation = statistics_covariance( ia1, ia2, array1(:), array2(:), undef ) &
1760  / ( statistics_stddev( ia1, array1(:), undef ) &
1761  * statistics_stddev( ia2, array2(:), undef ) )
1762  else
1763  statistics_correlation = undef
1764  end if
1765  else
1766  ! simple correlation coefficient
1767  statistics_correlation = statistics_covariance( ia1, ia2, array1(:), array2(:) ) &
1768  / ( statistics_stddev( ia1, array1(:) ) &
1769  * statistics_stddev( ia2, array2(:) ) )
1770  end if
1771 
1772  return

References statistics_covariance().

Referenced by statistics_lag_correlation(), statistics_partial_correlation(), and statistics_regression().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ statistics_regression()

real function, public scale_statistics::statistics_regression ( integer, intent(in)  IA1,
integer, intent(in)  IA2,
real(rp), dimension(ia1), intent(in)  ARRAY1,
real(rp), dimension(ia2), intent(in)  ARRAY2,
real(rp), intent(in), optional  UNDEF 
)

Definition at line 1786 of file scale_statistics.F90.

1786  implicit none
1787 
1788  ! arguments
1789  integer, intent(in) :: IA1, IA2
1790  real(RP), intent(in) :: ARRAY1(IA1)
1791  real(RP), intent(in) :: ARRAY2(IA2)
1792 
1793  real(RP), intent(in), optional :: UNDEF
1794 
1795  ! function result
1796  real :: STATISTICS_regression
1797  !---------------------------------------------------------------------------
1798 
1799  if( present( undef ) ) then
1800  if( abs( statistics_correlation( ia1, ia2, array1(:), array2(:), undef ) - undef ) > eps .and. &
1801  abs( statistics_stddev( ia1, array1(:), undef ) - undef ) > eps .and. &
1802  abs( statistics_stddev( ia2, array2(:), undef ) - undef ) > eps .and. &
1803  abs( statistics_stddev( ia2, array2(:), undef ) ) > 0.0_rp ) then
1804  ! regression coefficient without undefined value
1805  statistics_regression = statistics_correlation( ia1, ia2, array1(:), array2(:), undef ) &
1806  * statistics_stddev( ia1, array1(:), undef ) &
1807  / statistics_stddev( ia2, array2(:), undef )
1808  else
1809  statistics_regression = undef
1810  end if
1811  else
1812  ! simple regression coefficient
1813  statistics_regression = statistics_correlation( ia1, ia2, array1(:), array2(:) ) &
1814  * statistics_stddev( ia1, array1(:) ) &
1815  / statistics_stddev( ia2, array2(:) )
1816  end if
1817 
1818  return

References statistics_correlation().

Here is the call graph for this function:

◆ statistics_lag_correlation()

real(rp) function, public scale_statistics::statistics_lag_correlation ( integer, intent(in)  IA1,
integer, intent(in)  IA2,
real(rp), dimension(ia1), intent(in)  ARRAY1,
real(rp), dimension(ia2), intent(in)  ARRAY2,
integer, intent(in)  LAG,
real(rp), intent(in), optional  UNDEF 
)

Definition at line 1833 of file scale_statistics.F90.

1833  implicit none
1834 
1835  ! arguments
1836  integer, intent(in) :: IA1, IA2
1837  real(RP), intent(in) :: ARRAY1(IA1)
1838  real(RP), intent(in) :: ARRAY2(IA2)
1839 
1840  integer, intent(in) :: LAG
1841 
1842  real(RP), intent(in), optional :: UNDEF
1843 
1844  ! function result
1845  real(RP) :: STATISTICS_lag_correlation
1846 
1847  ! works
1848  real(RP) :: tmp(IA2)
1849  !---------------------------------------------------------------------------
1850 
1851  tmp(:) = eoshift( array2(:), lag, undef )
1852 
1853  statistics_lag_correlation = statistics_correlation( ia1, ia2, array1(:), tmp(:), undef )
1854 
1855  return

References statistics_correlation().

Here is the call graph for this function:

◆ statistics_partial_correlation()

real(rp) function, public scale_statistics::statistics_partial_correlation ( integer, intent(in)  IA1,
integer, intent(in)  IA2,
integer, intent(in)  IA3,
real(rp), dimension(ia1), intent(in)  ARRAY1,
real(rp), dimension(ia2), intent(in)  ARRAY2,
real(rp), dimension(ia3), intent(in)  ARRAY3,
real(rp), intent(in), optional  UNDEF 
)

Definition at line 1870 of file scale_statistics.F90.

1870  implicit none
1871 
1872  ! arguments
1873  integer, intent(in) :: IA1, IA2, IA3
1874  real(RP), intent(in) :: ARRAY1(IA1)
1875  real(RP), intent(in) :: ARRAY2(IA2)
1876  real(RP), intent(in) :: ARRAY3(IA3)
1877 
1878  real(RP), intent(in), optional :: UNDEF
1879 
1880  ! function result
1881  real(RP) :: STATISTICS_partial_correlation
1882 
1883  ! works
1884  real(RP) :: corr12
1885  real(RP) :: corr13
1886  real(RP) :: corr23
1887  !---------------------------------------------------------------------------
1888 
1889  if( present( undef ) ) then
1890  corr12 = statistics_correlation( ia1, ia2, array1(:), array2(:), undef )
1891  corr13 = statistics_correlation( ia1, ia3, array1(:), array3(:), undef )
1892  corr23 = statistics_correlation( ia2, ia3, array2(:), array3(:), undef )
1893 
1894  if( corr12 == undef .or. &
1895  corr13 == undef .or. &
1896  corr23 == undef ) then
1897  ! return undefined value if correlation can't calculate
1898  statistics_partial_correlation = undef
1899  else if( corr12 >= 1.0_rp .or. &
1900  corr13 >= 1.0_rp .or. &
1901  corr23 >= 1.0_rp ) then
1902  ! return undefined value if correlation is 1
1903  statistics_partial_correlation = undef
1904  else
1905  ! partial correlation without undefined value
1906  statistics_partial_correlation = ( corr12 - ( corr13 * corr23 ) ) &
1907  / sqrt( ( 1.0_rp - corr13**2 ) * ( 1.0_rp - corr23**2 ) )
1908  endif
1909  else
1910  corr12 = statistics_correlation( ia1, ia2, array1(:), array2(:) )
1911  corr13 = statistics_correlation( ia1, ia3, array1(:), array3(:) )
1912  corr23 = statistics_correlation( ia2, ia3, array2(:), array3(:) )
1913 
1914  ! simple partial correlation
1915  statistics_partial_correlation = ( corr12 - ( corr13 * corr23 ) ) &
1916  / sqrt( ( 1.0_rp - corr13**2 ) * ( 1.0_rp - corr23**2 ) )
1917  end if
1918 
1919  return

References statistics_correlation().

Here is the call graph for this function:

◆ statistics_undef_replace()

real(rp) function, dimension(ia), public scale_statistics::statistics_undef_replace ( integer, intent(in)  IA,
real(rp), dimension(ia), intent(in)  ARRAY,
real(rp), intent(in)  UNDEF1,
real(rp), intent(in)  UNDEF2 
)

Definition at line 1928 of file scale_statistics.F90.

1928  implicit none
1929 
1930  ! arguments
1931  integer, intent(in) :: IA
1932  real(RP), intent(in) :: ARRAY(IA)
1933  real(RP), intent(in) :: UNDEF1
1934  real(RP), intent(in) :: UNDEF2
1935 
1936  ! function result
1937  real(RP) :: STATISTICS_undef_replace(IA)
1938 
1939  ! works
1940  integer :: i
1941  !---------------------------------------------------------------------------
1942 
1943  do i = 1, ia
1944  if( abs( array(i) - undef2 ) > eps ) then
1945  statistics_undef_replace(i) = array(i)
1946  else
1947  statistics_undef_replace(i) = undef1
1948  end if
1949  end do
1950 
1951  return

◆ statistics_undef_embed()

real(rp) function, dimension(ia1), public scale_statistics::statistics_undef_embed ( integer, intent(in)  IA1,
integer, intent(in)  IA2,
real(rp), dimension(ia1), intent(in)  ARRAY1,
real(rp), dimension(ia2), intent(in)  ARRAY2,
real(rp), intent(in), optional  UNDEF 
)

Definition at line 1960 of file scale_statistics.F90.

1960  implicit none
1961 
1962  ! arguments
1963  integer, intent(in) :: IA1, IA2
1964  real(RP), intent(in) :: ARRAY1(IA1)
1965  real(RP), intent(in) :: ARRAY2(IA2)
1966 
1967  real(RP), intent(in), optional :: UNDEF
1968 
1969  ! function result
1970  real(RP) :: STATISTICS_undef_embed(IA1)
1971 
1972  ! works
1973  integer :: i
1974  !---------------------------------------------------------------------------
1975 
1976  ! check the length of input arrays
1977  if( ia1 /= ia2 ) then
1978  write(*,*) "Error: different array size !!"
1979  stop
1980  end if
1981 
1982  do i = 1, ia1
1983  if( abs( array2(i) - undef ) > eps ) then
1984  statistics_undef_embed(i) = array1(i)
1985  else
1986  statistics_undef_embed(i) = undef
1987  end if
1988  end do

Referenced by statistics_covariance().

Here is the caller graph for this function:

◆ statistics_undef_arraysize()

integer function, public scale_statistics::statistics_undef_arraysize ( integer, intent(in)  IA1,
integer, intent(in)  IA2,
real(rp), dimension(ia1), intent(in)  ARRAY1,
real(rp), dimension(ia2), intent(in)  ARRAY2,
real(rp), intent(in), optional  UNDEF 
)

Definition at line 1997 of file scale_statistics.F90.

1997  implicit none
1998 
1999  ! arguments
2000  integer, intent(in) :: IA1, IA2
2001  real(RP), intent(in) :: ARRAY1(IA1)
2002  real(RP), intent(in) :: ARRAY2(IA2)
2003 
2004  real(RP), intent(in), optional :: UNDEF
2005 
2006  ! function result
2007  integer :: STATISTICS_undef_arraysize
2008 
2009  ! works
2010  integer :: i, cnt
2011  !---------------------------------------------------------------------------
2012 
2013  ! check the length of input arrays
2014  if( ia1 /= ia2 ) then
2015  write(*,*) "Error: different array size !!"
2016  stop
2017  end if
2018 
2019  if( present( undef ) ) then
2020  cnt = 0
2021  do i = 1, ia1
2022  if( abs( array1(i) - undef ) > eps .and. &
2023  abs( array2(i) - undef ) > eps ) then
2024  cnt = cnt + 1
2025  end if
2026  end do
2027  ! array size excepting number of undefined values
2028  statistics_undef_arraysize = cnt
2029  else
2030  ! simple array size
2031  statistics_undef_arraysize = ia1
2032  end if
2033 
2034  return

Referenced by statistics_covariance().

Here is the caller graph for this function:

Variable Documentation

◆ statistics_checktotal

logical, public scale_statistics::statistics_checktotal = .false.
scale_comm_cartesc::comm_setup
subroutine, public comm_setup
Setup.
Definition: scale_comm_cartesC.F90:200
scale_comm_cartesc::comm_datatype
integer, public comm_datatype
datatype of variable
Definition: scale_comm_cartesC.F90:105
scale_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:350
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:35
scale_prc::prc_myrank
integer, public prc_myrank
process num in local communicator
Definition: scale_prc.F90:91
scale_index::va
integer, public va
Definition: scale_index.F90:35
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_precision::rp
integer, parameter, public rp
Definition: scale_precision.F90:41
scale_atmos_grid_cartesc_index::ie
integer, public ie
end point of inner domain: x, local
Definition: scale_atmos_grid_cartesC_index.F90:54
scale_const
module CONSTANT
Definition: scale_const.F90:11
scale_atmos_grid_cartesc_index::ia
integer, public ia
Definition: scale_atmos_grid_cartesC_index.F90:48
scale_const::const_huge
real(rp), public const_huge
huge number
Definition: scale_const.F90:37
scale_atmos_grid_cartesc_index::is
integer, public is
start point of inner domain: x, local
Definition: scale_atmos_grid_cartesC_index.F90:53
scale_prc::prc_nprocs
integer, public prc_nprocs
myrank in local communicator
Definition: scale_prc.F90:90
scale_comm_cartesc
module COMMUNICATION
Definition: scale_comm_cartesC.F90:11
scale_atmos_grid_cartesc_index::js
integer, public js
start point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:55
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:43
scale_atmos_grid_cartesc_index::je
integer, public je
end point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:56