1
1
program test_moment
2
2
use stdlib_experimental_error, only: assert
3
3
use stdlib_experimental_kinds, only: sp, dp, int32, int64
4
- use stdlib_experimental_stats, only: mean, moment, var
4
+ use stdlib_experimental_stats, only: moment
5
5
implicit none
6
6
7
7
@@ -19,20 +19,19 @@ program test_moment
19
19
cmplx (1.26401_sp , 0.00000_sp ),&
20
20
cmplx (0.00000_sp , 0.88833_sp ),&
21
21
cmplx (1.14352_sp , 0.00000_sp )]
22
- complex (dp) :: cd1(5 ) = [ cmplx (0.57706_dp , 0.00000_dp ),&
23
- cmplx (0.00000_dp , 1.44065_dp ),&
24
- cmplx (1.26401_dp , 0.00000_dp ),&
25
- cmplx (0.00000_dp , 0.88833_dp ),&
26
- cmplx (1.14352_dp , 0.00000_dp )]
27
22
complex (sp) :: cs(5 ,3 )
28
- complex (dp) :: cd(5 ,3 )
29
23
30
24
31
25
call test_sp(real (d1,sp), real (d,sp))
32
26
call test_dp(d1, d)
33
27
call test_int32(int (d1, int32), int (d, int32))
34
28
call test_int64(int (d1, int64), int (d, int64))
35
29
30
+ cs(:,1 ) = cs1
31
+ cs(:,2 ) = cs1* 3_sp
32
+ cs(:,3 ) = cs1* 1.5_sp
33
+ call test_csp(cs1, cs)
34
+
36
35
37
36
contains
38
37
subroutine test_sp (x1 , x2 )
@@ -161,17 +160,17 @@ subroutine test_sp(x1, x2)
161
160
162
161
print * ,' test_sp_3dim_mask_array' , order
163
162
call assert( abs (moment(x3, order, x3 < 11 ) - 7.7370242214532876_dp ) < sptol)
164
- call assert( all ( abs ( var (x3, 1 , x3 < 45 , corrected = .false. ) - &
163
+ call assert( all ( abs ( moment (x3, order, 1 , x3 < 45 ) - &
165
164
reshape ([5._sp , 5._sp , 1.25_sp , 20._sp , 20._sp , 5._sp ,&
166
165
80._sp , 80._sp , 32._sp / 3 .],&
167
166
[size (x3, 2 ), size (x3, 3 )])) < sptol ))
168
- call assert( all ( abs ( var (x3, 2 , x3 < 45 , corrected = .false. ) - &
167
+ call assert( all ( abs ( moment (x3, order, 2 , x3 < 45 ) - &
169
168
reshape ([ 38._sp / 3 ., 86._sp / 9 ., 62._sp / 9 ., 14._sp / 3 ., 152._sp / 3 .,&
170
169
344._sp / 9 ., 248._sp / 9 ., 168._sp / 9 ., 1824._sp / 9 .,&
171
170
1376._sp / 9 ., 992._sp / 9 ., 4._sp &
172
171
],&
173
172
[size (x3, 1 ), size (x3, 3 )])) < sptol ))
174
- call assert( all ( abs ( var (x3, 3 , x3 < 45 , corrected = .false. ) - &
173
+ call assert( all ( abs ( moment (x3, order, 3 , x3 < 45 ) - &
175
174
reshape ([14._sp / 9 ., 14._sp , 350._sp / 9 ., 686._sp / 9 ., 56._sp / 9 .,&
176
175
224._sp / 9 ., 56._sp , 896._sp / 9 ., 126._sp , 1400._sp / 9 .,&
177
176
1694._sp / 9 ., 36._sp &
@@ -307,17 +306,17 @@ subroutine test_dp(x1, x2)
307
306
308
307
print * ,' test_dp_3dim_mask_array' , order
309
308
call assert( abs (moment(x3, order, x3 < 11 ) - 7.7370242214532876_dp ) < dptol)
310
- call assert( all ( abs ( var (x3, 1 , x3 < 45 , corrected = .false. ) - &
309
+ call assert( all ( abs ( moment (x3, order, 1 , x3 < 45 ) - &
311
310
reshape ([5._dp , 5._dp , 1.25_dp , 20._dp , 20._dp , 5._dp ,&
312
311
80._dp , 80._dp , 32._dp / 3 .],&
313
312
[size (x3, 2 ), size (x3, 3 )])) < dptol ))
314
- call assert( all ( abs ( var (x3, 2 , x3 < 45 , corrected = .false. ) - &
313
+ call assert( all ( abs ( moment (x3, order, 2 , x3 < 45 ) - &
315
314
reshape ([ 38._dp / 3 ., 86._dp / 9 ., 62._dp / 9 ., 14._dp / 3 ., 152._dp / 3 .,&
316
315
344._dp / 9 ., 248._dp / 9 ., 168._dp / 9 ., 1824._dp / 9 .,&
317
316
1376._dp / 9 ., 992._dp / 9 ., 4._dp &
318
317
],&
319
318
[size (x3, 1 ), size (x3, 3 )])) < dptol ))
320
- call assert( all ( abs ( var (x3, 3 , x3 < 45 , corrected = .false. ) - &
319
+ call assert( all ( abs ( moment (x3, order, 3 , x3 < 45 ) - &
321
320
reshape ([14._dp / 9 ., 14._dp , 350._dp / 9 ., 686._dp / 9 ., 56._dp / 9 .,&
322
321
224._dp / 9 ., 56._dp , 896._dp / 9 ., 126._dp , 1400._dp / 9 .,&
323
322
1694._dp / 9 ., 36._dp &
@@ -453,17 +452,17 @@ subroutine test_int32(x1, x2)
453
452
454
453
print * ,' test_dp_3dim_mask_array' , order
455
454
call assert( abs (moment(x3, order, x3 < 11 ) - 7.7370242214532876_dp ) < dptol)
456
- call assert( all ( abs ( var (x3, 1 , x3 < 45 , corrected = .false. ) - &
455
+ call assert( all ( abs ( moment (x3, order, 1 , x3 < 45 ) - &
457
456
reshape ([5._dp , 5._dp , 1.25_dp , 20._dp , 20._dp , 5._dp ,&
458
457
80._dp , 80._dp , 32._dp / 3 .],&
459
458
[size (x3, 2 ), size (x3, 3 )])) < dptol ))
460
- call assert( all ( abs ( var (x3, 2 , x3 < 45 , corrected = .false. ) - &
459
+ call assert( all ( abs ( moment (x3, order, 2 , x3 < 45 ) - &
461
460
reshape ([ 38._dp / 3 ., 86._dp / 9 ., 62._dp / 9 ., 14._dp / 3 ., 152._dp / 3 .,&
462
461
344._dp / 9 ., 248._dp / 9 ., 168._dp / 9 ., 1824._dp / 9 .,&
463
462
1376._dp / 9 ., 992._dp / 9 ., 4._dp &
464
463
],&
465
464
[size (x3, 1 ), size (x3, 3 )])) < dptol ))
466
- call assert( all ( abs ( var (x3, 3 , x3 < 45 , corrected = .false. ) - &
465
+ call assert( all ( abs ( moment (x3, order, 3 , x3 < 45 ) - &
467
466
reshape ([14._dp / 9 ., 14._dp , 350._dp / 9 ., 686._dp / 9 ., 56._dp / 9 .,&
468
467
224._dp / 9 ., 56._dp , 896._dp / 9 ., 126._dp , 1400._dp / 9 .,&
469
468
1694._dp / 9 ., 36._dp &
@@ -599,22 +598,109 @@ subroutine test_int64(x1, x2)
599
598
600
599
print * ,' test_dp_3dim_mask_array' , order
601
600
call assert( abs (moment(x3, order, x3 < 11 ) - 7.7370242214532876_dp ) < dptol)
602
- call assert( all ( abs ( var (x3, 1 , x3 < 45 , corrected = .false. ) - &
601
+ call assert( all ( abs ( moment (x3, order, 1 , x3 < 45 ) - &
603
602
reshape ([5._dp , 5._dp , 1.25_dp , 20._dp , 20._dp , 5._dp ,&
604
603
80._dp , 80._dp , 32._dp / 3 .],&
605
604
[size (x3, 2 ), size (x3, 3 )])) < dptol ))
606
- call assert( all ( abs ( var (x3, 2 , x3 < 45 , corrected = .false. ) - &
605
+ call assert( all ( abs ( moment (x3, order, 2 , x3 < 45 ) - &
607
606
reshape ([ 38._dp / 3 ., 86._dp / 9 ., 62._dp / 9 ., 14._dp / 3 ., 152._dp / 3 .,&
608
607
344._dp / 9 ., 248._dp / 9 ., 168._dp / 9 ., 1824._dp / 9 .,&
609
608
1376._dp / 9 ., 992._dp / 9 ., 4._dp &
610
609
],&
611
610
[size (x3, 1 ), size (x3, 3 )])) < dptol ))
612
- call assert( all ( abs ( var (x3, 3 , x3 < 45 , corrected = .false. ) - &
611
+ call assert( all ( abs ( moment (x3, order, 3 , x3 < 45 ) - &
613
612
reshape ([14._dp / 9 ., 14._dp , 350._dp / 9 ., 686._dp / 9 ., 56._dp / 9 .,&
614
613
224._dp / 9 ., 56._dp , 896._dp / 9 ., 126._dp , 1400._dp / 9 .,&
615
614
1694._dp / 9 ., 36._dp &
616
615
], [size (x3,1 ), size (x3,2 )] ))&
617
616
< dptol ))
618
617
619
618
end subroutine
619
+
620
+ subroutine test_csp (x1 , x2 )
621
+ complex (sp), intent (in ) :: x1(:), x2(:, :)
622
+
623
+ integer :: order
624
+ complex (sp), allocatable :: x3(:, :, :)
625
+
626
+ order = 1
627
+
628
+ ! 1dim
629
+ print * ,' test_sp_1dim' , order
630
+ call assert( abs (moment(x1, order)) < sptol)
631
+ call assert( abs (moment(x1, order, dim= 1 )) < sptol)
632
+
633
+ print * ,' test_sp_1dim_mask' , order
634
+ call assert( isnan(abs (moment(x1, order, .false. ))))
635
+ call assert( isnan(abs (moment(x1, order, 1 , .false. ))))
636
+
637
+ print * ,' test_sp_1dim_mask_array' , order
638
+ call assert( abs (moment(x1, order, aimag (x1) == 0 )) < sptol)
639
+ call assert( abs (moment(x1, order, 1 , aimag (x1) == 0 )) < sptol)
640
+
641
+ ! 2dim
642
+ print * ,' test_sp_2dim' , order
643
+ call assert( abs (moment(x2, order)) < sptol)
644
+ call assert( all ( abs ( moment(x2, order, 1 )) < sptol))
645
+ call assert( all ( abs ( moment(x2, order, 2 )) < sptol))
646
+
647
+ print * ,' test_sp_2dim_mask' , order
648
+ call assert( isnan(abs (moment(x2, order, .false. ))))
649
+ call assert( any (isnan(abs (moment(x2, order, 1 , .false. )))))
650
+ call assert( any (isnan(abs (moment(x2, order, 2 , .false. )))))
651
+
652
+ print * ,' test_sp_2dim_mask_array' , order
653
+ call assert( abs (moment(x2, order, aimag (x2) == 0 )) < sptol)
654
+ call assert( all ( abs ( moment(x2, order, 1 , aimag (x2) == 0 )) < sptol))
655
+ call assert( any (isnan( abs ( moment(x2, order, 2 , aimag (x2) == 0 )))))
656
+
657
+ order = 2
658
+
659
+ ! 1dim
660
+ print * ,' test_sp_1dim' , order
661
+ call assert( abs (moment(x1, order) - (- 6.459422410E-02 ,- 0.556084037 )) < sptol)
662
+ call assert( abs (moment(x1, order, dim= 1 ) - &
663
+ (- 6.459422410E-02 ,- 0.556084037 )) < sptol)
664
+
665
+ print * ,' test_sp_1dim_mask' , order
666
+ call assert( isnan(abs (moment(x1, order, .false. ))))
667
+ call assert( isnan(abs (moment(x1, order, 1 , .false. ))))
668
+
669
+ print * ,' test_sp_1dim_mask_array' , order
670
+ call assert( abs (moment(x1, order, aimag (x1) == 0 ) - &
671
+ (8.969944715E-02 ,0.00000000 )) < sptol)
672
+ call assert( abs (moment(x1, order, 1 , aimag (x1) == 0 ) - &
673
+ (8.969944715E-02 ,0.00000000 )) < sptol)
674
+
675
+ ! 2dim
676
+ print * ,' test_sp_2dim' , order
677
+ call assert( abs (moment(x2, order) - (- 0.163121477 ,- 1.86906016 )) < sptol)
678
+ call assert( all ( abs ( moment(x2, order, 1 ) - &
679
+ [(- 6.459422410E-02 ,- 0.556084037 ),&
680
+ (- 0.581347823 ,- 5.00475645 ),&
681
+ (- 0.145336956 ,- 1.25118911 )]&
682
+ ) < sptol))
683
+ call assert( all ( abs ( moment(x2, order, 2 ) - &
684
+ [(0.240498722 ,0.00000000 ),&
685
+ (- 1.49895227 ,0.00000000 ),&
686
+ (1.15390968 ,0.00000000 ),&
687
+ (- 0.569927275 ,0.00000000 ),&
688
+ (0.944405317 ,0.00000000 )]&
689
+ ) < sptol))
690
+
691
+ print * ,' test_sp_2dim_mask' , order
692
+ call assert( isnan(abs (moment(x2, order, .false. ))))
693
+ call assert( any (isnan(abs (moment(x2, order, 1 , .false. )))))
694
+ call assert( any (isnan(abs (moment(x2, order, 2 , .false. )))))
695
+
696
+ print * ,' test_sp_2dim_mask_array' , order
697
+ call assert( abs (moment(x2, order, aimag (x2) == 0 )- &
698
+ (1.08109438 ,0.00000000 )) < sptol)
699
+ call assert( all ( abs ( moment(x2, order, 1 , aimag (x2)==0 ) - &
700
+ [(8.969944715E-02 ,0.00000000 ),&
701
+ (0.807295084 ,0.00000000 ),&
702
+ (0.201823771 ,0.00000000 )]&
703
+ ) < sptol))
704
+
705
+ end subroutine
620
706
end program
0 commit comments