function dirichlet_sample_test ( ) %*****************************************************************************80 % %% DIRICHLET_SAMPLE_TEST tests DIRICHLET_SAMPLE. % % Licensing: % % This code is distributed under the GNU LGPL license. % % Modified: % % 13 April 2016 % % Author: % % John Burkardt % n = 3; nsample = 1000; seed = 123456789; fprintf ( 1, '\n' ); fprintf ( 1, 'DIRICHLET_SAMPLE_TEST\n' ); fprintf ( 1, ' DIRICHLET_SAMPLE samples the Dirichlet distribution;\n' ); fprintf ( 1, ' DIRICHLET_MEAN computes the Dirichlet mean;\n' ); fprintf ( 1, ' DIRICHLET_VARIANCE computes the Dirichlet variance.\n' ); a(1:n) = [ 0.250, 0.500, 1.250 ]; check = dirichlet_check ( n, a ); if ( ~ check ); fprintf ( 1, '\n' ); fprintf ( 1, 'DIRICHLET_SAMPLE_TEST - Fatal error!\n' ); fprintf ( 1, ' The parameters are not legal.\n' ); return end fprintf ( 1, '\n' ); fprintf ( 1, ' Number of components N = %6d\n', n ); r8vec_print ( n, a, ' PDF parameters A:' ); mean(1:n) = dirichlet_mean ( n, a ); variance(1:n) = dirichlet_variance ( n, a ); fprintf ( 1, '\n' ); fprintf ( 1, ' PDF mean, variance:\n' ); fprintf ( 1, '\n' ); for i = 1 : n fprintf ( 1, ' %6d %14f %14f\n', i, mean(i), variance(i) ); end m2(1:n,1:n) = dirichlet_moment2 ( n, a ); r8mat_print ( n, n, m2, ' Second moment matrix:' ); for i = 1 : nsample [ x(1:n,i), seed ] = dirichlet_sample ( n, a, seed ); end xmax(1:n) = r8row_max ( n, nsample, x ); xmin(1:n) = r8row_min ( n, nsample, x ); mean(1:n) = r8row_mean ( n, nsample, x ); variance(1:n) = r8row_variance ( n, nsample, x ); fprintf ( 1, '\n' ); fprintf ( 1, ' Sample size = %d', nsample ) fprintf ( 1, '\n' ); fprintf ( 1, ' Observed Min, Max, Mean, Variance:\n' ); fprintf ( 1, '\n' ); for i = 1 : n fprintf ( 1, ' %6d %14g %14g %14g %14g\n', ... i, xmin(i), xmax(i), mean(i), variance(i) ); end return end