/* -*- Mode:C++; c-basic-offset:8; tab-width:8; indent-tabs-mode:t -*- */
/*
* Copyright (c) Xerox Corporation 1997. All rights reserved.
*
* 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
*
* Linking this file statically or dynamically with other modules is making
* a combined work based on this file. Thus, the terms and conditions of
* the GNU General Public License cover the whole combination.
*
* In addition, as a special exception, the copyright holders of this file
* give you permission to combine this file with free software programs or
* libraries that are released under the GNU LGPL and with code included in
* the standard release of ns-2 under the Apache 2.0 license or under
* otherwise-compatible licenses with advertising requirements (or modified
* versions of such code, with unchanged license). You may copy and
* distribute such a system following the terms of the GNU GPL for this
* file and the licenses of the other code concerned, provided that you
* include the source code of that other code when and as the GNU GPL
* requires distribution of source code.
*
* Note that people who make modified versions of this file are not
* obligated to grant this special exception for their modified versions;
* it is their choice whether to do so. The GNU General Public License
* gives permission to release a modified version without this exception;
* this exception also makes it possible to release a modified version
* which carries forward this exception.
*/
#ifndef lint
static const char rcsid[] =
"@(#) $Header: /nfs/jade/vint/CVSROOT/ns-2/tools/ranvar.cc,v 1.19 2005/08/26 05:05:31 tomh Exp $ (Xerox)";
#endif
#include <stdio.h>
#include "ranvar.h"
RandomVariable::RandomVariable()
{
rng_ = RNG::defaultrng();
}
int RandomVariable::command(int argc, const char*const* argv)
{
Tcl& tcl = Tcl::instance();
if (argc == 2) {
if (strcmp(argv[1], "value") == 0) {
tcl.resultf("%6e", value());
return(TCL_OK);
}
}
if (argc == 3) {
if (strcmp(argv[1], "use-rng") == 0) {
rng_ = (RNG*)TclObject::lookup(argv[2]);
if (rng_ == 0) {
tcl.resultf("no such RNG %s", argv[2]);
return(TCL_ERROR);
}
return(TCL_OK);
}
}
return(TclObject::command(argc, argv));
}
// Added by Debojyoti Dutta 12 October 2000
// This allows us to seed a randomvariable with
// our own RNG object. This command is called from
// expoo.cc and pareto.cc
int RandomVariable::seed(char *x){
Tcl& tcl = Tcl::instance();
rng_ = (RNG*)TclObject::lookup(x);
if (rng_ == 0) {
tcl.resultf("no such RNG %s", x);
return(TCL_ERROR);
}
return(TCL_OK);
}
static class UniformRandomVariableClass : public TclClass {
public:
UniformRandomVariableClass() : TclClass("RandomVariable/Uniform"){}
TclObject* create(int, const char*const*) {
return(new UniformRandomVariable());
}
} class_uniformranvar;
UniformRandomVariable::UniformRandomVariable()
{
bind("min_", &min_);
bind("max_", &max_);
}
UniformRandomVariable::UniformRandomVariable(double min, double max)
{
min_ = min;
max_ = max;
}
double UniformRandomVariable::value()
{
return(rng_->uniform(min_, max_));
}
static class ExponentialRandomVariableClass : public TclClass {
public:
ExponentialRandomVariableClass() : TclClass("RandomVariable/Exponential") {}
TclObject* create(int, const char*const*) {
return(new ExponentialRandomVariable());
}
} class_exponentialranvar;
ExponentialRandomVariable::ExponentialRandomVariable()
{
bind("avg_", &avg_);
}
ExponentialRandomVariable::ExponentialRandomVariable(double avg)
{
avg_ = avg;
}
double ExponentialRandomVariable::value()
{
return(rng_->exponential(avg_));
}
static class ParetoRandomVariableClass : public TclClass {
public:
ParetoRandomVariableClass() : TclClass("RandomVariable/Pareto") {}
TclObject* create(int, const char*const*) {
return(new ParetoRandomVariable());
}
} class_paretoranvar;
ParetoRandomVariable::ParetoRandomVariable()
{
bind("avg_", &avg_);
bind("shape_", &shape_);
}
ParetoRandomVariable::ParetoRandomVariable(double avg, double shape)
{
avg_ = avg;
shape_ = shape;
}
double ParetoRandomVariable::value()
{
/* yuck, user wants to specify shape and avg, but the
* computation here is simpler if we know the 'scale'
* parameter. right thing is to probably do away with
* the use of 'bind' and provide an API such that we
* can update the scale everytime the user updates shape
* or avg.
*/
return(rng_->pareto(avg_ * (shape_ -1)/shape_, shape_));
}
/* Pareto distribution of the second kind, aka. Lomax distribution */
static class ParetoIIRandomVariableClass : public TclClass {
public:
ParetoIIRandomVariableClass() : TclClass("RandomVariable/ParetoII") {}
TclObject* create(int, const char*const*) {
return(new ParetoIIRandomVariable());
}
} class_paretoIIranvar;
ParetoIIRandomVariable::ParetoIIRandomVariable()
{
bind("avg_", &avg_);
bind("shape_", &shape_);
}
ParetoIIRandomVariable::ParetoIIRandomVariable(double avg, double shape)
{
avg_ = avg;
shape_ = shape;
}
double ParetoIIRandomVariable::value()
{
return(rng_->paretoII(avg_ * (shape_ - 1), shape_));
}
static class NormalRandomVariableClass : public TclClass {
public:
NormalRandomVariableClass() : TclClass("RandomVariable/Normal") {}
TclObject* create(int, const char*const*) {
return(new NormalRandomVariable());
}
} class_normalranvar;
NormalRandomVariable::NormalRandomVariable()
{
bind("avg_", &avg_);
bind("std_", &std_);
}
double NormalRandomVariable::value()
{
return(rng_->normal(avg_, std_));
}
static class LogNormalRandomVariableClass : public TclClass {
public:
LogNormalRandomVariableClass() : TclClass("RandomVariable/LogNormal") {}
TclObject* create(int, const char*const*) {
return(new LogNormalRandomVariable());
}
} class_lognormalranvar;
LogNormalRandomVariable::LogNormalRandomVariable()
{
bind("avg_", &avg_);
bind("std_", &std_);
}
double LogNormalRandomVariable::value()
{
return(rng_->lognormal(avg_, std_));
}
static class ConstantRandomVariableClass : public TclClass {
public:
ConstantRandomVariableClass() : TclClass("RandomVariable/Constant"){}
TclObject* create(int, const char*const*) {
return(new ConstantRandomVariable());
}
} class_constantranvar;
ConstantRandomVariable::ConstantRandomVariable()
{
bind("val_", &val_);
}
ConstantRandomVariable::ConstantRandomVariable(double val)
{
val_ = val;
}
double ConstantRandomVariable::value()
{
return(val_);
}
/* Hyperexponential distribution code adapted from code provided
* by Ion Stoica.
*/
static class HyperExponentialRandomVariableClass : public TclClass {
public:
HyperExponentialRandomVariableClass() :
TclClass("RandomVariable/HyperExponential") {}
TclObject* create(int, const char*const*) {
return(new HyperExponentialRandomVariable());
}
} class_hyperexponentialranvar;
HyperExponentialRandomVariable::HyperExponentialRandomVariable()
{
bind("avg_", &avg_);
bind("cov_", &cov_);
alpha_ = .95;
}
HyperExponentialRandomVariable::HyperExponentialRandomVariable(double avg, double cov)
{
alpha_ = .95;
avg_ = avg;
cov_ = cov;
}
double HyperExponentialRandomVariable::value()
{
double temp, res;
double u = Random::uniform();
temp = sqrt((cov_ * cov_ - 1.0)/(2.0 * alpha_ * (1.0 - alpha_)));
if (u < alpha_)
res = Random::exponential(avg_ - temp * (1.0 - alpha_) * avg_);
else
res = Random::exponential(avg_ + temp * (alpha_) * avg_);
return(res);
}
/*
// Empirical Random Variable:
// CDF input from file with the following column
// 1. Possible values in a distrubutions
// 2. Number of occurances for those values
// 3. The CDF for those value
// code provided by Giao Nguyen
*/
static class EmpiricalRandomVariableClass : public TclClass {
public:
EmpiricalRandomVariableClass() : TclClass("RandomVariable/Empirical"){}
TclObject* create(int, const char*const*) {
return(new EmpiricalRandomVariable());
}
} class_empiricalranvar;
EmpiricalRandomVariable::EmpiricalRandomVariable() : minCDF_(0), maxCDF_(1), maxEntry_(32), table_(0)
{
bind("minCDF_", &minCDF_);
bind("maxCDF_", &maxCDF_);
bind("interpolation_", &interpolation_);
bind("maxEntry_", &maxEntry_);
}
int EmpiricalRandomVariable::command(int argc, const char*const* argv)
{
Tcl& tcl = Tcl::instance();
if (argc == 3) {
if (strcmp(argv[1], "loadCDF") == 0) {
if (loadCDF(argv[2]) == 0) {
tcl.resultf("%s loadCDF %s: invalid file",
name(), argv[2]);
return (TCL_ERROR);
}
return (TCL_OK);
}
}
return RandomVariable::command(argc, argv);
}
int EmpiricalRandomVariable::loadCDF(const char* filename)
{
FILE* fp;
char line[256];
CDFentry* e;
fp = fopen(filename, "r");
if (fp == 0)
return 0;
if (table_ == 0)
table_ = new CDFentry[maxEntry_];
for (numEntry_=0; fgets(line, 256, fp); numEntry_++) {
if (numEntry_ >= maxEntry_) { // resize the CDF table
maxEntry_ *= 2;
e = new CDFentry[maxEntry_];
for (int i=numEntry_-1; i >= 0; i--)
e[i] = table_[i];
delete table_;
table_ = e;
}
e = &table_[numEntry_];
// Use * and l together raises a warning
sscanf(line, "%lf %*f %lf", &e->val_, &e->cdf_);
}
fclose(fp);
return numEntry_;
}
double EmpiricalRandomVariable::value()
{
if (numEntry_ <= 0)
return 0;
double u = rng_->uniform(minCDF_, maxCDF_);
int mid = lookup(u);
if (mid && interpolation_ && u < table_[mid].cdf_)
return interpolate(u, table_[mid-1].cdf_, table_[mid-1].val_,
table_[mid].cdf_, table_[mid].val_);
return table_[mid].val_;
}
double EmpiricalRandomVariable::interpolate(double x, double x1, double y1, double x2, double y2)
{
double value = y1 + (x - x1) * (y2 - y1) / (x2 - x1);
if (interpolation_ == INTER_INTEGRAL) // round up
return ceil(value);
return value;
}
int EmpiricalRandomVariable::lookup(double u)
{
// always return an index whose value is >= u
int lo, hi, mid;
if (u <= table_[0].cdf_)
return 0;
for (lo=1, hi=numEntry_-1; lo < hi; ) {
mid = (lo + hi) / 2;
if (u > table_[mid].cdf_)
lo = mid + 1;
else hi = mid;
}
return lo;
}
syntax highlighted by Code2HTML, v. 0.9.1