Sorry your browser is not supported!

#include "stdafx.h"
#include "stdio.h"
#include <math.h>
#include "Math.h"
#include "stdlib.h"
#include "time.h"
#include "vector"
#include "..Headersglobstruct.h"

GlobStruct* g_pGlob = NULL;

#define NEGATIVE	(int) -1
#define ZERO		(int) 0
#define POSITIVE	(int) 1

char hexs[]={"0123456789ABCDEF"};

#define _BYTE(x)				((x) & 0xFF)

#define int2sll(X)	(((sll) (X)) << 32)
#define sll2int(X)	((long) ((X) >> 32))
#define sllint(X)	((X) & (ull) 0xffffffff00000000)
#define sllfrac(X)	((X) & (ull) 0x00000000ffffffff)
#define sllneg(X)	(-(X))
#define _slladd(X,Y)	((X) + (Y))
#define _sllsub(X,Y)	((X) - (Y))

/* Constants (converted from double) */
#define CONST_0		(sll) 0x0000000000000000
#define CONST_1		(sll) 0x0000000100000000
#define CONST_2		(sll) 0x0000000200000000
#define CONST_3		(sll) 0x0000000300000000
#define CONST_4		(sll) 0x0000000400000000
#define CONST_10	(sll) 0x0000000a00000000
#define CONST_1_2	(sll) 0x0000000080000000
#define CONST_1_3	(sll) 0x0000000055555555
#define CONST_1_4	(sll) 0x0000000040000000
#define CONST_1_5	(sll) 0x0000000033333333
#define CONST_1_6	(sll) 0x000000002aaaaaaa
#define CONST_1_7	(sll) 0x0000000024924924
#define CONST_1_8	(sll) 0x0000000020000000
#define CONST_1_9	(sll) 0x000000001c71c71c
#define CONST_1_10	(sll) 0x0000000019999999
#define CONST_1_11	(sll) 0x000000001745d174
#define CONST_1_12	(sll) 0x0000000015555555
#define CONST_1_20	(sll) 0x000000000ccccccc
#define CONST_1_30	(sll) 0x0000000008888888
#define CONST_1_42	(sll) 0x0000000006186186
#define CONST_1_56	(sll) 0x0000000004924924
#define CONST_1_72	(sll) 0x00000000038e38e3
#define CONST_1_90	(sll) 0x0000000002d82d82
#define CONST_1_110	(sll) 0x000000000253c825
#define CONST_1_132	(sll) 0x0000000001f07c1f
#define CONST_1_156	(sll) 0x0000000001a41a41
#define CONST_E		(sll) 0x00000002b7e15162
#define CONST_1_E	(sll) 0x000000005e2d58d8
#define CONST_SQRTE	(sll) 0x00000001a61298e1
#define CONST_1_SQRTE	(sll) 0x000000009b4597e3
#define CONST_LOG2_E	(sll) 0x0000000171547652
#define CONST_LOG10_E	(sll) 0x000000006f2dec54
#define CONST_LN2	(sll) 0x00000000b17217f7
#define CONST_LN10	(sll) 0x000000024d763776
#define CONST_PI	(sll) 0x00000003243f6a88
#define CONST_PI_2	(sll) 0x00000001921fb544
#define CONST_PI_4	(sll) 0x00000000c90fdaa2
#define CONST_1_PI	(sll) 0x00000000517cc1b7
#define CONST_2_PI	(sll) 0x00000000a2f9836e
#define CONST_2_SQRTPI	(sll) 0x0000000120dd7504
#define CONST_SQRT2	(sll) 0x000000016a09e667
#define CONST_1_SQRT2	(sll) 0x00000000b504f333

unsigned int seed;

BOOL APIENTRY DllMain( HANDLE hModule, 
                       DWORD  ul_reason_for_call, 
                       LPVOID lpReserved
					 )
{
	srand((unsigned) time(NULL));
	seed=time(NULL);
    switch (ul_reason_for_call)
	{
		case DLL_PROCESS_ATTACH:
		case DLL_THREAD_ATTACH:
		case DLL_THREAD_DETACH:
		case DLL_PROCESS_DETACH:
			break;
    }
    return TRUE;
}

void Constructor ( void )
{
	// Create memory here
}

void Destructor ( void )
{
	// Free memory here
}

void ReceiveCoreDataPtr ( LPVOID pCore )
{
	// Get Core Data Pointer here
	g_pGlob = (GlobStruct*)pCore;
}

// This is an example of an exported variable
__int64 hextointLI(LPSTR hex)
{
register int pos;
__int64 base,total;

	total=0;
	if ((hex) && strlen(hex)>0)
	{
		base=1;		
		for (pos=0; pos<(int) strlen(hex); pos++)
		{
		register char one;
		register int loop;

			one=toupper(*(hex+(strlen(hex)-1-pos)));
			for (loop=0; loop<(int) strlen(hexs); loop++)
			{
				if (one==hexs[loop])
				{
					// A valid number found
					total+=(loop*base);
					break;
				}
			}

			base<<=4;
		}
	}

	return (total);

}

int hextoint(LPSTR hex)
{
register int pos,base,total;

	total=0;
	if ((hex) && strlen(hex)>0)
	{
		base=1;		
		for (pos=0; pos<(int) strlen(hex); pos++)
		{
		register char one;
		register int loop;

			one=toupper(*(hex+(strlen(hex)-1-pos)));
			for (loop=0; loop<(int) strlen(hexs); loop++)
			{
				if (one==hexs[loop])
				{
					// A valid number found
					total+=(loop*base);
					break;
				}
			}

			base<<=4;
		}
	}

	return (total);
}	

double PI(void)
{
	return ((double) 3.14159265358979);
}

double PHI(void)
{
	return ((double) 1.6180339887);
}

double PHI2(void)
{
	return ((double) 0.6180339);
}

unsigned __int64 factorial(unsigned __int64 f)
{
register unsigned __int64 total;

	total=(f==0 || f==1 ? 1 : f);
	while (f>1)
	{
		total*=(f-1);
		f--;
	}
	return (total);
}

DWORD valF(LPSTR string)
{
float value;

	value=(float) atof(string);
	return (*(DWORD *) &value);
}

DWORD valF(LPSTR string,DWORD dp)
{
float value;

	value=(float) valD(string,dp);
	return (*(DWORD *) &value);
}

double valD(LPSTR string)
{
	return (atof(string));
}

double valD(LPSTR string,DWORD dp)
{
char *buffer;
int decimal,sign;
char *buffer2;
double value;

	value=0.0;
	buffer = _fcvt(atof(string), dp, &decimal, &sign );
	ZeroMemory(&buffer2,sizeof(buffer2));
	if ((buffer2=(char *) calloc(1,1024))!=NULL)
	{
		if (sign!=0)
		{
			strcat(buffer2,"-");
		}

		if (decimal>0)
		{
			strncat(buffer2,buffer,decimal);
			if (dp)
			{
				strcat(buffer2,".");
			}
		}
		else
		{
			if (dp)
			{
				strcat(buffer2,".");
			}
			decimal=0;
		}
		
		strcat(buffer2,(char *) (buffer+decimal));
		value=atof(buffer2);
	}

	
	free(buffer2);
	return (value);
}

int sgn(int value)
{
	return (value<0 ? NEGATIVE : 
			value>0 ? POSITIVE : ZERO);
}

int sgn(LONGLONG value)
{
	return (value<0 ? NEGATIVE : 
			value>0 ? POSITIVE : ZERO);
}

int sgn(float value)
{
	return (value<0.0 ? NEGATIVE : 
			value>0.0 ? POSITIVE : ZERO);
}

int sgn(double value)
{
	return (value<0.0 ? NEGATIVE : 
			value>0.0 ? POSITIVE : ZERO);
}

double besselFO0(double value)
{
	return (_j0(value));
}

double besselFO1(double value)
{
	return (_j1(value));
}

double besselFON(int value1,double value2)
{
	return (_jn(value1,value2));
}

double besselSO0(double value)
{
	return (_y0(value));
}

double besselSO1(double value)
{
	return (_y1(value));
}

double besselSON(int value1,double value2)
{
	return (_yn(value1,value2));
}

double _log(double value)
{
	return (log(value));
}

double _log10(double value)
{
	return (log10(value));
}

DWORD rrand(DWORD min,DWORD max)
{
register DWORD range;

	srand(seed*(unsigned int) time(NULL)); 
	range=max-min+1;
	seed=rand();	
	return ((seed % range + min));
}

double rndf(double min,double max)
{
register double value;

	srand(seed*(unsigned int) time(NULL)); 	
	seed=rand();
	value=(double) ((double) seed/(double) RAND_MAX);
	return ((double) (value*(max-min)+min));
}

DWORD strf(DWORD pOldString,double value)
{
LPSTR buffer=NULL;
const int size=256;

	if(pOldString) g_pGlob->CreateDeleteString ( (DWORD*)&pOldString, 0 );
	g_pGlob->CreateDeleteString ( (DWORD*)&buffer, size);

	ZeroMemory(buffer,size);
	sprintf(buffer,"%.15lf",value);
	return ((DWORD) buffer);
}

DWORD bytes16ToBottomDWORD(DWORD part2,DWORD part1)
{
	return ((_BYTE(part2)<<8) | (_BYTE(part1)));
}

DWORD bytes16ToTopDWORD(DWORD part2,DWORD part1)
{
	return ((_BYTE(part2)<<24) | (_BYTE(part1)<<16));
}

DWORD bytes32ToDWORD(DWORD part4,DWORD part3,DWORD part2,DWORD part1)
{
	return ((_BYTE(part4)<<24) | (_BYTE(part3)<<16) | (_BYTE(part2)<<8) | (_BYTE(part1)));
}

DWORD DWORDToByte(DWORD value,DWORD bank)
{
DWORD val;

	val=(bank==0 ? value : 
		 bank==1 ? value>>8 : 
		 bank==2 ? value>>16 : value>>24);
	return (_BYTE(val));
}

DWORD fcmp(float value1,float value2,float epsilon)
{
	return (fabs(value1-value2)<=fabs(value1)*epsilon ? (DWORD) true : (DWORD) false);
}

DWORD dfcmp(double value1,double value2,double epsilon)
{
	return (fabs(value1-value2)<=fabs(value1)*epsilon ? (DWORD) true : (DWORD) false);
}

sll slladd(sll x, sll y)
{
	return (x + y);
}

sll sllsub(sll x, sll y)
{
	return (x - y);
}

/*
 * Let a = A * 2^32 + a_hi * 2^0 + a_lo * 2^(-32)
 * Let b = B * 2^32 + b_hi * 2^0 + b_lo * 2^(-32)
 *
 * Where:
 *   *_hi is the integer part
 *   *_lo the fractional part
 *   A and B are the sign (0 for positive, -1 for negative).
 *
 * a * b = (A * 2^32 + a_hi * 2^0 + a_lo * 2^-32)
 *       * (B * 2^32 + b_hi * 2^0 + b_lo * 2^-32)
 *
 * Expanding the terms, we get:
 *
 *	 = A * B * 2^64 + A * b_h * 2^32 + A * b_l * 2^0
 *	 + a_h * B * 2^32 + a_h * b_h * 2^0 + a_h * b_l * 2^-32
 *	 + a_l * B * 2^0 + a_l * b_h * 2^-32 + a_l * b_l * 2^-64
 *
 * Grouping by powers of 2, we get:
 *
 *	 = A * B * 2^64
 *	 Meaningless overflow from sign extension - ignore
 *
 *	 + (A * b_h + a_h * B) * 2^32
 *	 Overflow which we can't handle - ignore
 *
 *	 + (A * b_l + a_h * b_h + a_l * B) * 2^0
 *	 We only need the low 32 bits of this term, as the rest is overflow
 *
 *	 + (a_h * b_l + a_l * b_h) * 2^-32
 *	 We need all 64 bits of this term
 *
 *	 +  a_l * b_l * 2^-64
 *	 We only need the high 32 bits of this term, as the rest is underflow
 *
 * Note that:
 *   a > 0 && b > 0: A =  0, B =  0 and the third term is a_h * b_h
 *   a < 0 && b > 0: A = -1, B =  0 and the third term is a_h * b_h - b_l
 *   a > 0 && b < 0: A =  0, B = -1 and the third term is a_h * b_h - a_l
 *   a < 0 && b < 0: A = -1, B = -1 and the third term is a_h * b_h - a_l - b_l
 */
/*
sll sllmul(sll left, sll right)
{
register sll retval;

	__asm__(
		"# sllmulnt"
		"	movl	%1, %%eaxnt"
		"	mull 	%3nt"
		"	movl	%%edx, %%ebxnt"
		"nt"
		"	movl	%2, %%eaxnt"
		"	mull 	%4nt"
		"	movl	%%eax, %%ecxnt"
		"nt"
		"	movl	%1, %%eaxnt"
		"	mull	%4nt"
		"	addl	%%eax, %%ebxnt"
		"	adcl	%%edx, %%ecxnt"
		"nt"
		"	movl	%2, %%eaxnt"
		"	mull	%3nt"
		"	addl	%%ebx, %%eaxnt"
		"	adcl	%%ecx, %%edxnt"
		"nt"
		"	btl	$31, %2nt"
		"	jnc	1fnt"
		"	subl	%3, %%edxnt"
		"1:	btl	$31, %4nt"
		"	jnc	1fnt"
		"	subl	%1, %%edxnt"
		"1:nt"
		: "=&A" (retval)
		: "m" (left), "m" (((unsigned *) &left)[1]),
		  "m" (right), "m" (((unsigned *) &right)[1])
		: "ebx", "ecx", "cc"
	);
	return retval;
}
*/
/* Plain C version: not optimal but portable. */
sll sllmul(sll a, sll b)
{
	sll a_lo, b_lo;
	sll a_hi, b_hi;
	sll x;

	a_lo = a;
	a_hi = (ull) a >> 32;
	b_lo = b;
	b_hi = (ull) b >> 32;

	x = ((ull) (a_hi * b_hi) << 32)
	  + (((ull) a_lo * b_lo) >> 32)
	  + (sll) a_lo * b_hi
	  + (sll) b_lo * a_hi;

	return x;
}

sll sllinv(sll v)
{
int sgn = 0;
sll u;
ull s = -1;

	/* Use positive numbers, or the approximation won't work */
	if (v < CONST_0) {
		v = sllneg(v);
		sgn = 1;
	}

	/* An approximation - must be larger than the actual value */
	for (u = v; u; u >>= 1)
		s >>= 1;

	/* Newton's Method */
	u = sllmul(s, _sllsub(CONST_2, sllmul(v, s)));
	u = sllmul(u, _sllsub(CONST_2, sllmul(v, u)));
	u = sllmul(u, _sllsub(CONST_2, sllmul(v, u)));
	u = sllmul(u, _sllsub(CONST_2, sllmul(v, u)));
	u = sllmul(u, _sllsub(CONST_2, sllmul(v, u)));
	u = sllmul(u, _sllsub(CONST_2, sllmul(v, u)));

	return ((sgn) ? sllneg(u): u);
}

sll slldiv(sll left, sll right)
{
	return sllmul(left, sllinv(right));
}

sll sllmul2(sll x)
{
	return x << 1;
}

sll sllmul4(sll x)
{
	return x << 2;
}

sll sllmul2n(sll x, int n)
{
sll y;

	y = x << n;
	return y;
}

sll slldiv2(sll x)
{
	return x >> 1;
}

sll slldiv4(sll x)
{
	return x >> 2;
}

sll slldiv2n(sll x, int n)
{
sll y;

	y = x >> n;
	return y;
}

/*
 * Unpack the IEEE floating point double format and put it in fixed point
 * sll format.
 */
sll dbl2sll(double dbl)
{
union {
	double d;
	unsigned u[2];
	ull ull;
	sll sll;
} in, retval;
register unsigned exp;

	/* Move into memory as args might be passed in regs */
	in.d = dbl;

	/* Leading 1 is assumed by IEEE */
	retval.u[1] = 0x40000000;

	/* Unpack the mantissa into the unsigned long */
	retval.u[1] |= (in.u[1] << 10) & 0x3ffffc00;
	retval.u[1] |= (in.u[0] >> 22) & 0x000003ff;
	retval.u[0] = in.u[0] << 10;

	/* Extract the exponent and align the decimals */
	exp = (in.u[1] >> 20) & 0x7ff;
	if (exp)
		retval.ull >>= 1053 - exp;
	else
		return 0L;

	/* Negate if negative flag set */
	if (in.u[1] & 0x80000000)
		retval.sll = -retval.sll;

	return retval.sll;
}

double sll2dbl(sll s)
{
union {
	double d;
	unsigned u[2];
	ull ull;
	sll sll;
} in, retval;
register unsigned int exp;
register unsigned int flag;

	if (s == 0)
		return 0.0;

	/* Move into memory as args might be passed in regs */
	in.sll = s;

	/* Handle the negative flag */
	if (in.sll < 1) {
		flag = 0x80000000;
		in.ull = sllneg(in.sll);
	} else
		flag = 0x00000000;

	/* Normalize */
	for (exp = 1053; in.ull && (in.u[1] & 0x80000000) == 0; exp--) {
		in.ull <<= 1;
	}
	in.ull <<= 1;
	exp++;
	in.ull >>= 12;
	retval.ull = in.ull;
	retval.u[1] |= flag | (exp << 20);

	return retval.d;
}

/*
 * Calculate cos x where -pi/4 <= x <= pi/4
 *
 * Description:
 *	cos x = 1 - x^2 / 2! + x^4 / 4! - ... + x^(2N) / (2N)!
 *	Note that (pi/4)^12 / 12! < 2^-32 which is the smallest possible number.
 */
sll _sllcos(sll x)
{
sll retval, x2;

	x2 = sllmul(x, x);
	/*
	 * cos x = t0 + t1 + t2 + t3 + t4 + t5 + t6
	 *
	 * f0 =  0! =  1
	 * f1 =  2! =  2 *  1 * f0 =   2 * f0
	 * f2 =  4! =  4 *  3 * f1 =  12 x f1
	 * f3 =  6! =  6 *  5 * f2 =  30 * f2
	 * f4 =  8! =  8 *  7 * f3 =  56 * f3
	 * f5 = 10! = 10 *  9 * f4 =  90 * f4
	 * f6 = 12! = 12 * 11 * f5 = 132 * f5
	 *
	 * t0 = 1
	 * t1 = -t0 * x2 /   2 = -t0 * x2 * CONST_1_2
	 * t2 = -t1 * x2 /  12 = -t1 * x2 * CONST_1_12
	 * t3 = -t2 * x2 /  30 = -t2 * x2 * CONST_1_30
	 * t4 = -t3 * x2 /  56 = -t3 * x2 * CONST_1_56
	 * t5 = -t4 * x2 /  90 = -t4 * x2 * CONST_1_90
	 * t6 = -t5 * x2 / 132 = -t5 * x2 * CONST_1_132
	 */
	retval = _sllsub(CONST_1, sllmul(x2, CONST_1_132));
	retval = _sllsub(CONST_1, sllmul(sllmul(x2, retval), CONST_1_90));
	retval = _sllsub(CONST_1, sllmul(sllmul(x2, retval), CONST_1_56));
	retval = _sllsub(CONST_1, sllmul(sllmul(x2, retval), CONST_1_30));
	retval = _sllsub(CONST_1, sllmul(sllmul(x2, retval), CONST_1_12));
	retval = _sllsub(CONST_1, slldiv2(sllmul(x2, retval)));
	return retval;
}

/*
 * Calculate sin x where -pi/4 <= x <= pi/4
 *
 * Description:
 *	sin x = x - x^3 / 3! + x^5 / 5! - ... + x^(2N+1) / (2N+1)!
 *	Note that (pi/4)^13 / 13! < 2^-32 which is the smallest possible number.
 */
sll _sllsin(sll x)
{
sll retval, x2;

	x2 = sllmul(x, x);
	/*
	 * sin x = t0 + t1 + t2 + t3 + t4 + t5 + t6
	 *
	 * f0 =  0! =  1
	 * f1 =  3! =  3 *  2 * f0 =   6 * f0
	 * f2 =  5! =  5 *  4 * f1 =  20 x f1
	 * f3 =  7! =  7 *  6 * f2 =  42 * f2
	 * f4 =  9! =  9 *  8 * f3 =  72 * f3
	 * f5 = 11! = 11 * 10 * f4 = 110 * f4
	 * f6 = 13! = 13 * 12 * f5 = 156 * f5
	 *
	 * t0 = 1
	 * t1 = -t0 * x2 /   6 = -t0 * x2 * CONST_1_6
	 * t2 = -t1 * x2 /  20 = -t1 * x2 * CONST_1_20
	 * t3 = -t2 * x2 /  42 = -t2 * x2 * CONST_1_42
	 * t4 = -t3 * x2 /  72 = -t3 * x2 * CONST_1_72
	 * t5 = -t4 * x2 / 110 = -t4 * x2 * CONST_1_110
	 * t6 = -t5 * x2 / 156 = -t5 * x2 * CONST_1_156
	 */
	retval = _sllsub(x, sllmul(x2, CONST_1_156));
	retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_110));
	retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_72));
	retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_42));
	retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_20));
	retval = _sllsub(x, sllmul(sllmul(x2, retval), CONST_1_6));
	return retval;
}

sll sllcos(sll x)
{
int i;
sll retval;

	/* Calculate cos (x - i * pi/2), where -pi/4 <= x - i * pi/2 <= pi/4  */
	i = sll2int(_slladd(sllmul(x, CONST_2_PI), CONST_1_2));
	x = _sllsub(x, sllmul(int2sll(i), CONST_PI_2));

	switch (i & 3) {
		default:
		case 0:
			retval = _sllcos(x);
			break;
		case 1:
			retval = sllneg(_sllsin(x));
			break;
		case 2:
			retval = sllneg(_sllcos(x));
			break;
		case 3:
			retval = _sllsin(x);
			break;
	}
	return retval;
}

sll sllsin(sll x)
{
int i;
sll retval;

	/* Calculate sin (x - n * pi/2), where -pi/4 <= x - i * pi/2 <= pi/4 */
	i = sll2int(_slladd(sllmul(x, CONST_2_PI), CONST_1_2));
	x = _sllsub(x, sllmul(int2sll(i), CONST_PI_2));

	switch (i & 3) {
		default:
		case 0:
			retval = _sllsin(x);
			break;
		case 1:
			retval = _sllcos(x);
			break;
		case 2:
			retval = sllneg(_sllsin(x));
			break;
		case 3:
			retval = sllneg(_sllcos(x));
			break;
	}
	return retval;
}

sll slltan(sll x)
{
int i;

	sll retval;
	i = sll2int(_slladd(sllmul(x, CONST_2_PI), CONST_1_2));
	x = _sllsub(x, sllmul(int2sll(i), CONST_PI_2));
	switch (i & 3) {
		default:
		case 0:
		case 2:
			retval = slldiv(_sllsin(x), _sllcos(x));
			break;
		case 1:
		case 3:
			retval = sllneg(slldiv(_sllcos(x), _sllsin(x)));
			break;
	}
	return retval;
}

/*
 * atan x = SUM[n=0,) (-1)^n * x^(2n + 1)/(2n + 1), |x| < 1
 *
 * Two term approximation
 *	a = x - x^3/3
 * Gives us
 *	atan x = a + ??
 * Let ?? = arctan ?
 *	atan x = a + arctan ?
 * Rearrange
 *	atan x - a = arctan ?
 * Apply tan to both sides
 *	tan (atan x - a) = tan arctan ?
 *	tan (atan x - a) = ?
 * Applying the standard formula
 *	tan (u - v) = (tan u - tan v) / (1 + tan u * tan v)
 * Gives us
 *	tan (atan x - a) = (tan atan x - tan a) / (1 + tan arctan x * tan a)
 * Let t = tan a
 *	tan (atan x - a) = (x - t) / (1 + x * t)
 * So finally
 *	arctan x = a + arctan ((tan x - t) / (1 + x * t))
 * And the typical worst case is x = 1.0 which converges in 3 iterations.
 */
sll _sllatan(sll x)
{
sll a, t, retval;

	/* First iteration */
	a = sllmul(x, _sllsub(CONST_1, sllmul(x, sllmul(x, CONST_1_3))));
	retval = a;

	/* Second iteration */
	t = slldiv(_sllsin(a), _sllcos(a));
	x = slldiv(_sllsub(x, t), _slladd(CONST_1, sllmul(t, x)));
	a = sllmul(x, _sllsub(CONST_1, sllmul(x, sllmul(x, CONST_1_3))));
	retval = _slladd(retval, a);

	/* Third  iteration */
	t = slldiv(_sllsin(a), _sllcos(a));
	x = slldiv(_sllsub(x, t), _slladd(CONST_1, sllmul(t, x)));
	a = sllmul(x, _sllsub(CONST_1, sllmul(x, sllmul(x, CONST_1_3))));
	return _slladd(retval, a);
}

sll sllatan(sll x)
{
sll retval;

	if (x < -sllneg(CONST_1))
		retval = sllneg(CONST_PI_2);
	else if (x > CONST_1)
		retval = CONST_PI_2;
	else
		return _sllatan(x);
	return _sllsub(retval, _sllatan(sllinv(x)));
}

/*
 * Calculate e^x where -0.5 <= x <= 0.5
 *
 * Description:
 *	e^x = x^0 / 0! + x^1 / 1! + ... + x^N / N!
 *	Note that 0.5^11 / 11! < 2^-32 which is the smallest possible number.
 */
sll _sllexp(sll x)
{
	sll retval;
	retval = _slladd(CONST_1, sllmul(0, sllmul(x, CONST_1_11)));
	retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_11)));
	retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_10)));
	retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_9)));
	retval = _slladd(CONST_1, sllmul(retval, slldiv2n(x, 3)));
	retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_7)));
	retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_6)));
	retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_5)));
	retval = _slladd(CONST_1, sllmul(retval, slldiv4(x)));
	retval = _slladd(CONST_1, sllmul(retval, sllmul(x, CONST_1_3)));
	retval = _slladd(CONST_1, sllmul(retval, slldiv2(x)));
	return retval;
}

/*
 * Calculate e^x where x is arbitrary
 */
sll sllexp(sll x)
{
int i;

	sll e, retval;

	e = CONST_E;

	/* -0.5 <= x <= 0.5  */
	i = sll2int(_slladd(x, CONST_1_2));
	retval = _sllexp(_sllsub(x, int2sll(i)));

	/* i >= 0 */
	if (i < 0) {
		i = -i;
		e = CONST_1_E;
	}

	/* Scale the result */
	for (;i; i >>= 1) {
		if (i & 1)
			retval = sllmul(retval, e);
		e = sllmul(e, e);
	}
	return retval;
}

/*
 * Calculate natural logarithm using Netwton-Raphson method
 */
sll slllog(sll x)
{
	sll x1, ln = 0;

	/* Scale: e^(-1/2) <= x <= e^(1/2) */
	while (x < CONST_1_SQRTE) {
		ln = _sllsub(ln, CONST_1);
		x = sllmul(x, CONST_E);
	}
	while (x > CONST_SQRTE) {
		ln = _slladd(ln, CONST_1);
		x = sllmul(x, CONST_1_E);
	}

	/* First iteration */
	x1 = sllmul(_sllsub(x, CONST_1), slldiv2(_sllsub(x, CONST_3)));
	ln = _sllsub(ln, x1);
	x = sllmul(x, _sllexp(x1));

	/* Second iteration */
	x1 = sllmul(_sllsub(x, CONST_1), slldiv2(_sllsub(x, CONST_3)));
	ln = _sllsub(ln, x1);
	x = sllmul(x, _sllexp(x1));

	/* Third iteration */
	x1 = sllmul(_sllsub(x, CONST_1), slldiv2(_sllsub(x, CONST_3)));
	ln = _sllsub(ln, x1);

	return ln;
}

/*
 * ln x^y = y * log x
 * e^(ln x^y) = e^(y * log x)
 * x^y = e^(y * ln x)
 */
sll sllpow(sll x, sll y)
{
	if (y == CONST_0)
		return CONST_1;
	return sllexp(sllmul(y, slllog(x)));
}

/*
 * Consider a parabola centered on the y-axis
 * 	y = a * x^2 + b
 * Has zeros (y = 0)  at
 *	a * x^2 + b = 0
 *	a * x^2 = -b
 *	x^2 = -b / a
 *	x = +- (-b / a)^(1 / 2)
 * Letting a = 1 and b = -X
 *	y = x^2 - X
 *	x = +- X^(1 / 2)
 * Which is convenient since we want to find the square root of X, and we can
 * use Newton's Method to find the zeros of any f(x)
 *	xn = x - f(x) / f'(x)
 * Applied Newton's Method to our parabola
 *	f(x) = x^2 - X
 *	xn = x - (x^2 - X) / (2 * x)
 *	xn = x - (x - X / x) / 2
 * To make this converge quickly, we scale X so that
 *	X = 4^N * z
 * Taking the roots of both sides
 *	X^(1 / 2) = (4^n * z)^(1 / 2)
 *	X^(1 / 2) = 2^n * z^(1 / 2)
 * Let N = 2^n
 *	x^(1 / 2) = N * z^(1 / 2)
 * We want this to converge to the positive root, so we must start at a point
 *	0 < start <= x^(1 / 2)
 * or
 *	x^(1/2) <= start <= infinity
 * since
 *	(1/2)^(1/2) = 0.707
 *	2^(1/2) = 1.414
 * A good choice is 1 which lies in the middle, and takes 4 iterations to
 * converge from either extreme.
 */
sll sllsqrt(sll x)
{
	sll n, xn;
       
	/* Start with a scaling factor of 1 */
	n = CONST_1;

	/* Quick solutions for the simple cases */
	if (x <= CONST_0 || x == CONST_1)
		return x;

	/* Scale x so that 0.5 <= x < 2 */
	while (x >= CONST_2) {
		x = slldiv4(x);
		n = sllmul2(n);
	}
	while (x < CONST_1_2) {
		x = sllmul4(x);
		n = slldiv2(n);
	}

	/* Simple solution if x = 4^n */
	if (x == CONST_1)
		return n;

	/* The starting point */
	xn = CONST_1;

	/* Four iterations will be enough */
	xn = _sllsub(xn, slldiv2(_sllsub(xn, slldiv(x, xn))));
	xn = _sllsub(xn, slldiv2(_sllsub(xn, slldiv(x, xn))));
	xn = _sllsub(xn, slldiv2(_sllsub(xn, slldiv(x, xn))));
	xn = _sllsub(xn, slldiv2(_sllsub(xn, slldiv(x, xn))));

	/* Scale the result */
	return sllmul(n, xn);
}

DWORD fixtostring(DWORD pOldString,ull value)
{
LPSTR buffer=NULL;
char buffer2[256];
double _value;

	if(pOldString) g_pGlob->CreateDeleteString ( (DWORD*)&pOldString, 0 );
	g_pGlob->CreateDeleteString ( (DWORD*)&buffer, sizeof(buffer2));

	ZeroMemory(buffer,sizeof(buffer2));
	_value=sll2dbl(value);
	sprintf(buffer2,"%.4lf",_value);
	strcpy(buffer,buffer2);
	return ((DWORD) buffer);
}