% File: wiener.w [CWEB source code] % Residence: http://jonsson.eu/programs/cweb/wiener/ % Created: November 10, 2011 [v.1.0] % Last change: November 10, 2011 [v.1.0] % Author: Fredrik Jonsson % Description: The CWEB source code for the WIENER simulator, generating a % sequence of data points corresponding to the Wiener process. % The program employs Donald Knuth's proposed random number % generator, as listed in and the Box-Muller transform. For % information on the CWEB programming language, see % http://www.literateprogramming.com. % Compilation: Compile this program by using the enclosed Makefile, or use % the blocks of the Makefile as listed in the documentation % (file wiener.ps or wiener.pdf). The C source code (as % generated from this CWEB code) conforms to the ANSI standard % for the C programming language (which is equivalent to the % ISO C90 standard for C). % % Copyright (C) 2011, Fredrik Jonsson % Non-commercial copying welcome. % \input epsf %% \input amstex \def\version{1.0} \def\lastrevdate{November 11, 2011} \font\eightcmr=cmr8 \font\tensc=cmcsc10 \font\eightcmssq=cmssq8 \font\eightcmssqi=cmssqi8 \font\twentycmcsc=cmcsc10 at 20 truept \def\wiener{{\eightcmr WIENER\spacefactor1000}} \def\poincare{{\eightcmr POINCARE\spacefactor1000}} \def\ANSI{{\eightcmr ANSI\spacefactor1000}} % The language standard we stick to \def\GNU{{\eightcmr GNU\spacefactor1000}} % GNU is Not Unix \def\GCC{{\eightcmr GCC\spacefactor1000}} % The GNU C-compiler \def\CEE{{\eightcmr C\spacefactor1000}} % The C programming language \def\ANSICEE{{\eightcmr ANSI~C\spacefactor1000}}% The ANSI-C standard language \def\CWEB{{\eightcmr CWEB\spacefactor1000}} % The CWEB programming language \def\UNIX{{\eightcmr UNIX\spacefactor1000}} \def\CYGWIN{{\eightcmr CYGWIN\spacefactor1000}} \def\GNUPLOT{{\eightcmr GNUPLOT\spacefactor1000}} \def\CTANGLE{{\eightcmr CTANGLE\spacefactor1000}} \def\CWEAVE{{\eightcmr CWEAVE\spacefactor1000}} \def\mod{\mathop{\rm mod}\nolimits} % The real part of a complex number \def\re{\mathop{\rm Re}\nolimits} % The real part of a complex number \def\im{\mathop{\rm Im}\nolimits} % The imaginary part of a complex number %% \def\bbr{\mathbb{R}} % The blackboard bold 'R' as in 'R^3' \def\backslash{\char'134} % The `\' character \def\vertbar{\char'174} % The `|' character \def\dollar{\char'044} % The `$' character \def\tilde{\char'176} % The `~' character \def\tothepower{\char'136} % The `^' character \def\onehalf{{\textstyle{{1}\over{2}}}} \def\threefourth{{\textstyle{{3}\over{4}}}} \def\endalg{\vrule height 1.4ex width .6ex depth .4ex} % Rectangular bullet \def\ie{i.\thinspace{e.}~\ignorespaces} \def\eg{e.\thinspace{g.}~\ignorespaces} \let\,\thinspace %% \def\subsection[#1]{\goodbreak\noindent{\it #1}\par\nobreak\noindent} % % Define a handy macro for listing the steps of an algorithm. % \newdimen\aitemindent \aitemindent=26pt \newdimen\aitemleftskip \aitemleftskip=36pt \def\aitem[#1]{\smallbreak\noindent\hbox to 10pt{}% \hbox to\aitemindent{\bf #1\hfill}% \hangindent\aitemleftskip\ignorespaces} % % Define a handy macro for bibliographic references. % \newdimen\refitemindent \refitemindent=18pt \def\refitem[#1]{\smallbreak\noindent% \hbox to\refitemindent{[#1]\hfill}% \hangindent\refitemindent\ignorespaces} \def\references{\medskip\goodbreak\noindent{\bf References}\nobreak\smallskip} % % Define a handy macro for nicely typeset variable descriptions. % \newdimen\varitemindent \varitemindent=100pt \def\varitem[#1]#2{\smallbreak\noindent\hbox to 20pt{}% \hbox to\varitemindent{#1\hfill}% \hangindent 120pt\ignorespaces#2\par} % % Define a handy macro for nicely typeset descriptions of command line options. % \newdimen\optitemindent \optitemindent=80pt \def\optitem[#1]#2{\smallbreak\noindent\hbox to\hsize{\hskip40pt{#1}\hfill} \hangindent\optitemindent\ignorespaces{\hskip\optitemindent{#2}}\smallskip} % % Define a handy macro for the list of program revisions. % \newdimen\citemindent \citemindent=80pt \newdimen\citemleftskip \citemleftskip=90pt \def\citem[#1]{\smallbreak\noindent\hbox to 10pt{}% \hbox to\citemindent{\bf #1\hfill}% \hangindent\citemleftskip\ignorespaces} \def\citindent{\smallbreak\noindent\hbox to 10pt{}% \hbox to\citemindent{\hfil}% \hangindent\citemleftskip\ignorespaces} % % Define a handy macro for the typesetting of figure captions. The double % '{{' and '}}' are here in place to prevent the increased \leftskip and % \rightskip to apply globally after the macro has been expanded. % \newdimen\figcapindent \figcapindent=40pt \def\figcaption[#1]#2{{\advance\leftskip by\figcapindent \advance\rightskip by\figcapindent \smallskip\noindent{\bf Figure #1.} {#2}\smallskip}} % % Define the \beginvrulealign and \endvrulealign macros as described in % Donald Knuth's The TeXbook, Appendix D: Dirty Tricks. These macros are % used in typesetting nicely looking tables. % \def\beginvrulealign{\setbox0=\vbox\bgroup} \def\endvrulealign{\egroup % now \box0 holds the entire alignment \setbox0=\vbox{\setbox2=\hbox{\vrule height\ht0 depth\dp0 width0pt} \unvbox0 \setbox0=\lastbox % now \box0 is the bottom row \nointerlineskip \copy0 % put it back \global\setbox1=\hbox{} % initialize box that will contain rules \setbox4=\hbox{\unhbox0 % now open up the bottom row \loop \skip0=\lastskip \unskip % remove tabskip glue \advance\skip0 by-.4pt % rules are .4 pt wide \divide\skip0 by 2 \global\setbox1=\hbox{\hskip\skip0\vrule\hskip\skip0 \unhbox2\unhbox1}% \setbox2=\lastbox % remove alignment entry \ifhbox2 \setbox2=\hbox{\kern\wd2}\repeat}}% \hbox{\rlap{\box0}\box1}} % superimpose the alignment on the rules % % Make sure that the METAPOST logo can be loaded even in plain TeX. % \ifx\MP\undefined \ifx\manfnt\undefined \font\manfnt=logo10 \fi \ifx\manfntsl\undefined \font\manfntsl=logosl10 \fi \def\MP{{\ifdim\fontdimen1\font>0pt \let\manfnt = \manfntsl \fi {\manfnt META}\-{\manfnt POST}}\spacefactor1000 }% \fi \ifx\METAPOST\undefined \let\METAPOST=\MP \fi \datethis @*Introduction. \vskip 60pt \centerline{\twentycmcsc Wiener} \vskip 20pt \centerline{A computer program for the generation of numerical data simulating a Wiener process.} \vskip 2pt \centerline{(Version \version\ of \lastrevdate)} \vskip 10pt \centerline{Written by Fredrik Jonsson} \vskip 10pt \centerline{\epsfxsize=88mm\epsfbox{fig3.1}} \vskip 10pt \noindent This \CWEB\footnote{${}^\dagger$}{For information on the \CWEB\ programming language by Donald E.~Knuth, as well as samples of \CWEB\ programs, see {\tt http://www-cs-faculty.stanford.edu/\~\ \kern -5pt knuth/cweb.html}. For general information on literate programming, see {\tt http://www.literateprogramming.com}.} computer program computes a series of floating-point numbers corresponding to a Wiener process in $D$ dimensions. It relies on the random number generator as proposed by Donald Knuth in {\sl The Art of Computer Programming, Volume 1 -- Fundamental Algorithms}, 3rd edition (Addison-Wesley, Boston, 1998), generating numbers which are fed into the Box--Muller transform to generate the normal distribution associated with the Wiener process. Besides providing a simulator of the Wiener process, the \wiener\ program can also be used in a ``lock-in'' mode with zero expectation value for each data point, providing a pretty good random number generator for large series of stochastic data, not relying on the (rather poor) generators available in standard \CEE\ libraries. The \wiener\ program does not solve any problem {\sl per se}, but is merely to be considered as a generator of statistical data to be used by other applications in modeling of physical, chemical or financial processes. \bigskip \noindent Copyright \copyright\ Fredrik Jonsson, 2011. All rights reserved. Non-commercial copying welcome. \vfill @*The Wiener process. ``In mathematics, the Wiener process is a continuous-time stochastic process named in honor of Norbert Wiener. It is often called standard Brownian motion, after Robert Brown. It is one of the best known L\'evy processes (c\`adl\`ag stochastic processes with stationary independent increments) and occurs frequently in pure and applied mathematics, economics and physics. The Wiener process plays an important role both in pure and applied mathematics. In pure mathematics, the Wiener process gave rise to the study of continuous time martingales. It is a key process in terms of which more complicated stochastic processes can be described. As such, it plays a vital role in stochastic calculus, diffusion processes and even potential theory. It is the driving process of Schramm--Loewner evolution. In applied mathematics, the Wiener process is used to represent the integral of a Gaussian white noise process, and so is useful as a model of noise in electronics engineering, instrument errors in filtering theory and unknown forces in control theory. The Wiener process has applications throughout the mathematical sciences. In physics it is used to study Brownian motion, the diffusion of minute particles suspended in fluid, and other types of diffusion via the Fokker--Planck and Langevin equations. It also forms the basis for the rigorous path integral formulation of quantum mechanics (by the Feynman--Kac formula, a solution to the Schr\"odinger equation can be represented in terms of the Wiener process) and the study of eternal inflation in physical cosmology. It is also prominent in the mathematical theory of finance, in particular the Black--Scholes option pricing model.'' \smallskip {\eightcmssq --\,Wikipedia, ``Wiener process'' (2011)} @ What the \wiener\ program does (and doesn't). The present \CWEB\ program does not solve any problems related to any of the processes described by models involving the Wiener process, but is merely an attempt to produce an as-good-as-possible result when simulating the Wiener process as such. In the \wiener\ program, special attention has been paid to the generation of random numbers, as this is a crucial and rather tricky problem when it comes to generating large non-recurring series of data. In the present program, the random number generator proposed by Donald Knuth\footnote{$\dagger$}{Donald E.~Knuth, {\sl The Art of Computer Programming, Volume 1 -- Fundamental Algorithms}, 3rd edition (Addison-Wesley, Boston, 1998).} has been employed, generating uniformly distributed numbers which are fed into the Box--Muller transform to generate the normal distribution associated with the Wiener process. Apart from being a pretty good and reliable generator of statistical data, to be used by other applications in modeling of physical, chemical or financial processes, the \wiener\ program does not solve any problems {\sl per se}. @ The CWEB programming language. For the reader who might not be familiar with the concept of the \CWEB\ programming language, the following citations hopefully will be useful. For further information, as well as freeware compilers for compiling \CWEB\ source code, see {\tt http://www.literateprogramming.com}. \bigskip {\narrower\narrower\narrower\narrower\eightcmssqi\noindent I believe that the time is ripe for significantly better documentation of programs, and that we can best achieve this by considering programs to be works of literature. Hence, my title: `Literate Programming.' Let us change our traditional attitude to the construction of programs: Instead of imagining that our main task is to instruct a computer what to do, let us concentrate rather on explaining to human beings what we want a computer to do. The practitioner of literate programming can be regarded as an essayist, whose main concern is with exposition and excellence of style. Such an author, with thesaurus in hand, chooses the names of variables carefully and explains what each variable means. He or she strives for a program that is comprehensible because its concepts have been introduced in an order that is best for human understanding, using a mixture of formal and informal methods that reinforce each other. \smallskip {\eightcmssq --\,Donald Knuth,} {\eightcmssqi The CWEB System of Structured Documentation} {\eightcmssq (Addison-Wesley, Massachusetts, 1994)} } \bigskip {\narrower\narrower\narrower\narrower\eightcmssqi\noindent The philosophy behind CWEB is that an experienced system programmer, who wants to provide the best possible documentation of his or her software products, needs two things simultaneously: a language like \TeX\ for formatting, and a language like C for programming. Neither type of language can provide the best documentation by itself; but when both are appropriately combined, we obtain a system that is much more useful than either language separately. The structure of a software program may be thought of as a `WEB' that is made up of many interconnected pieces. To document such a program we want to explain each individual part of the web and how it relates to its neighbors. The typographic tools provided by \TeX\ give us an opportunity to explain the local structure of each part by making that structure visible, and the programming tools provided by languages like C make it possible for us to specify the algorithms formally and unambiguously. By combining the two, we can develop a style of programming that maximizes our ability to perceive the structure of a complex piece of software, and at the same time the documented programs can be mechanically translated into a working software system that matches the documentation. Besides providing a documentation tool, CWEB enhances the C language by providing the ability to permute pieces of the program text, so that a large system can be understood entirely in terms of small sections and their local interrelationships. The CTANGLE program is so named because it takes a given web and moves the sections from their web structure into the order required by C; the advantage of programming in CWEB is that the algorithms can be expressed in ``untangled'' form, with each section explained separately. The CWEAVE program is so named because it takes a given web and intertwines the \TeX\ and C portions contained in each section, then it knits the whole fabric into a structured document. \smallskip {\eightcmssq --\,Donald Knuth, ``Literate Programming'', in} {\eightcmssqi Literate Programming} {\eightcmssq (CSLI Lecture Notes, Stanford, 1992)} } @ Revision history of the program. \medskip \citem[2011-11-11]{[v.1.0]} {\tt }\hfill\break First properly working version of the \wiener\ program, written in \CWEB\ and (\ANSI-conformant) \CEE. @*Compiling the source code. The program is written in \CWEB, generating \ANSICEE\ (ISO C90) conforming source code and documentation as plain \TeX-source, and is to be compiled using the sequences as outlined in the {\tt Makefile} listed below. \bigskip {\obeylines\obeyspaces\tt \# \# Makefile designed for use with ctangle, cweave, gcc, and plain TeX. \# \# Copyright (C) 2002-2011, Fredrik Jonsson \# \# The CTANGLE program converts a CWEB source document into a C program \# which may be compiled in the usual way. The output file includes \#line \# specifications so that debugging can be done in terms of the CWEB source \# file. \# \# The CWEAVE program converts the same CWEB file into a TeX file that may \# be formatted and printed in the usual way. It takes appropriate care of \# typographic details like page layout and the use of indentation, italics, \# boldface, etc., and it supplies extensive cross-index information that it \# gathers automatically. \# \# CWEB allows you to prepare a single document containing all the informa- \# tion that is needed both to produce a compilable C program and to produce \# a well-formatted document describing the program in as much detail as the \# writer may desire. The user of CWEB ought to be familiar with TeX as well \# as C. \# PROJECT = wiener FIGURES = fig1.eps fig2.eps fig3.eps ~ ~ CTANGLE = ctangle CWEAVE = cweave CC = gcc CCOPTS = -O2 -Wall -ansi -std=iso9899:1990 -pedantic LNOPTS = -lm TEX = tex DVIPS = dvips DVIPSOPT = -ta4 -D1200 PS2PDF = ps2pdf METAPOST = mpost ~ ~ all: \$(PROJECT) \$(FIGURES) \$(PROJECT).pdf ~ ~ \$(PROJECT): \$(PROJECT).o ~ \$(CC) \$(CCOPTS) -o \$(PROJECT) \$(PROJECT).o \$(LNOPTS) ~ ~ \$(PROJECT).o: \$(PROJECT).c ~ \$(CC) \$(CCOPTS) -c \$(PROJECT).c ~ ~ \$(PROJECT).c: \$(PROJECT).w ~ \$(CTANGLE) \$(PROJECT) ~ ~ \$(PROJECT).pdf: \$(PROJECT).ps ~ \$(PS2PDF) \$(PROJECT).ps \$(PROJECT).pdf ~ ~ \$(PROJECT).ps: \$(PROJECT).dvi ~ \$(DVIPS) \$(DVIPSOPT) \$(PROJECT).dvi -o \$(PROJECT).ps ~ ~ \$(PROJECT).dvi: \$(PROJECT).tex ~ \$(TEX) \$(PROJECT).tex ~ ~ \$(PROJECT).tex: \$(PROJECT).w ~ \$(CWEAVE) \$(PROJECT) ~ ~ \# \# Generate the Encapsulated PostScript image fig1.eps for the documentation. \# This is a 2D scatter plot of the uniformly distributed pseudo-random numbers \# prior to having been fed into the Box-Muller transform. \# fig1.eps: Makefile \$(PROJECT).w ~ wiener --uniform -D 2 -M 10000 > fig1.dat; ~ @echo "input graph;{\backslash} ~ ~ def mpdot = btex{\backslash} ~ ~ ~ {{\backslash}vrule height 0.5pt width 1.0pt depth 0.5pt}{\backslash} ~ ~ etex enddef;{\backslash} ~ ~ beginfig(1);{\backslash} ~ ~ draw begingraph(86mm,86mm);{\backslash} ~ ~ setrange(0,0,1,1);{\backslash} ~ ~ pickup pencircle scaled .5pt;{\backslash} ~ ~ gdraw {\backslash}"fig1.dat{\backslash}" plot mpdot;{\backslash} ~ ~ pickup pencircle scaled .25pt;{\backslash} ~ ~ autogrid(itick bot,itick lft);{\backslash} ~ ~ glabel.bot(btex \$\$ x\$\$ etex,OUT);{\backslash} ~ ~ glabel.lft(btex \$\$ y\$\$ etex,OUT);{\backslash} ~ ~ endgraph; endfig; end" > fig1.mp ~ ~ \$(METAPOST) fig1.mp ~ \$(TEX) -jobname=fig1 "{\backslash}input epsf{\backslash}nopagenumbers{\backslash} ~ ~ {\backslash}centerline{{\backslash}epsfbox{fig1.1}}{\backslash}bye" ~ \$(DVIPS) -D1200 -E fig1.dvi -o fig1.eps ~ ~ \# \# Generate the Encapsulated PostScript image fig2.eps for the documentation. \# This is a 2D scatter plot of the normally distributed pseudo-random numbers \# resulting from the Box-Muller transform. \# fig2.eps: Makefile \$(PROJECT).w ~ wiener --normal -D 2 -M 10000 > fig2.dat; ~ @echo "input graph;{\backslash} ~ ~ def mpdot = btex{\backslash} ~ ~ ~ {{\backslash}vrule height 0.5pt width 1.0pt depth 0.5pt}{\backslash} ~ ~ etex enddef;{\backslash} ~ ~ beginfig(1);{\backslash} ~ ~ draw begingraph(86mm,86mm);{\backslash} ~ ~ setrange(whatever,whatever,whatever,whatever);{\backslash} ~ ~ pickup pencircle scaled .5pt;{\backslash} ~ ~ gdraw {\backslash}"fig2.dat{\backslash}" plot mpdot;{\backslash} ~ ~ pickup pencircle scaled .25pt;{\backslash} ~ ~ autogrid(itick bot,itick lft);{\backslash} ~ ~ glabel.bot(btex \$\$ x\$\$ etex,OUT);{\backslash} ~ ~ glabel.lft(btex \$\$ y\$\$ etex,OUT);{\backslash} ~ ~ endgraph; endfig; end" > fig2.mp ~ ~ \$(METAPOST) fig2.mp ~ \$(TEX) -jobname=fig2 "{\backslash}input epsf{\backslash}nopagenumbers{\backslash} ~ ~ {\backslash}centerline{{\backslash}epsfbox{fig2.1}}{\backslash}bye" ~ \$(DVIPS) -D1200 -E fig2.dvi -o fig2.eps ~ ~ \# \# Generate the Encapsulated PostScript image fig3.eps for the documentation. \# This is a 2D graph showing the resulting simulated Wiener process. \# fig3.eps: Makefile \$(PROJECT).w ~ wiener -D 2 -M 10000 > fig3.dat; ~ @echo "input graph;{\backslash} ~ ~ beginfig(1);{\backslash} ~ ~ draw begingraph(86mm,86mm);{\backslash} ~ ~ setrange(whatever,whatever,whatever,whatever);{\backslash} ~ ~ pickup pencircle scaled .5pt;{\backslash} ~ ~ gdraw {\backslash}"fig3.dat{\backslash}";{\backslash} ~ ~ pickup pencircle scaled .25pt;{\backslash} ~ ~ autogrid(itick bot,itick lft);{\backslash} ~ ~ glabel.bot(btex \$\$ x\$\$ etex,OUT);{\backslash} ~ ~ glabel.lft(btex \$\$ y\$\$ etex,OUT);{\backslash} ~ ~ endgraph; endfig; end" > fig3.mp ~ ~ \$(METAPOST) fig3.mp ~ \$(TEX) -jobname=fig3 "{\backslash}input epsf{\backslash}nopagenumbers{\backslash} ~ ~ {\backslash}centerline{{\backslash}epsfbox{fig3.1}}{\backslash}bye" ~ \$(DVIPS) -D1200 -E fig3.dvi -o fig3.eps ~ ~ clean: ~ -rm -Rf \$(PROJECT) *~ *.c *.o *.exe *.dat *.pdf *.mp *.trj *.mpx ~ -rm -Rf *.tex *.aux *.log *.toc *.idx *.scn *.dvi *.ps *.1 *.eps ~ ~ archive: ~ make -ik clean ~ tar --gzip --directory=../ -cf ../\$(PROJECT).tgz \$(PROJECT) } \bigskip \noindent This {\tt Makefile} essentially executes two major calls. First, the \CTANGLE\ program parses the \CWEB\ source document {\tt wiener.w} to extract a \CEE\ source file {\tt wiener.c} which may be compiled in the usual way using any \ANSICEE\ conformant compiler. The output source file includes {\tt \#line} specifications so that any debugging can be done conveniently in terms of the original \CWEB\ source file. Second, the \CWEAVE\ program parses the same \CWEB\ source file to extract a plain \TeX\ source file {\tt wiener.tex} which may be compiled in the usual way. It takes appropriate care of typographic details like page layout and the use of indentation, italics, boldface, and so on, and it supplies extensive cross-index information that it gathers automatically. After having executed {\tt make} in the same catalogue where the files {\tt wiener.w} and {\tt Makefile} are located, one is left with an executable file {\tt wiener}, being the ready-to-use compiled program, and a PostScript file {\tt wiener.ps} which contains the full documentation of the program, that is to say the document you currently are reading. On platforms running any operating system by Microsoft, the executable file will instead automatically be called {\tt wiener.exe}. This convention also applies to programs compiled under the \UNIX-like environment \CYGWIN. @*Running the program. The program is entirely controlled by the command line options supplied when invoking the program. The syntax for executing the program is\par \medskip \line{\hskip 40pt{\tt wiener [options]}\hfill}\par \medskip \noindent where {\tt options} include the following, given in their long (`{\tt --}') as well as their short (`{\tt -}') forms: \medskip \optitem[{\tt --help}, {\tt -h}]{Display a brief help message and exit clean.} \optitem[{\tt --verbose}, {\tt -v}]{Toggle verbose mode. Default: off.} \optitem[{\tt --num\_samples}, {\tt -M} $M$]{Generate $M$ samples of data. Here $M$ should always be an even number, greater than the long lag |KK|. If an odd number is specified, the program will automatically adjust this to the next higher integer. Default: $M=|KK|=100$.} \optitem[{\tt --dimension}, {\tt -D} $D$]{Specifies the dimension $D$ of the Wiener process, that is to say generating a set of $D$ numbers for each of the $M$ data points in the seqence. Default: $D=1$.} \optitem[{\tt --seed}, {\tt -s} $S$]{Define a custom seed number $S$ for the initialization of the random number generator. Default: $S=|DEFAULT_SEED|=310952$.} \optitem[{\tt --uniform}, {\tt -u}]{Instead of generating a sequence of data corresponding to a Wiener process, lock the program to simply generate a uniform distribution of $D$-dimensional points, with each element distributed over the interval $[0,1]$.} \optitem[{\tt --normal}, {\tt -n}]{Instead of generating a sequence of data corresponding to a Wiener process, lock the program to simply generate a normal distribution of $D$-dimensional points, with each element distributed with zero expectation value and unit variance.} \noindent One may look upon the two last options as verification options, generating data suitable for spectral tests on the quality of the generator of pseudo-random numbers. @ Plotting the results using \GNUPLOT. The data sets generated by \wiener\ may be displayed and saved as Encapsulated PostScript images using, say, \GNUPLOT\ or \METAPOST. While I personally prefer MetaPost, mainly due to the possibility of incorporating the same typygraphic elements as in \TeX, many people consider \GNUPLOT\ to be easier in operation. In order to save a scatter graph as Encapsulated PostScript using \GNUPLOT, you may include the following calls in, say, a {\tt Makefile} or a shell script: \bigskip {\obeylines\obeyspaces\tt ~ wiener -D 2 -M 10000 > figure.dat; ~ echo "set term postscript eps;\backslash ~ set output \\"figure.eps\\";\backslash ~ set size square;\backslash ~ set nolabel;\backslash ~ plot \\"figure.dat\\" with dots notitle;\backslash ~ quit" \vertbar\ gnuplot } @ Plotting the results using \METAPOST. Another choice is to go for the \METAPOST\ way. This is illustrated with the following blocks, taken directly from the enclosed Makefile and generating the figures which can be seen in the section relating to the generation of normally distributed variables (routine |normdist|): \bigskip {\obeylines\obeyspaces\tt PROJECT = wiener TEX = tex DVIPS = dvips METAPOST = mpost ~ ~ \# \# Generate the Encapsulated PostScript image fig1.eps for the documentation. \# This is a 2D scatter plot of the uniformly distributed pseudo-random numbers \# prior to having been fed into the Box-Muller transform. \# fig1.eps: Makefile \$(PROJECT).w ~ wiener --uniform -D 2 -M 10000 > fig1.dat; ~ @echo "input graph;{\backslash} ~ ~ def mpdot = btex{\backslash} ~ ~ ~ {{\backslash}vrule height 0.5pt width 1.0pt depth 0.5pt}{\backslash} ~ ~ etex enddef;{\backslash} ~ ~ beginfig(1);{\backslash} ~ ~ draw begingraph(86mm,86mm);{\backslash} ~ ~ setrange(0,0,1,1);{\backslash} ~ ~ pickup pencircle scaled .5pt;{\backslash} ~ ~ gdraw {\backslash}"fig1.dat{\backslash}" plot mpdot;{\backslash} ~ ~ pickup pencircle scaled .25pt;{\backslash} ~ ~ autogrid(itick bot,itick lft);{\backslash} ~ ~ glabel.bot(btex \$\$ x\$\$ etex,OUT);{\backslash} ~ ~ glabel.lft(btex \$\$ y\$\$ etex,OUT);{\backslash} ~ ~ endgraph; endfig; end" > fig1.mp ~ ~ \$(METAPOST) fig1.mp ~ \$(TEX) -jobname=fig1 "{\backslash}input epsf{\backslash}nopagenumbers{\backslash} ~ ~ {\backslash}centerline{{\backslash}epsfbox{fig1.1}}{\backslash}bye" ~ \$(DVIPS) -D1200 -E fig1.dvi -o fig1.eps ~ ~ \# \# Generate the Encapsulated PostScript image fig2.eps for the documentation. \# This is a 2D scatter plot of the normally distributed pseudo-random numbers \# resulting from the Box-Muller transform. \# fig2.eps: Makefile \$(PROJECT).w ~ wiener --normal -D 2 -M 10000 > fig2.dat; ~ @echo "input graph;{\backslash} ~ ~ def mpdot = btex{\backslash} ~ ~ ~ {{\backslash}vrule height 0.5pt width 1.0pt depth 0.5pt}{\backslash} ~ ~ etex enddef;{\backslash} ~ ~ beginfig(1);{\backslash} ~ ~ draw begingraph(86mm,86mm);{\backslash} ~ ~ setrange(whatever,whatever,whatever,whatever);{\backslash} ~ ~ pickup pencircle scaled .5pt;{\backslash} ~ ~ gdraw {\backslash}"fig2.dat{\backslash}" plot mpdot;{\backslash} ~ ~ pickup pencircle scaled .25pt;{\backslash} ~ ~ autogrid(itick bot,itick lft);{\backslash} ~ ~ glabel.bot(btex \$\$ x\$\$ etex,OUT);{\backslash} ~ ~ glabel.lft(btex \$\$ y\$\$ etex,OUT);{\backslash} ~ ~ endgraph; endfig; end" > fig2.mp ~ ~ \$(METAPOST) fig2.mp ~ \$(TEX) -jobname=fig2 "{\backslash}input epsf{\backslash}nopagenumbers{\backslash} ~ ~ {\backslash}centerline{{\backslash}epsfbox{fig2.1}}{\backslash}bye" ~ \$(DVIPS) -D1200 -E fig2.dvi -o fig2.eps ~ ~ \# \# Generate the Encapsulated PostScript image fig3.eps for the documentation. \# This is a 2D graph showing the resulting simulated Wiener process. \# fig3.eps: Makefile \$(PROJECT).w ~ wiener -D 2 -M 10000 > fig3.dat; ~ @echo "input graph;{\backslash} ~ ~ beginfig(1);{\backslash} ~ ~ draw begingraph(86mm,86mm);{\backslash} ~ ~ setrange(whatever,whatever,whatever,whatever);{\backslash} ~ ~ pickup pencircle scaled .5pt;{\backslash} ~ ~ gdraw {\backslash}"fig3.dat{\backslash}";{\backslash} ~ ~ pickup pencircle scaled .25pt;{\backslash} ~ ~ autogrid(itick bot,itick lft);{\backslash} ~ ~ glabel.bot(btex \$\$ x\$\$ etex,OUT);{\backslash} ~ ~ glabel.lft(btex \$\$ y\$\$ etex,OUT);{\backslash} ~ ~ endgraph; endfig; end" > fig3.mp ~ ~ \$(METAPOST) fig3.mp ~ \$(TEX) -jobname=fig3 "{\backslash}input epsf{\backslash}nopagenumbers{\backslash} ~ ~ {\backslash}centerline{{\backslash}epsfbox{fig3.1}}{\backslash}bye" ~ \$(DVIPS) -D1200 -E fig3.dvi -o fig3.eps ~ ~ } @*The main program. Here follows the general outline of the main program. @c @@; @@; @@; @@; int main(int argc, char *argv[]) { @@; @@; @@; @@; @@; @@; return(SUCCESS); } @ Library dependencies. The standard \ANSICEE\ libraries included in this program are:\par \varitem[{\tt math.h}]{For access to common mathematical functions.} \varitem[{\tt stdio.h}]{For file access and any block involving |fprintf|.} \varitem[{\tt stdlib.h}]{For access to the |exit| function.} \varitem[{\tt string.h}]{For string manipulation, |strcpy|, |strcmp| etc.} \varitem[{\tt ctype.h}]{For access to the |isalnum| function.} @= #include #include #include #include #include @ Global definitions. These are the global definitions present in the \wiener\ program:\par \varitem[{\tt VERSION}]{The current program revision number.} \varitem[{\tt COPYRIGHT}]{The copyright banner.} \varitem[{\tt SUCCESS}]{The return code for successful program termination.} \varitem[{\tt FAILURE}]{The return code for unsuccessful program termination.} \varitem[|QUALITY|]{The recommended ``quality level'' for high-resolution use, according to Knuth. Used by the |ranf_array| routine.} \varitem[{\tt KK}]{The ``long lag'' used by routines |ranf_array| and |ranf_start|.} \varitem[{\tt LL}]{The ``short lag'' used by routines |ranf_array| and |ranf_start|.} \varitem[|DEFAULT_SEED|]{The default seed to use when initializing the random number generator, using the routine |ranf_start|. The seed can be hand-picked using the {\tt --seed} command-line option.} \varitem[{\tt cmd\_match(s,l,c)}]{Check if the string |s| and/or |l| matches the string |c|. This short-hand macro is used when parsing the command line for options.} @= #define VERSION "1.0" #define COPYRIGHT "Copyright (C) 2011, Fredrik Jonsson " #define SUCCESS (0) #define FAILURE (1) #define QUALITY (1009) #define KK (100) #define LL (37) #define DEFAULT_SEED (310952) #define cmd_match(s,l,c) ((!strcmp((c),(s)))||(!strcmp((c),(l)))) #define MODE_WIENER_PROCESS (0) #define MODE_LOCKED_UNIFORM_DISTRIBUTION (1) #define MODE_LOCKED_NORMAL_DISTRIBUTION (2) @ Declaration of global variables. Usually, the only global variables allowed in my programs are |optarg|, which is a pointer to the the string of characters that specified the call from the command line, and |progname|, which simply is a pointer to the string containing the name of the program, as it was invoked from the command line. However, as Donald Knuth has a {\sl faiblesse} for global variables, I have for the sake of consistency with the routines related to the random number generator kept his definitions in a global scope. @= extern char *optarg; char *progname; double ranf_arr_buf[QUALITY]; double ranf_arr_dummy=-1.0, ranf_arr_started=-1.0; double *ranf_arr_ptr=&ranf_arr_dummy; /* the next random fraction, or -1 */ @*Declaration of routines. These routines exclusively concern the generation of random numbers and the generation of normally distributed data points in the Wiener process. @= @@; @@; @@; @@; @@; @@; @@; @@; @@; @ Generation of uniformly distributed random numbers. We here avoid relying on the standard functions available in \CEE,\footnote{$\dagger$}{ Here only included as a reference, the (primitive) standard \CEE\ way of generating random integer numbers, in this case initializing with a simple {\tt srand(time(NULL))} and generating random numbers between 0 and |RAND_MAX|, is {\obeylines\obeyspaces\tt int rand\_stdc() { srand(time(NULL)); /* Initialize random number generator */ return(rand()); } } \noindent I have personally {\sl not} (yet) checked this approach using the spectral tests recommended by Knuth.} but rather take resort to the algorithm devised by Donald Knuth in {\sl The Art of Computer Programming, Volume 1 -- Fundamental Algorithms}, 3rd edition, Section 3.6. (Addison--Wesley, Boston, 1998).% \footnote{$\ddagger$}{ The credit for this random number generator, which employs double floating-point precision rather than the original integer version, goes entirely to Donald Knuth and the persons who contributed. The original source code %and full list of credits are available at {\tt http://www-cs-faculty.stanford.edu/{\tilde}knuth/programs/rng-double.c}. The current routine takes into account changes as explained in the errata to the 2nd edition, see {\tt http://www-cs-faculty.stanford.edu/{\tilde}knuth/taocp.html} in the changes to Volume 2 on pages 171 and following.} The variables to the routine are as follows: \varitem[{\tt double aa[]}]{On return, this array contains |n| new random numbers, following the recurrence $X_j=(X_{j-100}-X_{j-37})\mod2^{30}.$ Not used on input.} \varitem[{\tt double n}]{The number |n| of random numbers to be generated. This array length must be at least |KK| elements.} \noindent The |mod_sum(x,y)| macro is here merely a shorthand for ``$(x+y)\mod1.0$.'' @= #define mod_sum(x,y) (((x)+(y))-(int)((x)+(y))) double ran_u[KK]; /* the generator state */ void ranf_array(double aa[], int n) /* put n new random fractions in aa */ { register int i,j; for (j=0;j= #define TT 70 /* guaranteed separation between streams */ #define is_odd(s) ((s)&1) void ranf_start(long seed) { register int t,s,j; double u[KK+KK-1]; double ulp=(1.0/(1L<<30))/(1L<<22); /* 2 to the -52 */ double ss=2.0*ulp*((seed&0x3fffffff)+2); for (j=0;j=1.0) ss-=1.0-2*ulp; /* cyclic shift of 51 bits */ } u[1]+=ulp; /* make u[1] (and only u[1]) "odd" */ for (s=seed&0x3fffffff,t=TT-1; t; ) { for (j=KK-1;j>0;j--) u[j+j]=u[j],u[j+j-1]=0.0; /* "square" */ for (j=KK+KK-2;j>=KK;j--) { u[j-(KK-LL)]=mod_sum(u[j-(KK-LL)],u[j]); u[j-KK]=mod_sum(u[j-KK],u[j]); } if (is_odd(s)) { /* "multiply by z" */ for (j=KK;j>0;j--) u[j]=u[j-1]; u[0]=u[KK]; /* shift the buffer cyclically */ u[LL]=mod_sum(u[LL],u[KK]); } if (s) s>>=1; else t--; } for (j=0;j=0? *ranf_arr_ptr++: ranf_arr_cycle()) double ranf_arr_cycle() { if (ranf_arr_ptr==&ranf_arr_dummy) ranf_start(314159L); /* the user forgot to initialize */ ranf_array(ranf_arr_buf,QUALITY); ranf_arr_buf[KK]=-1; ranf_arr_ptr=ranf_arr_buf+1; return ranf_arr_buf[0]; } @ Generation of normally distributed variables. Accepting uniformly distributed random numbers as input, the Box--Muller transform creates a set of normally distributed numbers. This transform originally appeared in G.~E. P.~Box and Mervin E.~Muller, {\sl A Note on the Generation of Random Normal Deviates}, Annals Math.~Stat. {\bf 29}, 610--611 (1958). In Donald Knuth's {\sl The Art of Computer Programming, Volume 1 -- Fundamental Algorithms}, 3rd edition (Addison--Wesley, Boston, 1998), the Box--Muller method is denoted {the polar method} and is described in detail in Section 3.4.1, Algorithm P, on page 122. \centerline{\epsfbox{fig1.1}} \figcaption[1]{Sample two-dimensional output from the |ranf_matrix| routine, in this case 10\,000 data points uniformly distributed over the domain $0\le\{x,y\}\le 1$. The data for this graph was generated by \wiener\ using {\tt wiener --uniform -D 2 -M 10000 > fig1.dat}. See the {\tt fig1.eps} block in the enclosed {\tt Makefile} for details on how \METAPOST\ was used in the generation of the encapsulated PostScript image of the graph.} \centerline{\epsfbox{fig2.1}} \figcaption[2]{The same data points as in Fig.~1, but after having applied the Box--Muller transform to yield a normal distribution of pseudo-random numbers. The data for this graph was generated by \wiener\ using {\tt wiener --normal -D 2 -M 10000 > fig2.dat}. See the {\tt fig2.eps} block in the enclosed {\tt Makefile} for details on how \METAPOST\ was used in the generation of the encapsulated PostScript image of the graph.} @= #define PI (3.14159265358979323846264338) void normdist(double **aa, int mm, int dd) { register int j, k; register double f, z; for (j=0;j0.0) { f=sqrt(-2*log(aa[k][j])); z=2.0*PI*aa[k+1][j]; aa[k][j]=f*cos(z); aa[k+1][j]=f*sin(z); } else { fprintf(stderr,"%s: Error: Zero element detected!\n",progname); fprintf(stderr,"%s: (row %d, column %d)\n",progname,k,j); } } } return; } #undef PI @ Routine for generation of numerical data describing a Wiener process. This is the core routine of the \wiener\ program. After having initialized the random number generator with the supplied seed value (calling |ranf_start(seed)|), the generation of a sequence of numbers describing a Wiener process starts with the generation of $M\times D$ random and uniformly distributed floating-point numbers ($M\times D/2$ pairs, assuming $M\times D$ to be an even number), from which $M\times D$ normally distributed floating-point numbers are computed using the Box--Muller transform\footnote{$\dagger$}{For example, see {\tt http://en.wikipedia.org/wiki/Box\%E2\%80\%93Muller\_transform}.} The actual computation of the random numbers (at first corresponding to a uniform distribution in the interval $[0,1]$) is done by the routine |ranf_matrix(aa,mm,dd)|, which fills an $[M\times D]$ array. \centerline{\epsfbox{fig3.1}} \figcaption[3]{The same data points as in Fig.~2, but after having chained the normally distributed points to form the simulated Wiener process. The data for this graph was generated by \wiener\ using {\tt wiener -D 2 -M 10000 > fig3.dat}. The trajectory starts with data point 1 at $(0,0)$ and end up with data point 10\,000 at approximately $(89.9,12.6)$ See the {\tt fig3.eps} block in the enclosed {\tt Makefile} for details on how \METAPOST\ was used in the generation of the encapsulated PostScript image of the graph.} The variables interfacing the |wiener| routine are as follows: \medskip \varitem[|aa|]{[Output] The $M\times D$ matrix $A$, on return containing the $M$ number of $D$-dimensional data points for the generated Wiener process.} \varitem[|mm|]{[Input] The $M$ number of data points to generate. This equals to the number of rows in the |aa| array, and must be at least |KK| elements.} \varitem[|dd|]{[Input] The dimension of each generated data point.} \varitem[|seed|]{[Input] The seed to use when initializing the random number generator, using the routine |ranf_start|.} \varitem[|mode|]{[Input] Determines if the sequence should be locked to simply generate a uniform distribution over the interval $[0,1]$ (|MODE_LOCKED_UNIFORM_DISTRIBUTION|) or a normal (Gaussian) distribution with expectation value zero and unit variance (|MODE_LOCKED_NORMAL_DISTRIBUTION|). Otherwise, the series of data will be generated to simulate a Wiener process, as is the main purpose of the \wiener\ program. One may look upon the two first modes as verification options, generating data suitable for spectral tests on the quality of the generator of pseudo-random numbers.} @= void wiener(double **aa, int mm, int dd, int seed, short mode) { register int j, k; ranf_start(seed); ranf_matrix(aa,mm,dd); /* Uniform distribution over $[0,1]$ */ if (mode==MODE_LOCKED_UNIFORM_DISTRIBUTION) return; normdist(aa, mm, dd); /* Normal distribution of unit variance around zero */ if (mode==MODE_LOCKED_NORMAL_DISTRIBUTION) return; for (j=0;j= double **dmatrix(long nrl, long nrh, long ncl, long nch) { long i, nrow=nrh-nrl+1,ncol=nch-ncl+1; double **m; m=(double **) malloc((size_t)((nrow+1)*sizeof(double*))); if (!m) { fprintf(stderr,"%s: Allocation failure 1 in dmatrix()\n",progname); exit(FAILURE); } m += 1; m -= nrl; m[nrl]=(double*) malloc((size_t)((nrow*ncol+1)*sizeof(double))); if (!m[nrl]) { fprintf(stderr,"%s: Allocation failure 2 in dmatrix()\n",progname); exit(FAILURE); } m[nrl] += 1; m[nrl] -= ncl; for(i=nrl+1;i<=nrh;i++) m[i]=m[i-1]+ncol; return m; } @ Memory de-allocation. The |free_dmatrix| routine {\sl releases} the memory occupied by the double floating-point precision matrix |v[nl..nh]|, as allocated by |dmatrix()|. @= void free_dmatrix(double **m, long nrl, long nrh, long ncl, long nch) { free((char*) (m[nrl]+ncl-1)); free((char*) (m+nrl-1)); } @ Displaying a help message at terminal output. @= void display_help_message(void) { fprintf(stderr,"Usage: %s M [options]\n",progname); fprintf(stderr,"Options:\n"); fprintf(stderr," -h, --help\n"); fprintf(stderr," Display this help message and exit clean.\n"); fprintf(stderr," -v, --verbose\n"); fprintf(stderr," Toggle verbose mode. Default: off.\n"); fprintf(stderr," -M, --num_samples \n"); fprintf(stderr," Generate M samples of data. Here M should always be\n"); fprintf(stderr," an even number, greater than the long lag KK=%d.\n",KK); fprintf(stderr," If an odd number is specified, the program will\n"); fprintf(stderr," automatically adjust this to the next higher\n"); fprintf(stderr," integer. Default: M=%d.\n",KK); fprintf(stderr," -D, --dimension \n"); fprintf(stderr," Generate D-dimensional samples of data. Default: D=1.\n"); fprintf(stderr," -s, --seed \n"); fprintf(stderr," Define a custom seed number for the initialization\n"); fprintf(stderr," of the random number generator. Default: seed=%d.\n", DEFAULT_SEED); fprintf(stderr," -u, --uniform\n"); fprintf(stderr," Instead of generating a sequence of D-dimensional\n"); fprintf(stderr," corresponding to a Wiener process, lock the program\n"); fprintf(stderr," to simply generate a uniform distribution of\n"); fprintf(stderr," D-dimensional points, with each element distributed\n"); fprintf(stderr," over the interval [0,1].\n"); fprintf(stderr," -n, --normal\n"); fprintf(stderr," Instead of generating a sequence of D-dimensional\n"); fprintf(stderr," corresponding to a Wiener process, lock the program\n"); fprintf(stderr," to simply generate a normal distribution of\n"); fprintf(stderr," D-dimensional points, with each element distributed\n"); fprintf(stderr," with zero expectation value and unit variance.\n"); } @ Checking for a valid path character. The |pathcharacter| routine takes one character |ch| as argument, and returns 1 (``true'') if the character is valid character of a path string, otherwise 0 (``false'') is returned. @= short pathcharacter(int ch) { return(isalnum(ch)||(ch=='.')||(ch=='/')||(ch=='\\')||(ch=='_')|| (ch=='-')||(ch=='+')); } @ Stripping path string from a file name. The |strip_away_path| routine takes a character string |filename| as argument, and returns a pointer to the same string but without any preceding path segments. This routine is, for example, useful for removing paths from program names as parsed from the command line. @= char *strip_away_path(char filename[]) { int j,k=0; while (pathcharacter(filename[k])) k++; j=(--k); /* this is the uppermost index of the full path+file string */ while (isalnum((int)(filename[j]))) j--; j++; /* this is the lowermost index of the stripped file name */ return(&filename[j]); } @ Declaration of local variables of the |main| program. In \CWEB\ one has the option of adding variables along the program, for example by locally adding temporary variables related to a given sub-block of code. However, the philosophy in the \wiener\ program is to keep all variables of the |main| section collected together, so as to simplify tasks as, for example, tracking down a given variable type definition. Exceptions to this rule are dummy variables merely used for the iteration over loops, not participating in any other tasks. The local variables of the program are as follows: \medskip \varitem[|aa|]{The $M\times D$ matrix $A$, containing the $M$ number of $D$-dimensional data points for the generated Wiener process.} \varitem[|mm|]{The $M$ number of data points to generate. This equals to the number of rows in the |aa| array (of dimension $[M\times D]$), and must be at least |KK| elements. The default initialization is $|mm|=|KK|$; however this may change depending on supplied command-line parameters.} \varitem[|dd|]{The dimension of each generated data point. Default value: 1.} \varitem[|seed|]{The seed to use when initializing the random number generator, using the routine |ranf_start|. The seed can be hand-picked using the {\tt --seed} command-line option. Default value: 310952.} \varitem[|mode|]{Determines if the sequence should be locked to simply generate a uniform distribution over the interval $[0,1]$ (|MODE_LOCKED_UNIFORM_DISTRIBUTION|) or a normal (Gaussian) distribution with expectation value zero and unit variance (|MODE_LOCKED_NORMAL_DISTRIBUTION|). Otherwise, the series of data will be generated to simulate a Wiener process, as is the main purpose of the \wiener\ program. One may look upon the two first modes as verification options, generating data suitable for spectral tests on the quality of the generator of pseudo-random numbers.} @= double **aa; unsigned long mm=KK; unsigned dd=1; int seed=DEFAULT_SEED; short mode=MODE_WIENER_PROCESS, verbose=0; int no_arg; register int j, k; @ Parse command line for parameters. We here use the possibility open in \CWEB\ to add {\tt getopt.h} to the inclusions of libraries, as this block is the only one making use of the definitions therein. @= progname=strip_away_path(argv[0]); no_arg=argc; while(--argc) { if(cmd_match("-h", "--help", argv[no_arg-argc])) { display_help_message(); exit(SUCCESS); } else if(cmd_match("-v", "--verbose", argv[no_arg-argc])) { verbose=(verbose?0:1); } else if(cmd_match("-M", "--num_samples", argv[no_arg-argc])) { --argc; if(!sscanf(argv[no_arg-argc],"%lu",&mm)) { fprintf(stderr,"%s: Error detected when parsing the number of " "samples.\n", progname); display_help_message(); exit(FAILURE); } if (mm= KK = %d.\n", progname, KK); exit(FAILURE); } if (is_odd(mm)) mm++; /* If odd, then make it even */ } else if(cmd_match("-D", "--dimension", argv[no_arg-argc])) { --argc; if(!sscanf(argv[no_arg-argc],"%ud",&dd)) { fprintf(stderr,"%s: Error detected when parsing dimension.\n", progname); display_help_message(); exit(FAILURE); } if (dd<1) { fprintf(stderr,"%s: Dimension D should be at least 1.\n",progname); exit(FAILURE); } } else if(cmd_match("-s", "--seed", argv[no_arg-argc])) { --argc; if(!sscanf(argv[no_arg-argc],"%d",&seed)) { fprintf(stderr,"%s: Error detected when parsing the seed of the " "initializer.\n", progname); display_help_message(); exit(FAILURE); } } else if(cmd_match("-u", "--uniform", argv[no_arg-argc])) { mode=MODE_LOCKED_UNIFORM_DISTRIBUTION; } else if(cmd_match("-n", "--normal", argv[no_arg-argc])) { mode=MODE_LOCKED_NORMAL_DISTRIBUTION; } else { fprintf(stderr,"%s: Sorry, I do not recognize option '%s'.\n", progname, argv[no_arg-argc]); display_help_message(); exit(FAILURE); } } if (verbose) fprintf(stdout,"This is %s v.%s. %s\n",progname,VERSION,COPYRIGHT); @ Allocate memory for a vector containing $M\times D$ elements. We here make a call to the operating system for the allocation of a sufficcient amount of memory to accommodate $M$ data points, each of dimension $D$. We here apply the common convention in \CEE, starting the indexing of the allocated array at zero. @= aa=dmatrix(0,mm-1,0,dd-1); @ Fill vector with $M$ number of $D$:tuples describing the Wiener process. @= wiener(aa,mm,dd,seed,mode); @ Print the generated vector at standard terminal output. @= for (k=0;k= free_dmatrix(aa,0,mm-1,0,dd-1); @*Index.