From 5df79c53745fde5d6c3340a2979b1429cd5892c1 Mon Sep 17 00:00:00 2001 From: Henrik Rydberg Date: Sat, 8 Oct 2011 20:30:28 +0200 Subject: Initial import of htcd system 1.0 Signed-off-by: Henrik Rydberg --- src/mt/random.cc | 95 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 95 insertions(+) create mode 100644 src/mt/random.cc (limited to 'src/mt/random.cc') diff --git a/src/mt/random.cc b/src/mt/random.cc new file mode 100644 index 0000000..87a4c8c --- /dev/null +++ b/src/mt/random.cc @@ -0,0 +1,95 @@ +/************************************************************************* + * + * HTCd - Copyright (C) 1998-2006 Henrik Rydberg + * + * This program is free software; you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation; either version 2 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program; if not, write to the Free Software + * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + */ + +#include +#include +#include + +////////////////////////////////////////////////////// +// +// This file implements some routines from NumRec +// + +const int NTAB=32; // NOTE: change in header if changing this + +const int IA=16807; +const int IM=2147483647; +const float AM=1/float(IM); +const int IQ=127773; +const int IR=2836; +const int NDIV=1+(IM-1)/NTAB; +const float EPS=1.2e-7; +const float RNMX=1-EPS; + +static int idum,iy,iv[32],flag; +static float extra; +static int inited; +static mutex_t mutex; + +unsigned long RandomInit(unsigned long seed) +{ + inited=1; + idum=seed?seed:1; + for(int j=NTAB+7;j>=0;j--) { + int k=idum/IQ; + idum=IA*(idum-k*IQ)-IR*k; + if(idum<0) idum+=IM; + if(jRNMX) temp=RNMX; + return temp; +} + +float Exponential() +{ + return -log(Uniform()); +} + +float Gaussian() +{ + MLOCK(mutex); + float v1,v2,rsq; + if(flag) { + flag=0; + return extra; + } else { + do { + v1=2*Uniform()-1; + v2=2*Uniform()-1; + rsq=v1*v1+v2*v2; + } while(rsq>=1||rsq<=0); + float fac=sqrt(-2*log(rsq)/rsq); + flag=1; extra=v1*fac; + return v2*fac; + } +} -- cgit v1.2.3