FE-Project
Loading...
Searching...
No Matches
scale_meshfield_statistics.F90
Go to the documentation of this file.
1
9#include "scaleFElib.h"
11
12 !-----------------------------------------------------------------------------
13 !
14 !++ used modules
15 !
16 use mpi ! TODO: to replace functions in scale_comm module
17 use scale_precision
18 use scale_io
19 use scale_prc, only: &
20 prc_local_comm_world
21 use scale_prof
22
26 use scale_meshfield_base, only: &
28
29 !-----------------------------------------------------------------------------
30 implicit none
31 private
32
33 !-----------------------------------------------------------------------------
34 !
35 !++ Public type & procedure
36 !
38
40 module procedure :: meshfield_statistics_detail_1d
41 module procedure :: meshfield_statistics_detail_2d
42 module procedure :: meshfield_statistics_detail_3d
43 end interface
45
47 module procedure :: meshfield_statistics_maxmin_1d
48 module procedure :: meshfield_statistics_maxmin_2d
49 module procedure :: meshfield_statistics_maxmin_3d
50 end interface
52
54 module procedure :: meshfield_statistics_total_1d
55 module procedure :: meshfield_statistics_total_2d
56 module procedure :: meshfield_statistics_total_3d
57 end interface
59
60 !-----------------------------------------------------------------------------
61 !
62 !++ Public parameters & variables
63 !
64
65 !-----------------------------------------------------------------------------
66 !
67 !++ Private type & procedure
68 !
69 type :: meshfield_statistics_base
70 logical :: use_globalcomm
71 integer :: comm_datatype
72 end type meshfield_statistics_base
73
74 !-----------------------------------------------------------------------------
75 !
76 !++ Private parameters & variables
77 !
78 type(MeshField_statistics_base), private :: base
79
80contains
81
82 !-----------------------------------------------------------------------------
84!OCL SERIAL
86 use scale_prc, only: &
87 prc_abort
88 implicit none
89
90 logical :: use_globalcomm = .false.
91 namelist / param_meshfield_statistics / &
92 use_globalcomm
93
94 integer :: ierr
95 !---------------------------------------------------------------------------
96
97 log_newline
98 log_info("MeshField_statistics_setup",*) 'Setup'
99
100 !--- read namelist
101 rewind(io_fid_conf)
102 read(io_fid_conf,nml=param_meshfield_statistics,iostat=ierr)
103 if( ierr < 0 ) then !--- missing
104 log_info("MeshField_statistics_setup",*) 'Not found namelist. Default used.'
105 elseif( ierr > 0 ) then !--- fatal error
106 log_error("MeshField_statistics_setup",*) 'Not appropriate names in namelist PARAM_MeshField_statistics. Check!'
107 call prc_abort
108 endif
109 log_nml(param_meshfield_statistics)
110
111 base%use_globalcomm = use_globalcomm
112
113 if ( base%use_globalcomm ) then
114 log_info_cont(*) '=> The total is calculated for the global domain.'
115 else
116 log_info_cont(*) '=> The total is calculated for the local domain.'
117 endif
118
119 if ( rp == kind(0.d0) ) then
120 base%comm_datatype = mpi_double_precision
121 elseif( rp == kind(0.0) ) then
122 base%comm_datatype = mpi_real
123 else
124 log_error("MeshField_statistics_setup",*) 'precision is not supportd'
125 call prc_abort
126 endif
127
128 return
129 end subroutine meshfield_statistics_setup
130
131 !-----------------------------------------------------------------------------
133!OCL SERIAL
134 subroutine meshfield_statistics_total_1d( field, & ! (in)
135 log_suppress, global, & ! (in)
136 mean, sum ) ! (out)
137
138 implicit none
139
140 class(meshfield1d), intent(in) :: field
141 logical, intent(in), optional :: log_suppress
142 logical, intent(in), optional :: global
143 real(RP), intent(out), optional :: mean
144 real(DP), intent(out), optional :: sum
145
146 real(DP) :: statval
147 real(DP) :: total
148 !---------------------------------------------------------------------------
149
150 call calculate_statval( field, field%mesh%lcmesh_list(:), & ! (in)
151 statval, total ) ! (out)
152
153 call statistics_total_core( &
154 field%varname, statval, total, & ! (in)
155 log_suppress, global, & ! (in)
156 mean, sum ) ! (out)
157
158 return
159 end subroutine meshfield_statistics_total_1d
160
161 !-----------------------------------------------------------------------------
163!OCL SERIAL
164 subroutine meshfield_statistics_total_2d( field, & ! (in)
165 log_suppress, global, & ! (in)
166 mean, sum ) ! (out)
167
168 implicit none
169
170 class(meshfield2d), intent(in) :: field
171 logical, intent(in), optional :: log_suppress
172 logical, intent(in), optional :: global
173 real(RP), intent(out), optional :: mean
174 real(DP), intent(out), optional :: sum
175
176 real(DP) :: statval
177 real(DP) :: total
178 !---------------------------------------------------------------------------
179
180 call calculate_statval( field, field%mesh%lcmesh_list(:), & ! (in)
181 statval, total ) ! (out)
182
183 call statistics_total_core( &
184 field%varname, statval, total, & ! (in)
185 log_suppress, global, & ! (in)
186 mean, sum ) ! (out)
187
188 return
189 end subroutine meshfield_statistics_total_2d
190
191 !-----------------------------------------------------------------------------
193!OCL SERIAL
194 subroutine meshfield_statistics_total_3d( field, & ! (in)
195 log_suppress, global, & ! (in)
196 mean, sum ) ! (out)
197
198 implicit none
199
200 class(meshfield3d), intent(in) :: field
201 logical, intent(in), optional :: log_suppress
202 logical, intent(in), optional :: global
203 real(RP), intent(out), optional :: mean
204 real(DP), intent(out), optional :: sum
205
206 real(DP) :: statval
207 real(DP) :: total
208 !---------------------------------------------------------------------------
209
210 call calculate_statval( field, field%mesh%lcmesh_list(:), & ! (in)
211 statval, total ) ! (out)
212
213 call statistics_total_core( &
214 field%varname, statval, total, & ! (in)
215 log_suppress, global, & ! (in)
216 mean, sum ) ! (out)
217
218 return
219 end subroutine meshfield_statistics_total_3d
220
221
222 !-----------------------------------------------------------------------------
224 !-----------------------------------------------------------------------------
225
227 subroutine meshfield_statistics_maxmin_1d( field_list, & ! (in)
228 log_suppress, global, & ! (in)
229 maxval, minval ) ! (out)
230
231 implicit none
232
233 class(meshfield1d), intent(in) :: field_list(:)
234 logical, intent(in), optional :: log_suppress
235 logical, intent(in), optional :: global
236 real(RP), intent(out), optional :: maxval(size(field_list))
237 real(RP), intent(out), optional :: minval(size(field_list))
238
239 real(RP) :: statval_l( size(field_list),2)
240 integer :: statidx_l(3,size(field_list),2)
241 integer :: VA
242 !---------------------------------------------------------------------------
243
244 va = size(field_list)
245 call search_maxmin_local( &
246 va, field_list, field_list(1)%mesh%lcmesh_list, & ! (in)
247 statval_l, statidx_l ) ! (out)
248
249 call statistics_detail_core( &
250 va, field_list, statval_l, statidx_l, .not. global, log_suppress, & ! (in)
251 maxval, minval ) ! (out)
252
253 return
254 end subroutine meshfield_statistics_maxmin_1d
255
256!OCL SERIAL
257 subroutine meshfield_statistics_detail_1d( field_list, local )
258 implicit none
259
260 class(meshfield1d), intent(in) :: field_list(:)
261 logical, intent(in), optional :: local
262
263 real(RP) :: statval_l( size(field_list),2)
264 integer :: statidx_l(3,size(field_list),2)
265 integer :: VA
266 !---------------------------------------------------------------------------
267
268 log_newline
269 log_info("MeshField_STATISTICS_detail_1D",*) 'Variable Statistics '
270
271 va = size(field_list)
272 call search_maxmin_local( &
273 va, field_list, field_list(1)%mesh%lcmesh_list, & ! (in)
274 statval_l, statidx_l ) ! (out)
275
276 call statistics_detail_core( &
277 va, field_list, statval_l, statidx_l, local ) ! (in)
278
279 log_newline
280
281 return
282 end subroutine meshfield_statistics_detail_1d
283
284 !-----------------------------------------------------------------------------
286 subroutine meshfield_statistics_maxmin_2d( field_list, & ! (in)
287 log_suppress, global, & ! (in)
288 maxval, minval ) ! (out)
289
290 implicit none
291
292 class(meshfield2d), intent(in) :: field_list(:)
293 logical, intent(in), optional :: log_suppress
294 logical, intent(in), optional :: global
295 real(RP), intent(out), optional :: maxval(size(field_list))
296 real(RP), intent(out), optional :: minval(size(field_list))
297
298 real(RP) :: statval_l( size(field_list),2)
299 integer :: statidx_l(3,size(field_list),2)
300 integer :: VA
301 !---------------------------------------------------------------------------
302
303 va = size(field_list)
304 call search_maxmin_local( &
305 va, field_list, field_list(1)%mesh%lcmesh_list, & ! (in)
306 statval_l, statidx_l ) ! (out)
307
308 call statistics_detail_core( &
309 va, field_list, statval_l, statidx_l, .not. global, log_suppress, & ! (in)
310 maxval, minval ) ! (out)
311
312 return
313 end subroutine meshfield_statistics_maxmin_2d
314
315!OCL SERIAL
316 subroutine meshfield_statistics_detail_2d( field_list, local )
317 implicit none
318
319 class(meshfield2d), intent(in) :: field_list(:)
320 logical, intent(in), optional :: local
321
322 real(RP) :: statval_l( size(field_list),2)
323 integer :: statidx_l(3,size(field_list),2)
324 integer :: VA
325 !---------------------------------------------------------------------------
326
327 log_newline
328 log_info("MeshField_STATISTICS_detail_2D",*) 'Variable Statistics '
329
330 va = size(field_list)
331 call search_maxmin_local( &
332 va, field_list, field_list(1)%mesh%lcmesh_list, & ! (in)
333 statval_l, statidx_l ) ! (out)
334
335 call statistics_detail_core( &
336 va, field_list, statval_l, statidx_l, local ) ! (in)
337
338 log_newline
339
340 return
341 end subroutine meshfield_statistics_detail_2d
342
343 !-----------------------------------------------------------------------------
345!OCL SERIAL
346 subroutine meshfield_statistics_maxmin_3d( field_list, & ! (in)
347 log_suppress, global, & ! (in)
348 maxval, minval ) ! (out)
349
350 implicit none
351
352 class(meshfield3d), intent(in) :: field_list(:)
353 logical, intent(in), optional :: log_suppress
354 logical, intent(in), optional :: global
355 real(RP), intent(out), optional :: maxval(size(field_list))
356 real(DP), intent(out), optional :: minval(size(field_list))
357
358 real(RP) :: statval_l( size(field_list),2)
359 integer :: statidx_l(3,size(field_list),2)
360 integer :: VA
361 !---------------------------------------------------------------------------
362
363 va = size(field_list)
364 call search_maxmin_local( &
365 va, field_list, field_list(1)%mesh%lcmesh_list, & ! (in)
366 statval_l, statidx_l ) ! (out)
367
368 call statistics_detail_core( &
369 va, field_list, statval_l, statidx_l, .not. global, log_suppress, & ! (in)
370 maxval, minval ) ! (out)
371
372 return
373 end subroutine meshfield_statistics_maxmin_3d
374
375!OCL SERIAL
376 subroutine meshfield_statistics_detail_3d( field_list, local )
377 implicit none
378
379 class(meshfield3d), intent(in) :: field_list(:)
380 logical, intent(in), optional :: local
381
382 real(RP) :: statval_l( size(field_list),2)
383 integer :: statidx_l(3,size(field_list),2)
384 integer :: VA
385 !---------------------------------------------------------------------------
386
387 log_newline
388 log_info("MeshField_STATISTICS_detail_3D",*) 'Variable Statistics '
389
390 va = size(field_list)
391 call search_maxmin_local( &
392 va, field_list, field_list(1)%mesh%lcmesh_list, & ! (in)
393 statval_l, statidx_l ) ! (out)
394
395 call statistics_detail_core( &
396 va, field_list, statval_l, statidx_l, local ) ! (in)
397
398 log_newline
399
400 return
401 end subroutine meshfield_statistics_detail_3d
402
403!-- private --------------------------------------------------------------------
404
405!OCL SERIAL
406 subroutine search_maxmin_local( VA, field_list, lcmesh_list, &
407 statval_l, statidx_l )
408
409 implicit none
410 integer, intent(in) :: VA
411 class(meshfieldbase), intent(in) :: field_list(VA)
412 class(localmeshbase), intent(in), target :: lcmesh_list(:)
413 real(RP), intent(out) :: statval_l ( VA,2)
414 integer, intent(out) :: statidx_l (3,VA,2)
415
416 type(localmeshbase), pointer :: lcmesh
417 class(localmeshfieldbase), pointer :: lcfield
418 class(elementbase), pointer :: refElem
419
420 integer :: v
421 integer :: n
422 integer :: ke, p
423 !---------------------------------------------------------------------------
424
425 do n=1, size(lcmesh_list)
426 lcmesh => lcmesh_list(n)
427 refelem => lcmesh%refElem
428 do v=1, va
429 call field_list(v)%GetLocalMeshField(n, lcfield)
430 if (n==1) then
431 statval_l( v,:) = lcfield%val(1,lcmesh%NeS)
432 statidx_l(1,v,:) = 1
433 statidx_l(2,v,:) = lcmesh%NeS
434 statidx_l(3,v,:) = n
435 end if
436 do ke=lcmesh%NeS, lcmesh%NeE
437 do p =1, refelem%Np
438 if ( lcfield%val(p,ke) > statval_l(v,1) ) then
439 statval_l( v,1) = lcfield%val(p,ke)
440 statidx_l(:,v,1) = (/ p, ke, n /)
441 end if
442 if ( lcfield%val(p,ke) < statval_l(v,2) ) then
443 statval_l( v,2) = lcfield%val(p,ke)
444 statidx_l(:,v,2) = (/ p, ke, n /)
445 end if
446 end do
447 end do
448 end do
449 end do
450
451 return
452 end subroutine search_maxmin_local
453
454!OCL SERIAL
455 subroutine statistics_detail_core( &
456 VA, field_list, statval_l, statidx_l, & ! (in)
457 local, log_supress, & ! (in)
458 maxval, minval ) ! (out)
459
460 use scale_prc, only: &
461 prc_nprocs
462 implicit none
463
464 integer, intent(in) :: VA
465 class(meshfieldbase), intent(in) :: field_list(VA)
466 real(DP), intent(in) :: statval_l( VA,2)
467 integer, intent(in) :: statidx_l(3,VA,2)
468 logical, intent(in), optional :: local
469 logical, intent(in), optional :: log_supress
470 real(RP), intent(out), optional :: maxval(VA)
471 real(RP), intent(out), optional :: minval(VA)
472
473 real(RP) :: statval ( VA,2,0:PRC_nprocs-1)
474 integer :: statidx (3,VA,2,0:PRC_nprocs-1)
475 real(RP) :: allstatval(VA,2)
476 integer :: allstatidx(VA,2)
477 logical :: do_globalcomm
478
479 integer :: v
480 integer :: p
481 integer :: ierr
482
483 logical :: supress_
484 !---------------------------------------------------------------------------
485
486 do_globalcomm = base%use_globalcomm
487 if ( present(local) ) do_globalcomm = ( .not. local )
488
489 if ( present(log_supress) ) then
490 supress_ = log_supress
491 else
492 supress_ = .false.
493 end if
494
495 if ( do_globalcomm ) then
496 call prof_rapstart('COMM_Bcast', 2)
497
498 call mpi_allgather( statval_l(:,:), &
499 va*2, &
500 base%comm_datatype, &
501 statval(:,:,:), &
502 va*2, &
503 base%comm_datatype, &
504 prc_local_comm_world, &
505 ierr )
506
507 call mpi_allgather( statidx_l(:,:,:), &
508 3*va*2, &
509 mpi_integer, &
510 statidx(:,:,:,:), &
511 3*va*2, &
512 mpi_integer, &
513 prc_local_comm_world, &
514 ierr )
515
516 call prof_rapend ('COMM_Bcast', 2)
517
518 do v = 1, va
519 allstatval(v,1) = statval(v,1,0)
520 allstatval(v,2) = statval(v,2,0)
521 allstatidx(v,:) = 0
522 do p = 1, prc_nprocs-1
523 if ( statval(v,1,p) > allstatval(v,1) ) then
524 allstatval(v,1) = statval(v,1,p)
525 allstatidx(v,1) = p
526 end if
527 if ( statval(v,2,p) < allstatval(v,2) ) then
528 allstatval(v,2) = statval(v,2,p)
529 allstatidx(v,2) = p
530 end if
531 end do
532 if ( .not. supress_ ) then
533 log_info_cont(*) '[', trim(field_list(v)%varname), ']'
534 log_info_cont('(1x,A,ES17.10,A,4(I5,A))') ' MAX =', &
535 allstatval(v,1), ' (rank=', &
536 allstatidx(v,1), '; ', &
537 statidx(1,v,1,allstatidx(v,1)),',', &
538 statidx(2,v,1,allstatidx(v,1)),',', &
539 statidx(3,v,1,allstatidx(v,1)),')'
540 log_info_cont('(1x,A,ES17.10,A,4(I5,A))') ' MIN =', &
541 allstatval(v,2), ' (rank=', &
542 allstatidx(v,2), '; ', &
543 statidx(1,v,2,allstatidx(v,2)),',', &
544 statidx(2,v,2,allstatidx(v,2)),',', &
545 statidx(3,v,2,allstatidx(v,2)),')'
546 end if
547
548 if ( present(maxval) ) maxval(v) = allstatval(v,1)
549 if ( present(minval) ) minval(v) = allstatval(v,2)
550 enddo
551 else
552 ! statistics on each node
553 do v = 1, va
554 if ( .not. supress_ ) then
555 log_info_cont(*) '[', trim(field_list(v)%varname), ']'
556 log_info_cont('(1x,A,ES17.10,A,3(I5,A))') 'MAX = ', &
557 statval_l( v,1),' (', &
558 statidx_l(1,v,1),',', &
559 statidx_l(2,v,1),',', &
560 statidx_l(3,v,1),')'
561 log_info_cont('(1x,A,ES17.10,A,3(I5,A))') 'MIN = ', &
562 statval_l( v,2),' (', &
563 statidx_l(1,v,2),',', &
564 statidx_l(2,v,2),',', &
565 statidx_l(3,v,2),')'
566 end if
567 if ( present(maxval) ) maxval(v) = statval_l(v,1)
568 if ( present(minval) ) minval(v) = statval_l(v,2)
569 enddo
570 endif
571
572 return
573 end subroutine statistics_detail_core
574
575 !-----
576
577!OCL SERIAL
578 subroutine calculate_statval( field, lcmesh_list, &
579 statval, total )
580
581 implicit none
582
583 class(meshfieldbase), intent(in) :: field
584 class(localmeshbase), intent(in), target :: lcmesh_list(:)
585 real(DP), intent(out) :: statval
586 real(DP), intent(out) :: total
587
588 type(localmeshbase), pointer :: lcmesh
589 class(localmeshfieldbase), pointer :: lcfield
590 class(elementbase), pointer :: refElem
591
592 real(DP) :: statval_lc
593 real(DP) :: total_lc
594 real(DP) :: weight
595
596 integer :: n
597 integer :: ke, p
598 !---------------------------------------------------------------------------
599
600 statval = 0.0_dp
601 total = 0.0_dp
602 do n=1, size(lcmesh_list)
603 lcmesh => lcmesh_list(n)
604 refelem => lcmesh%refElem
605 call field%GetLocalMeshField(n, lcfield)
606
607 total_lc = 0.0_dp
608 statval_lc = 0.0_dp
609 !$omp parallel do private(p, weight) reduction(+:total_lc, statval_lc)
610 do ke=lcmesh%NeS, lcmesh%NeE
611 do p=1, refelem%Np
612 weight = lcmesh%J(p,ke) * lcmesh%Gsqrt(p,ke) * refelem%IntWeight_lgl(p)
613 total_lc = total_lc + weight
614 statval_lc = statval_lc + weight * lcfield%val(p,ke)
615 end do
616 end do
617 total = total + total_lc
618 statval = statval + statval_lc
619 end do
620
621 return
622 end subroutine calculate_statval
623
624!OCL SERIAL
625 subroutine statistics_total_core( &
626 varname, statval, total, &
627 log_suppress, &
628 global, &
629 mean, sum )
630 use scale_prc, only: &
631 prc_myrank, &
632 prc_abort
633 use scale_const, only: &
634 eps => const_eps, &
635 undef => const_undef
636
637 implicit none
638
639 character(len=*), intent(in) :: varname
640 real(DP), intent(in) :: statval
641 real(DP), intent(in) :: total
642 logical, intent(in), optional :: log_suppress
643 logical, intent(in), optional :: global
644 real(RP), intent(out), optional :: mean
645 real(DP), intent(out), optional :: sum
646
647 real(DP) :: sendbuf(2), recvbuf(2)
648 real(DP) :: sum_, mean_
649
650 logical :: suppress_, global_
651 integer :: ierr
652 !---------------------------------------------------------------------------
653
654 if ( .NOT. ( statval > -1.0_dp .OR. statval < 1.0_dp ) ) then ! must be NaN
655 log_error("MeshField_STATISTICS_total",*) 'NaN is detected for ', trim(varname), ' in rank ', prc_myrank
656 call prc_abort
657 endif
658
659 if ( present(log_suppress) ) then
660 suppress_ = log_suppress
661 else
662 suppress_ = .false.
663 end if
664
665 if ( present(global) ) then
666 global_ = global
667 else
668 global_ = base%use_globalcomm
669 end if
670
671 if ( global_ ) then
672 call prof_rapstart('COMM_Allreduce', 2)
673 sendbuf(1) = statval
674 sendbuf(2) = total
675 ! All reduce
676 call mpi_allreduce( sendbuf(:), recvbuf(:), &
677 2, &
678 mpi_double_precision, &
679 mpi_sum, &
680 prc_local_comm_world, &
681 ierr )
682 call prof_rapend ('COMM_Allreduce', 2)
683
684 if ( recvbuf(2) < eps ) then
685 sum_ = undef
686 mean_ = undef
687 else
688 sum_ = recvbuf(1)
689 mean_ = recvbuf(1) / recvbuf(2)
690 end if
691 ! statistics over the all node
692 if ( .not. suppress_ ) then ! if varname is empty, suppress output
693 log_info("MeshField_STATISTICS_total_3D",'(1x,A,A24,A,ES24.17)') &
694 '[', trim(varname), '] MEAN(global) = ', mean_
695 endif
696 else
697 if ( total < eps ) then
698 sum_ = undef
699 mean_ = undef
700 else
701 sum_ = statval
702 mean_ = statval / total
703 end if
704
705 ! statistics on each node
706 if ( .not. suppress_ ) then ! if varname is empty, suppress output
707 log_info("MeshField_STATISTICS_total_3D",'(1x,A,A24,A,ES24.17)') &
708 '[', trim(varname), '] MEAN(local) = ', mean_
709 endif
710 endif
711
712 if ( present(mean) ) mean = mean_
713 if ( present(sum ) ) sum = sum_
714
715 return
716 end subroutine statistics_total_core
717
module FElib / Element / Base
module FElib / Mesh / Local, Base
module FElib / Data / base
module FElib / Data / Statistics
subroutine, public meshfield_statistics_setup()
Setup.