RanLip Maple toolbox 

School of Engineering and IT, Deakin University, Australia
Copyright Gleb Beliakov, 2006
Australia
gbeliako@gmail.com 

Installation instructions 

To use this worksheet you need to have MapleRanLip.dll library in the curent directory or on the path 

Introduction 

This worksheet illustrates how RanLip package can be used with Maple (under MS Windows).   

RanLip is a programming library to generate non-uniformly distributed random vectors with Lipschitz continuous densities. 

 

Please see the user manual supplied with this worksheet. 

This worksheet is provided as a template for user's own distributions. It shows how to declare the distribution density and  

generate random vectors. 

It should be used in conjunction with the supplied dll (MapleRanLip.dll), and the user manual (ranlipmaple.pdf). 

The typical calling sequence is as follows: 

1. execute define_external procedures in the next subsection 

2. declare your density function (and also set the boundaries of the domain) 

3. call the relevant  routines with correct arguments
 

Please change the RANLIPPATH line below to the location of this package on your computer! 

Initialization: Declarations of RanLip procedures found in the dll 

>	restart; with(plots):

 

Execute the commands below to declare the external subroutines from RanLip library. 

>

RANLIPPATH:="mapleranlip.dll":   # if you need to provide full path, use, e.g.,RANLIPPATH:="c:/work/maple/mapleranlip.dll": PrepareHatFunction := define_external('MWRAP_PrepareHatFunction','MAPLE', LIB=RANLIPPATH): InitRanLip := define_external('MWRAP_InitRanLip','MAPLE', LIB=RANLIPPATH): PrepareHatFunctionRanLip := define_external('MWRAP_PrepareHatFunctionRanLip','MAPLE', LIB=RANLIPPATH): PrepareHatFunctionAuto := define_external('MWRAP_PrepareHatFunctionAuto','MAPLE', LIB=RANLIPPATH): RandomVec := define_external('MWRAP_RandomVecRanLip','MAPLE', LIB=RANLIPPATH): GetSeed := define_external('MWRAP_GetSeedRanLip','MAPLE', LIB=RANLIPPATH): Lipschitz := define_external('MWRAP_LipschitzRanLip','MAPLE', LIB=RANLIPPATH): Count_total := define_external('MWRAP_Count_totalRanLip','MAPLE', LIB=RANLIPPATH): Count_error := define_external('MWRAP_Count_errorRanLip','MAPLE', LIB=RANLIPPATH): SavePartition := define_external('MWRAP_SavePartitionRanlip','MAPLE', LIB=RANLIPPATH): LoadPartition := define_external('MWRAP_LoadPartitionRanLip','MAPLE', LIB=RANLIPPATH): Seed := define_external('MWRAP_SeedRanLip','MAPLE', LIB=RANLIPPATH): SetDistFunction := define_external('MWRAP_SetDistFunctionRanLip','MAPLE', LIB=RANLIPPATH): FreeMem := define_external('MWRAP_FreeMemRanLip','MAPLE', LIB=RANLIPPATH):

 

 

Section 1: Declaration of the density function  

We need to declare the distribution density function and also box constraints lo_i <= x_i <= up_i. 

Declaration of the density function f(x) 

Declare f as a procedure. The objective function f takes n arguments, the components of vector x. It returns a real value.
The density does not need to be normalised. 

>	f := proc( x, y)  # user's density function        exp(-(y-x^2)^2 - (x^2+y^2)/2);    end proc;

 

 

>	# test the density function:  evaluate f(x), plot it f(0.2,0.6); plot3d(f(x,y),x=-3..3,y=-2..3, axes=boxed,orientation=[25,50]);

 

 

 

Declaring the bounds and space to receive generated random vectors 

Boundaries, and also the dimension of x and declaration of xrand.  

>	Lo:=Vector(1..2,datatype=float[8],[-2,-3]):   Up:=Vector(1..2,datatype=float[8],[3,3]): dim:=2; xrand:=Vector(1..2,datatype=float[8]);      # will receive generated random vector

 

 

 

Executing RanLip methods 

Construct the hat function, using a grid 50 x 50, and using 16 x 16 function evaluations in each cell.  Lipschitz constant is 10. 

>	PrepareHatFunction(2,Lo, Up,f, 50,16,10);

 

Section 2: Generation of random vectors 

Generate and plot random vectors with density f. 

Generate a few random vectors 

>	RandomVec(xrand); xrand; for i from 1 to 5 do    RandomVec(xrand); print(xrand); od:

 

 

 

 

 

 

 

Generating many random vectors and plotting them 

First prepare space to hold random vectors 

>	M:=3000: XR:=Matrix(1..M,1..2,datatype=float[8],order=C_order); RandomVec(XR,M);

 

 

>	with(plots): pointplot({seq([XR[i,1],XR[i,2]],i=1..M)},symbol=point);

 

 

Get some statistics from the process: how many random numbers generated, how may errors (error is the hat function too small - increase Lipschitz constant) 

>	Count_total(); Count_error();

 

 

 

Section 3: Saving and retrieving the hat function 

>	SavePartition(`partition.txt`);

 

Hat function successfully saved

 

The next command clears the hat function. No generation of random vectors is possible after it, until we load it or prepare a new hat function . 

>	FreeMem(); RandomVec(xrand); print(xrand); # meaningless result RandomVec(xrand); print(xrand);

 

 

 

>	LoadPartition(`partition.txt`); SetDistFunction(f); # do not forget this command and also define the vector xrand xrand:=Vector(1..2,datatype=float[8]);

 

Hat function  loaded 

 

Now we can generate random vectors again (we keep the same distribution function) 

>	Seed(25); for i from 1 to 5 do    RandomVec(xrand); print(xrand); od:

 

 

 

 

 

 

>	`The seed was `; GetSeed();

 

 

 

>

 

Section 4: Another example, putting it all together 

>	g := proc( x, y)  # user's density function       local r,r1:       r:=(x^2+y^2)^(0.5):       r1:= (x+0.2)^2+(y+0.2)^2:      (r-1)^2*exp(-r1/3);    end proc:

 

>	g(x,y); g(0,0);

 

 

 

>	plot3d(g(x,y),x=-4..4,y=-4..4, axes=boxed,orientation=[25,50],grid=[50,50]);

 

 

just define the boundaries and prepare arrays 

>	dim:=2: Lo:=Vector(1..dim,datatype=float[8],[-4,-4]):   Up:=Vector(1..dim,datatype=float[8],[4,4]): xrand:=Vector(1..dim,datatype=float[8]):

 

generate the hat function 

>	PrepareHatFunctionAuto(dim,Lo, Up,g, 50,16); `Computed Lipschitz constant is `;  Lipschitz();

 

 

 

generate random vectors 

>	RandomVec(xrand); xrand;

 

 

generate and plot many random vectors 

>	M:=10000: XR:=Matrix(1..M,1..dim,datatype=float[8],order=C_order): RandomVec(XR,M); pointplot({seq([XR[i,1],XR[i,2]],i=1..M)},symbol=point,axes=boxed); `random numbers generated =`;Count_total(); `errors=`;Count_error();

 

 

 

 

 

 

>

 

References 

Please consult http://www.deakin.edu.au/~gleb for updates. 

We appreciate your comments about your experiences using RanLip and this Maple toolbox. Please send your feedback to 

gbeliako@gmail.com 

This project was carried out with the assistance of Ilya Khriapin, which is greatefully acknowledged. 

Legal Notice: The copyright for this application is owned by the author(s). Neither Maplesoft nor the author are responsible for any errors contained within and are not liable for any damages resulting from the use of this material. This application is intended for non-commercial, non-profit use only. Contact the author for permission if you wish to use this application in for-profit activities.