App-boxmuller

 view release on metacpan or  search on metacpan

scripts/boxmuller  view on Meta::CPAN


   sub LLN ( $ ) ; * LLN = * looks_like_number ; # 関数名が長すぎるので、短くした。
   sub printErr( $ ){ print STDERR BRIGHT_RED "Option -$_[0] should have a numeric specification.\n" ; exit 1 }
   $mu = $o{m} ? LLN $o{m} ? $o{m} : printErr "m" : 0 ;  #m : 平均
   $sd = $o{d} ? LLN $o{d} ? $o{d} : printErr "d" : $o{v} ? LLN $o{v} ? sqrt  $o{v} : printErr "v" : 1 ;  #sd:標準偏差 
}

sub main ( ) {  #  乱数の出力
   sub getrand  ;
   sub boxmuller ( $$ ) ;  #ボックスミュラー法によるガウス乱数の作成
   * getrand = * boxmuller ; 
   * getrand = * lognormal if $o{L} ;  # 対数正規分布の指定があった場合。

   while ( $count < $upto ) { 
   	   sleep $o{'@'} if defined $o{'@'} ;
       my $x = getrand $mu, $sd ; 
       $x = sprintf "% .$o{'.'}f" , $x if defined $o{'.'} ; # <-- May be efficientized. 
       $s1 += $x ; 
       $s2 += $x ** 2 ; 
       $count ++ ; # 出力個数を計数するので個々を選んだ。
       print "$count\t" if $o{':'} ;  # <-- Maybe effiecientized by other code structure.
       print "$x\n" ; 
   }

  sub boxmuller ( $$ ) {  #  ガウス乱数の生成
      state $piW = atan2(1,1)* 8  ;  #6.28=3.14*2
      state $z = undef ; 
      do { my $t = $z ; undef $z ; return $t * $_[1] + $_[0] }  if defined $z ; 
      
      my $r1 = 1 - rand ; # ; $r1 = rand until $r1 ; 値が0になることを阻止。
      my $r2 = $piW * rand ;
      my $r1R = sqrt ( -2 * log $r1 ); 
      my $r2S = sin $r2 ; 
      my $r2C = cos $r2 ; 
      ( my $t , $z ) = ( $r1R * $r2S , $r1R * $r2C ) ; 
      return $t * $_[1] + $_[0] ; 
  }

  sub lognormal ( $$ ) { 
  	return exp boxmuller $_[0], $_[1] ;
  }
}

sub SecondInfo( ) {   #  処理したことについての二次情報を出力
    return if $o{1} ;
    use FindBin qw [ $Script ] ; 
    my $cmd = "$Script -m $mu -d $sd" ; 
    $cmd .= ' -L' if $o{L} ;
    print STDERR 
       CYAN "printed lines: ", BRIGHT_CYAN $count ,
       CYAN " , used random seed: " , BRIGHT_CYAN  $o{s} ,
       CYAN " , elapsed seconds: " , BRIGHT_CYAN  tv_interval ($time0) ,
       RESET "\n" , 
       CYAN "sum = " , BRIGHT_CYAN  sprintf("%g", $s1 ) ,
       CYAN " , squared sum = " , BRIGHT_CYAN  sprintf( "%g" , $s2 ) , 
       CYAN " ($cmd) " , RESET "\n" ;
 }

## ヘルプとバージョン情報
BEGIN {
  our $VERSION = 0.24 ;
  $Getopt::Std::STANDARD_HELP_VERSION = 1 ; 
  grep { m/--help/} @ARGV and *VERSION_MESSAGE = sub {} ; 
    # 最初は 0.21 を目安とする。
    # 1.00 以上とする必要条件は英語版のヘルプをきちんと出すこと。
    # 2.00 以上とする必要条件はテストコードが含むこと。
    # 0.22 : -g inf を指定可能とした。
    # 0.23 : 説明文を分かり安くした。
    # 0.24 : 英語のマニュアルをPOD化する。
}  

sub HELP_MESSAGE {
    use FindBin qw[ $Script $Bin ] ;
    sub EnvJ ( ) { $ENV{LANG} =~ m/^ja_JP/ ? 1 : 0 } ; # # ja_JP.UTF-8 
    sub en( ) { grep ( /^en(g(i(sh?)?)?)?/i , @ARGV ) ? 1 : 0 } # English という文字列を先頭から2文字以上を含むか 
    sub ja( ) { grep ( /^jp$|^ja(p(a(n?)?)?)?/i , @ARGV ) ? 1 : 0 } # jp または japan という文字列を先頭から2文字以上を含むか 
    sub opt( ) { grep (/^opt(i(o(ns?)?)?)?$/i, @ARGV ) ? 1 : 0 } # options という文字列を先頭から3文字以上含むから
    sub noPOD ( ) { grep (/^no-?pod\b/i, @ARGV) ? 1 : 0 } # POD を使わないと言う指定がされているかどうか
    my $jd = "JapaneseManual" ;
    my $flagE = ! ja && ( en || ! EnvJ ) ; # 英語にするかどうかのフラグ

    exec "perldoc $0" if $flagE &&  ! opt && ! noPOD   ; 
    $ARGV[1] //= '' ;
    open my $FH , '<' , $0 ;
    while(<$FH>){
        s/\Qboxmuller\E/$Script/gi ;
        s/\$Bin/$Bin/gi ;
        if ( s/^=head1\b\s*// .. s/^=cut\b\s*// ) { 
            if ( s/^=begin\s+$jd\b\s*// .. s/^=end\s+$jd\b\s*// xor $flagE ) {
                print $_ if ! opt || m/^\s+\-/  ; 
            }
        } 
    }
    close $FH ;
    exit 0 ;
}

=encoding utf8 

=head1 NAME

boxmuller

=head1 VERSION

0.24 -- 2018-07-03

=head1 SYNOPSIS

boxmuller [B<-m> mean] [B<-v> variance | B<-d> standard_deviation] 
[B<-g> how_many_you_want] [B<-.> digits_after_decimal_point] [B<-s> random_seed] 
[B<-L>(log normal)] [B<-@> seconds] [B<-1>] [B<-:>]

boxmuller [B<--help> [ja|en] [opt] [nopod]] [B<--version>]

=head1 DESCRIPTION

Generates Gaussian random variables by Box-Muller method.
The used random seed and the sums of the generated numbers and the square of them are also 
provided to STDERR.

=head1 OPTION

=over 4

=item -m N 

Population B<Mean (average)>. Default value is 0.


=item -d N 

Population B<Standard Deviation>. Default value is 1.

=item -v N 

Population B<Variance>. Default value is 1. If -d is given, -v would be nullified.

=item -. N 

The digits to be shown after the decimal point. Without this specification 
14 digits after the decimal point may be shown.

=item -g N 

How many numbers to be produced. Default value is 6. "Inf" (Infinity) can be given.

=item -s N 

Random B<seed> specification. The residual divided by 2**32 is essential.

=item -L 

Outputs variables from the B<log normal distribution> instead of the normal distribution.

=item -1

Only output the random number without other secondary information.

=item -: 

Attaches serial number beginning from 1. 

=item -@ N 

Waiting time in B<seconds> for each output line, that can be spedicifed 6 digits after the decimal points



( run in 3.072 seconds using v1.01-cache-2.11-cpan-39bf76dae61 )