...
 
......@@ -3,8 +3,9 @@ this is <CHANGELOG>
Recent development in Seitosh (bug fixes, new features, etc)
------------------------------------------------------------
01.11.2019: libtfxx
01.11.2019: libtfxx, randomseries
- commit 27491f4: let RNGgaussian support selection of rng type
- commit d29bb3d: provide new program synt/misc/randomseries.cc
28.06.2019: any2matlab
- commit 7efe5ac: fix variable declaration to make compiler happy
......
......@@ -147,6 +147,58 @@ namespace tfxx {
}
} // RNGgaussian::UTException::UTException
/* ====================================================================== */
const char* RNGgaussian::comment_gsl_rng_ranlux=
{
"The recommended random number generator type is \"ranlxd2\"\n"
"\n"
"Quotation from \n"
"https://www.gnu.org/software/gsl/doc/html/index.html\n"
"https://www.gnu.org/software/gsl/doc/html/rng.html#random-number-generator-algorithms\n"
"\n"
"The following generators are recommended for use in simulation. They have\n"
"extremely long periods, low correlation and pass most statistical tests.\n"
"For the most reliable source of uncorrelated numbers, the second-generation\n"
"RANLUX generators have the strongest proof of randomness.\n"
"\n"
"gsl_rng_ranlxd1\n"
"gsl_rng_ranlxd2\n"
"\n"
" These generators produce double precision output (48 bits) from the\n"
" RANLXS generator. The library provides two luxury levels ranlxd1 and\n"
" ranlxd2, in increasing order of strength.\n"
"\n"
"gsl_rng_ranlux\n"
"gsl_rng_ranlux389\n"
"\n"
" The ranlux generator is an implementation of the original algorithm\n"
" developed by Luscher. It uses a lagged-fibonacci-with-skipping\n"
" algorithm to produce “luxury random numbers”. It is a 24-bit\n"
" generator, originally designed for single-precision IEEE floating\n"
" point numbers. This implementation is based on integer arithmetic,\n"
" while the second-generation versions RANLXS and RANLXD described\n"
" above provide floating-point implementations which will be faster\n"
" on many platforms. The period of the generator is about 10^{171}.\n"
" The algorithm has mathematically proven properties and it can\n"
" provide truly decorrelated numbers at a known level of randomness.\n"
" The default level of decorrelation recommended by Luscher is\n"
" provided by gsl_rng_ranlux, while gsl_rng_ranlux389 gives the\n"
" highest level of randomness, with all 24 bits decorrelated. Both\n"
" types of generator use 24 words of state per generator.\n"
"\n"
" For more information see,\n"
"\n"
" * M. Luscher, “A portable high-quality random number generator for\n"
" lattice field theory calculations”, Computer Physics Communications,\n"
" 79 (1994) 100–110.\n"
" DOI: 10.1016/0010-4655(94)90232-1.\n"
" * F. James, “RANLUX: A Fortran implementation of the high-quality\n"
" pseudo-random number generator of Luscher”, Computer Physics\n"
" Communications, 79 (1994) 111–114\n"
" DOI: 10.1016/0010-4655(94)90233-X\n"
}; // RNGgaussian::comment_gsl_rng_ranlux
} // namespace numeric
} // namespace tfxx
......
......@@ -100,6 +100,8 @@ namespace tfxx {
double mean() const { return Mmean; }
//! print list of random number generators to stream.
static void rng_list_types(std::ostream& os);
//! comment on GSL RNGs
static const char* comment_gsl_rng_ranlux;
private:
//! store standard deviation and mean value
double Mstd, Mmean;
......
......@@ -59,11 +59,14 @@ int main(int iargc, char* argv[])
"usage: rngtest [-l] [-v] [-t type] [-seed v] [-mean v]\n"
" [-std v] [-n n] [run]\n"
" or: rngtest --help|-h" "\n"
" or: rngtest --xhelp" "\n"
};
// define full help text
char help_text[]=
{
"--xhelp print comments regarding GSL RNG implementation\n"
"\n"
"run run the program\n"
"-l print list of random number generators and exit\n"
"-v be verbose\n"
......@@ -103,6 +106,8 @@ int main(int iargc, char* argv[])
{"std",arg_yes,"1."},
// 7: number of samples to be printed
{"n",arg_yes,"10"},
// 8: print comment text
{"xhelp",arg_no,"-"},
{NULL}
};
......@@ -117,10 +122,15 @@ int main(int iargc, char* argv[])
Commandline cmdline(iargc, argv, options);
// help requested? print full help text...
if (cmdline.optset(0))
if (cmdline.optset(0) || cmdline.optset(8))
{
cerr << usage_text << endl;
cerr << help_text << endl;
if (cmdline.optset(8))
{
cerr << endl;
cerr << tfxx::numeric::RNGgaussian::comment_gsl_rng_ranlux << endl;
}
exit(0);
}
......
......@@ -35,7 +35,8 @@
PROGRAMS=dispfield kette planefield coupledoscillators \
dispfieldx lamb lambx modeinterference \
siggen tesiff teswf hamres siggenx phasedsignals
siggen tesiff teswf hamres siggenx phasedsignals \
randomseries
# gez does not compile correctly
.PHONY: all
......@@ -109,6 +110,10 @@ siggenx: %x: %.o
modeinterference: modeinterference.o
$(FC) -o $@ $< -ltf $(TF_LINK_PGPLOT) $(LDFLAGS) -L$(LOCLIBDIR)
randomseries: randomseries.o
$(CXX) -o $@ $< -ltfxx -ldatrwxx -lsffxx -lgsexx -ltime++ -laff \
-lgsl -lgslcblas $(LDFLAGS) -L$(LOCLIBDIR)
.PHONY: list-libraries
list-libraries:
grep ' -l' Makefile | tr ' ' '\n' | egrep '^-l' |sort | uniq
......
/*! \file randomseries.cc
* \brief provide time series of incoherent random numbers
*
* ----------------------------------------------------------------------------
*
* \author Thomas Forbriger
* \date 31/10/2019
*
* provide time series of incoherent random numbers
*
* Copyright (c) 2019 by Thomas Forbriger (BFO Schiltach)
*
* REVISIONS and CHANGES
* - 31/10/2019 V1.0 Thomas Forbriger
*
* ============================================================================
*/
#define RANDOMSERIES_VERSION \
"RANDOMSERIES V1.0 provide time series of incoherent random numbers"
#include <iostream>
#include <fstream>
#include <string>
#include <sstream>
#include <tfxx/commandline.h>
#include <tfxx/seitosh.h>
#include <tfxx/error.h>
#include <tfxx/rng.h>
#include <aff/series.h>
#include <aff/iterator.h>
#include <datrwxx/writeany.h>
using std::cout;
using std::cerr;
using std::endl;
/* ---------------------------------------------------------------------- */
struct Options {
bool verbose, overwrite;
std::string rngtype, filename, filetype;
int nseries, nsamples;
double std, mean, dt;
long unsigned int seed;
}; // struct Options
/* ---------------------------------------------------------------------- */
int main(int iargc, char* argv[])
{
// define usage information
char usage_text[]=
{
RANDOMSERIES_VERSION "\n"
"usage: randomseries filename [-v] [-t type] [-o]" "\n"
" [-nseries n] [-nsamples n] [-dt v]" "\n"
" [-std v] [-mean v] [-seed n]" "\n"
" or: randomseries --help|-h" "\n"
" or: randomseries --xhelp[=type]\n"
};
// define full help text
char help_text[]=
{
"-help print description\n"
"-xhelp[=type] print detailed information\n"
" if 'type' is not specfied, or equals 'all', print\n"
" information regarding all supported file formats\n"
" if 'type' is specified, just print text for file\n"
" format 'type'\n"
" if 'type' is 'gslrng' print recommendations regarding\n"
" random number generator type\n"
"\n"
"filename name of output file\n"
"-v be verbose\n"
"-o overwrite existing output file\n"
"-t type select output file type\n"
"-nseries n produce 'n' time series\n"
"-nsamples n produce 'n' samples per series\n"
"-dt v set sampling interval to 'v' seconds\n"
"-seed n initialize random number generator with 'n'\n"
"-std v set standard deviation to 'v'\n"
"-mean v set mean value of samples to 'v'\n"
"\n"
"The program uses the GSL (GNU Scientific Library) random number\n"
"generators. See https://www.gnu.org/software/gsl/doc/html/rng.html\n"
"\n"
"The programs main purpose is to provide several time series\n"
"which are incoherent (as far as possible). All samples are taken\n"
"sequentially from the random number generator with the generator\n"
"being seeded only at the very beginning. If several (incoherent)\n"
"time series should be computed by individual calls to a random\n"
"number generator program, this would require careful seeding\n"
"for each program invocation to make sure that different sequences\n"
"of samples are indeed independent or at least incoherent.\n"
};
// define commandline options
using namespace tfxx::cmdline;
static Declare options[]=
{
// 0: print help
{"help",arg_no,"-"},
// 1: verbose mode
{"v",arg_no,"-"},
// 2: overwrite existing output file
{"o",arg_no,"-"},
// 3: select output file type
{"t",arg_yes,"sff"},
// 4: select output file type
{"nseries",arg_yes,"1"},
// 5: select output file type
{"nsamples",arg_yes,"10000"},
// 6: select output file type
{"seed",arg_yes,"0"},
// 7: select output file type
{"std",arg_yes,"1."},
// 8: select output file type
{"mean",arg_yes,"0."},
// 9: output file format
{"xhelp",arg_opt,"all"},
// 10: output file format
{"rngtype",arg_yes,"default"},
// 11: sampling interval
{"dt",arg_yes,"1."},
{NULL}
};
// no arguments? print usage...
if (iargc<2)
{
cerr << usage_text << endl;
cerr << endl << tfxx::seitosh::repository_reference << endl;
exit(0);
}
// collect options from commandline
Commandline cmdline(iargc, argv, options);
// help requested? print full help text...
if (cmdline.optset(0))
{
cerr << usage_text << endl;
cerr << help_text << endl;
cerr << endl;
cerr << "supported output data file types:" << endl;
datrw::supported_output_data_types(cerr);
cerr << endl << tfxx::seitosh::repository_reference << endl;
exit(0);
}
// help on file format details requested?
if (cmdline.optset(9))
{
cerr << usage_text << endl;
cerr << endl;
if (cmdline.string_arg(9) == "all")
{
datrw::online_help(cerr);
}
else if (cmdline.string_arg(9) == "gslrng")
{
cerr << tfxx::numeric::RNGgaussian::comment_gsl_rng_ranlux << endl;
cerr << endl;
cerr << "Available random number generator types:" << endl;
tfxx::numeric::RNGgaussian::rng_list_types(cerr);
}
else
{
datrw::online_help(cmdline.string_arg(9), cerr);
}
exit(0);
}
TFXX_assert(cmdline.extra(), "Missing file name");
// extract command line options
Options opt;
opt.filename=cmdline.next();
opt.verbose=cmdline.optset(1);
opt.overwrite=cmdline.optset(2);
opt.filetype=cmdline.string_arg(3);
opt.nseries=cmdline.int_arg(4);
opt.nsamples=cmdline.int_arg(5);
opt.seed=cmdline.int_arg(6);
opt.std=cmdline.double_arg(7);
opt.mean=cmdline.double_arg(8);
opt.rngtype=cmdline.string_arg(10);
opt.dt=cmdline.double_arg(11);
// check parameters
TFXX_assert(opt.nseries>0, "number of series must be positive");
TFXX_assert(opt.nsamples>0, "number of samples must be positive");
TFXX_assert(opt.seed>=0, "number of samples must be non-negative");
TFXX_assert(opt.dt>0, "sampling interval must be positive");
if (opt.verbose)
{
cout << RANDOMSERIES_VERSION << endl;
}
// open output file
if (opt.verbose)
{
cout << "open output file " << opt.filename
<< " with format " << opt.filetype << endl;
}
if (!opt.overwrite) { datrw::abort_if_exists(opt.filename); }
std::ofstream ofs(opt.filename.c_str(),
datrw::oanystream::openmode(opt.filetype));
datrw::oanystream os(ofs, opt.filetype);
// create random number generator
tfxx::numeric::RNGgaussian rng(opt.std, opt.mean, opt.rngtype.c_str());
if (opt.verbose)
{
cout << "RNG parameters:\n"
<< " type: " << rng.type() << "\n"
<< " seed: " << rng.seed() << "\n"
<< " mean: " << rng.mean() << "\n"
<< " std: " << rng.std() << endl;
cout << "produce " << opt.nseries << " traces with "
<< opt.nsamples << " samples each" << endl;
cout << "set sampling interval to " << opt.dt << " seconds" << endl;
}
sff::WID2 wid2;
wid2.nsamples=opt.nsamples;
wid2.instype=opt.rngtype;
wid2.station="RNG";
wid2.dt=opt.dt;
aff::Series<double> series(opt.nsamples);
for (int iseries=0; iseries<opt.nseries; iseries++)
{
if (opt.verbose)
{
cout << "produce series #" << iseries+1 << endl;
}
aff::Iterator<aff::Series<double> > CS(series);
while (CS.valid())
{
*CS = rng.value();
++CS;
}
std::ostringstream oss;
oss << iseries+1;
wid2.channel=oss.str();
os << wid2;
os << series;
}
}
/* ----- END OF randomseries.cc ----- */