1 package MathUtil;
   2 #
   3 # $RCSfile: MathUtil.pm,v $
   4 # $Date: 2012/03/17 21:35:26 $
   5 # $Revision: 1.21 $
   6 #
   7 # Author: Manish Sud <msud@san.rr.com>
   8 #
   9 # Copyright (C) 2004-2012 Manish Sud. All rights reserved.
  10 #
  11 # This file is part of MayaChemTools.
  12 #
  13 # MayaChemTools is free software; you can redistribute it and/or modify it under
  14 # the terms of the GNU Lesser General Public License as published by the Free
  15 # Software Foundation; either version 3 of the License, or (at your option) any
  16 # later version.
  17 #
  18 # MayaChemTools is distributed in the hope that it will be useful, but without
  19 # any warranty; without even the implied warranty of merchantability of fitness
  20 # for a particular purpose.  See the GNU Lesser General Public License for more
  21 # details.
  22 #
  23 # You should have received a copy of the GNU Lesser General Public License
  24 # along with MayaChemTools; if not, see <http://www.gnu.org/licenses/> or
  25 # write to the Free Software Foundation Inc., 59 Temple Place, Suite 330,
  26 # Boston, MA, 02111-1307, USA.
  27 #
  28 
  29 use strict;
  30 use Exporter;
  31 use Constants;
  32 use Math::Trig ();
  33 use POSIX ();
  34 
  35 use vars qw(@ISA @EXPORT @EXPORT_OK %EXPORT_TAGS);
  36 
  37 @ISA = qw(Exporter);
  38 @EXPORT = qw(acos asin atan tan ceil floor log10 min max srandom random round);
  39 @EXPORT_OK = qw();
  40 
  41 %EXPORT_TAGS = (all  => [@EXPORT, @EXPORT_OK]
  42                );
  43 
  44 
  45 # Return next largest integer...
  46 sub ceil ($) {
  47   my($Value) = @_;
  48 
  49   return POSIX::ceil($Value);
  50 }
  51 
  52 # Return previous smallest integer...
  53 sub floor ($) {
  54   my($Value) = @_;
  55 
  56   return POSIX::floor($Value);
  57 }
  58 
  59 # Calculate log value using base 10...
  60 sub log10 ($) {
  61   my($Value) = @_;
  62 
  63   return CORE::log($Value)/CORE::log(10);
  64 }
  65 
  66 # Return the smaller of two numbers...
  67 sub min ($$) {
  68   my($Value1, $Value2) = @_;
  69 
  70   return ($Value1 <= $Value2) ? $Value1 : $Value2;
  71 }
  72 
  73 # Return the larger of two numbers...
  74 sub max ($$) {
  75   my($Value1, $Value2) = @_;
  76 
  77   return ($Value1 >= $Value2) ? $Value1 : $Value2;
  78 }
  79 
  80 # The random number generator implemented in MayaChemTools is a variant of linear
  81 # congruential generator (LCG) as described by Miller et al. [ Ref 120 ]. It is
  82 # also referred to as Lehmer random number generator or Park-Miller random number
  83 # generator.
  84 #
  85 # Unlike Perl's core random number generator function rand, the random number
  86 # generator implemented in MayaChemTools generates consistent random values
  87 # across different platforms - Windows, CygWin, Linux, Unix - for a specific random
  88 # seed.
  89 #
  90 
  91 # $RandomModulus = 2**31 - 1;
  92 # $RandomMultiplier = 16807;
  93 # $RandomQuotient = $RandomModulus / $RandomMultiplier;
  94 # $RandomRemainder = $RandomModulus % $RandomMultiplier
  95 #
  96 # $MaxRandomSeed = 2*31 -2
  97 #
  98 my($MaxRandomSeed, $RandomSeed, $RandomModulus, $RandomMultiplier, $RandomQuotient, $RandomRemainder);
  99 
 100 $MaxRandomSeed = 2147483646;
 101 $RandomSeed = 123456789;
 102 
 103 $RandomModulus = 2147483647;
 104 $RandomMultiplier = 16807;
 105 $RandomQuotient = 127773;
 106 $RandomRemainder = 2836;
 107 
 108 # Set random number seed...
 109 #
 110 # The intial value of random number seed is recommeded to be an integer between 1
 111 # and 2**31 - 2 [Ref 120] which translates to be 1 and 2147483646
 112 #
 113 sub srandom ($) {
 114   my($Seed) = @_;
 115 
 116   if ($Seed <= 0 ) {
 117     die "Error: srandom: Specified seed value must be greater than 0...";
 118   }
 119 
 120   $RandomSeed = ($Seed > $MaxRandomSeed) ? ($Seed % $MaxRandomSeed) : $Seed;
 121 
 122   return $RandomSeed;
 123 }
 124 
 125 # Retrun a random number between 0 and less than 1 or specified size...
 126 #
 127 sub random (;$) {
 128   my($Size) = @_;
 129   my($Value, $LowValue, $HighValue);
 130 
 131   $Size = defined $Size ? $Size : 1.0;
 132 
 133   $HighValue = $RandomSeed / $RandomQuotient;
 134   $LowValue = $RandomSeed % $RandomQuotient;
 135 
 136   $Value = $RandomMultiplier * $LowValue - $RandomRemainder * $HighValue;
 137 
 138   $RandomSeed = ($Value > 0) ? $Value : ($Value + $RandomModulus);
 139 
 140   return ($RandomSeed / $RandomModulus) * $Size;
 141 }
 142 
 143 # Round a integer/real number to:
 144 # . A nearest integer
 145 # . Specified number of decimal places
 146 #
 147 sub round ($;$) {
 148   my($Value, $DecimalPlaces) = @_;
 149   my($RoundedValue);
 150 
 151   if (defined($DecimalPlaces) && $DecimalPlaces > 0) {
 152     $RoundedValue = sprintf "%.${DecimalPlaces}f", $Value;
 153   }
 154   else {
 155     if ($Value < 0) {
 156       $RoundedValue = int($Value - 0.5);
 157     }
 158     else {
 159       $RoundedValue = int($Value + 0.5);
 160     }
 161   }
 162   return $RoundedValue;
 163 }
 164 
 165 # Return tangent of an angle expressed in radians.
 166 sub tan {
 167   my($Value) = @_;
 168 
 169   return (CORE::sin($Value)/CORE::cos($Value));
 170 }
 171 
 172 # Return inverse sine of an angle expressed in radians.
 173 #
 174 # For a right angle triangle defined by sides X and Y in a unit circle, Pythagorean theorem implies
 175 # X**2 + Y**2 = 1 and sin value corresponds to Y. So asin is equivalent to atan2(Y, sqrt(1-Y**2)).
 176 # However, taking sqrt of negative numbers is problematic; Math::Trig::asin handles it using complex
 177 # numbers.
 178 #
 179 sub asin ($) {
 180   my($Value) = @_;
 181 
 182   return Math::Trig::asin($Value);
 183 }
 184 
 185 # Return inverse cosine of an angle expressed in radians.
 186 #
 187 # For a right angle triangle defined by sides X and Y in a unit circle, Pythagorean theorem implies
 188 # X**2 + Y**2 = 1 and cos value corresponds to X. So asin is equivalent to atan2(sqrt(1-X**2), X)
 189 # However, taking sqrt of negative numbers is problematic; Math::Trig::acos handles it using complex
 190 # numbers.
 191 #
 192 sub acos ($) {
 193   my($Value) = @_;
 194 
 195   return Math::Trig::acos($Value);
 196 }
 197