test portability of PRNG

testlwrandom.c

/**********************************************************************
 *
 * PostGIS - Spatial Types for PostgreSQL
 * http://postgis.net
 *
 * PostGIS 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.
 *
 * PostGIS 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 PostGIS.  If not, see <http://www.gnu.org/licenses/>.
 *
 **********************************************************************
 *
 * Copyright 2019 Mike Taves
 *
 **********************************************************************/

#include <stdlib.h>
#include <stdio.h>
#include <stdint.h>

#include "lwrandom.h"


/* See asa183_prb.c and asa183_prb_output.txt by John Burkardt
 * https://people.sc.fsu.edu/~jburkardt/c_src/asa183/asa183.html
 */

void test10 ( )
{
	int i;
	double r;

	printf ( "\n" );
	printf ( "TEST10\n" );
	printf ( "  R8_UNI computes pseudorandom values.\n" );
	printf ( "  Two seeds, S1 and S2, are used.\n" );

	_lwrandom_set_seeds(12345, 34567);

	printf ( "\n" );
	printf ( "      R                     S1            S2\n" );
	printf ( "\n" );
	printf ( "                  %12d  %12d\n", _lwrandom_get_seed(1), _lwrandom_get_seed(2) );

	for ( i = 0; i < 10; i++ )
	{
		r = lwrandom_uniform();
		printf ( "  %14f  %12d  %12d\n", r, _lwrandom_get_seed(1), _lwrandom_get_seed(2) );
	}
	return;
}

void test11 ( )
{
	int i;
	int n = 100000;
	double *u;
	double u_avg;
	double u_var;

	u = ( double * ) malloc ( n * sizeof ( double ) );

	printf ( "\n" );
	printf ( "TEST11\n" );
	printf ( "  Examine the average and variance of a\n" );
	printf ( "  sequence generated by R8_UNI.\n" );

	_lwrandom_set_seeds(12345, 34567);

	printf ( "\n" );
	printf ( "  Now compute %d elements.\n", n );

	u_avg = 0.0;
	for ( i = 0; i < n; i++ )
	{
		u[i] = lwrandom_uniform();
		u_avg += u[i];
	}

	u_avg = u_avg / ( float ) ( n );

	u_var = 0.0;
	for ( i = 0; i < n; i++ )
	{
		u_var += ( u[i] - u_avg ) * ( u[i] - u_avg );
	}

	u_var = u_var / ( float ) ( n - 1 );

	printf ( "\n" );
	printf ( "  Average value = %f\n", u_avg );
	printf ( "  Expecting       %f\n", 0.5 );;

	printf ( "\n" );
	printf ( "  Variance =      %f\n", u_var );
	printf ( "  Expecting       %f\n", 1.0 / 12.0 );

	free ( u );

	return;
}

void test12 ( )
{
	int i;
	double r;
	int32_t s1;
	int32_t s2;
	int32_t s1_save;
	int32_t s2_save;

	printf ( "\n" );
	printf ( "TEST12\n" );
	printf ( "  Show how the seeds used by R8_UNI,\n" );
	printf ( "  which change on each step, can be reset to\n" );
	printf ( "  restore any part of the sequence.\n" );

	s1_save = 12345;
	s2_save = 34567;

	_lwrandom_set_seeds(s1_save, s2_save);

	printf ( "\n" );
	printf ( "  Begin sequence with following seeds\n" );
	printf ( "\n" );
	printf ( "  S1 = %d\n", _lwrandom_get_seed(1) );
	printf ( "  S2 = %d\n", _lwrandom_get_seed(2) );
	printf ( "\n" );
	printf ( "         I      R                     S1            S2\n" );
	printf ( "\n" );

	for ( i = 1; i <= 10; i++ )
	{
		r = lwrandom_uniform();
		printf ( "  %8d  %14f  %12d  %12d\n", i, r, _lwrandom_get_seed(1), _lwrandom_get_seed(2) );

		if ( i == 5 )
		{
			s1_save = _lwrandom_get_seed(1);
			s2_save = _lwrandom_get_seed(2);
		}
	}
	_lwrandom_set_seeds(s1_save, s2_save);

	printf ( "\n" );
	printf ( "  Restart the sequence, using the seeds\n" );
	printf ( "  produced after step 5:\n" );
	printf ( "\n" );
	printf ( "  S1 = %d\n", _lwrandom_get_seed(1) );
	printf ( "  S2 = %d\n", _lwrandom_get_seed(2) );
	printf ( "\n" );
	printf ( "         I      R                     S1            S2\n" );
	printf ( "\n" );

	for ( i = 1; i <= 10; i++ )
	{
		r = lwrandom_uniform();
		printf ( "  %8d  %14f  %12d  %12d\n", i, r, _lwrandom_get_seed(1), _lwrandom_get_seed(2) );
	}
	return;
}

/* Compare with expected output excerpts from:
* https://people.sc.fsu.edu/~jburkardt/c_src/asa183/asa183_prb_output.txt
*/
int main()
{
	test10();
	test11();
	test12();
	return 0;
}