/* COMPUTE NATURAL LOGARITHM OF 3/2 */
/* 03/09/13 (dkc) */
/* */
/*****************************************************************************/
#include <stdio.h>
#include <math.h>
unsigned int carry(unsigned int a, unsigned int b, unsigned int sum);
unsigned int lmbd(unsigned int mode, unsigned int a);
void div128_64(unsigned int a0, unsigned int a1, unsigned int a2,
unsigned int a3, unsigned int *quotient, unsigned int d2,
unsigned int d3);
void setn(unsigned int *a, unsigned int b, unsigned int n);
void lshiftn(unsigned int *a, unsigned int *b, unsigned int shift, unsigned int n);
void addn(unsigned int *a, unsigned int *b, unsigned int n);
void copyn(unsigned int *a, unsigned int *b, unsigned int n);
void div512_32(unsigned int a0, unsigned int a1, unsigned int a2,
unsigned int a3, unsigned int a4, unsigned int a5,
unsigned int a6, unsigned int a7, unsigned int a8,
unsigned int a9, unsigned int a10, unsigned int a11,
unsigned int a12, unsigned int a13, unsigned int a14,
unsigned int a15, unsigned int *quotient, unsigned int d15);
int main ()
{
unsigned int i,j,n;
unsigned int A[16],B[16],C[16],D[4],E[4],T[16],U[16],L[16],S[16];
double ln2,ln3,temp;
FILE *Outfp;
Outfp = fopen("output.dat","w");
ln2=log(2.0);
ln3=log(3.0);
printf("log(2)=%e \n",ln2);
printf("log(3)=%e \n",ln3);
printf("\n");
n=4;
setn(A,0,n);
A[0]=0x20000000;
setn(C,0,n);
setn(T,1,n);
setn(U,1,n);
for (i=1; i<1000000; i++) {
if (i!=(i/2)*2) {
copyn(T,S,n);
addn(T,T,n);
addn(S,T,n);
if (T[1]!=0) {
printf("iterations=%d \n",i);
break;
}
setn(D,0,n);
for (j=0; j<i; j++)
addn(T,D,n);
}
else {
copyn(U,S,n);
addn(U,U,n);
addn(S,U,n);
if (U[1]!=0) {
printf("iterations=%d \n",i);
break;
}
setn(E,0,n);
for (j=0; j<i; j++)
addn(U,E,n);
}
if (i!=(i/2)*2) {
div128_64(A[0],A[1],A[2],A[3],B,D[2],D[3]);
div128_64(B[0],B[1],B[2],B[3],B,U[2],U[3]);
}
else {
div128_64(A[0],A[1],A[2],A[3],B,E[2],E[3]);
div128_64(B[0],B[1],B[2],B[3],B,T[2],T[3]);
}
addn(B,C,n);
}
printf("log(3/2)= \n");
for (i=0; i<4; i++) {
printf(" %#010x \n",C[i]);
fprintf(Outfp," %#010x \n",C[i]);
}
L[0]=0x162e42fe;
L[1]=0xfa39ef35;
L[2]=0x793c7673;
L[3]=0x007e5ed5;
//
temp=(double)(L[0])/65536.0/8192.0;
printf("\n");
printf("log(2)=%e \n",temp);
addn(C,L,n);
printf("\n");
printf("log(3)= \n");
for (i=0; i<4; i++) {
printf(" %#010x \n",L[i]);
}
temp=(double)(L[0])/65536.0/8192.0;
printf("log(3)=%e \n",temp);
//
// check values
//
n=16;
printf("\n");
printf("compare results to OEIS (to about four words) \n");
// log(3)=1.09861 228866 810969 139524 523692 252570 464749
setn(A,464749,n);
//
// next word
//
setn(B,252570,n);
lshiftn(B,T,19,n);
addn(T,A,n);
lshiftn(B,T,18,n);
addn(T,A,n);
lshiftn(B,T,17,n);
addn(T,A,n);
lshiftn(B,T,16,n);
addn(T,A,n);
lshiftn(B,T,14,n);
addn(T,A,n);
lshiftn(B,T,9,n);
addn(T,A,n);
lshiftn(B,T,6,n);
addn(T,A,n);
//
// next word
//
setn(C,0,n);
setn(B,523692,n);
lshiftn(B,T,19,n);
addn(T,C,n);
lshiftn(B,T,18,n);
addn(T,C,n);
lshiftn(B,T,17,n);
addn(T,C,n);
lshiftn(B,T,16,n);
addn(T,C,n);
lshiftn(B,T,14,n);
addn(T,C,n);
lshiftn(B,T,9,n);
addn(T,C,n);
lshiftn(B,T,6,n);
addn(T,C,n);
//
lshiftn(C,T,19,n);
addn(T,A,n);
lshiftn(C,T,18,n);
addn(T,A,n);
lshiftn(C,T,17,n);
addn(T,A,n);
lshiftn(C,T,16,n);
addn(T,A,n);
lshiftn(C,T,14,n);
addn(T,A,n);
lshiftn(C,T,9,n);
addn(T,A,n);
lshiftn(C,T,6,n);
addn(T,A,n);
//
// next word
//
setn(C,0,n);
setn(B,139524,n);
lshiftn(B,T,19,n);
addn(T,C,n);
lshiftn(B,T,18,n);
addn(T,C,n);
lshiftn(B,T,17,n);
addn(T,C,n);
lshiftn(B,T,16,n);
addn(T,C,n);
lshiftn(B,T,14,n);
addn(T,C,n);
lshiftn(B,T,9,n);
addn(T,C,n);
lshiftn(B,T,6,n);
addn(T,C,n);
//
setn(U,0,n);
lshiftn(C,T,19,n);
addn(T,U,n);
lshiftn(C,T,18,n);
addn(T,U,n);
lshiftn(C,T,17,n);
addn(T,U,n);
lshiftn(C,T,16,n);
addn(T,U,n);
lshiftn(C,T,14,n);
addn(T,U,n);
lshiftn(C,T,9,n);
addn(T,U,n);
lshiftn(C,T,6,n);
addn(T,U,n);
//
lshiftn(U,T,19,n);
addn(T,A,n);
lshiftn(U,T,18,n);
addn(T,A,n);
lshiftn(U,T,17,n);
addn(T,A,n);
lshiftn(U,T,16,n);
addn(T,A,n);
lshiftn(U,T,14,n);
addn(T,A,n);
lshiftn(U,T,9,n);
addn(T,A,n);
lshiftn(U,T,6,n);
addn(T,A,n);
//
// next word
//
setn(C,0,n);
setn(B,810969,n);
lshiftn(B,T,19,n);
addn(T,C,n);
lshiftn(B,T,18,n);
addn(T,C,n);
lshiftn(B,T,17,n);
addn(T,C,n);
lshiftn(B,T,16,n);
addn(T,C,n);
lshiftn(B,T,14,n);
addn(T,C,n);
lshiftn(B,T,9,n);
addn(T,C,n);
lshiftn(B,T,6,n);
addn(T,C,n);
//
setn(U,0,n);
lshiftn(C,T,19,n);
addn(T,U,n);
lshiftn(C,T,18,n);
addn(T,U,n);
lshiftn(C,T,17,n);
addn(T,U,n);
lshiftn(C,T,16,n);
addn(T,U,n);
lshiftn(C,T,14,n);
addn(T,U,n);
lshiftn(C,T,9,n);
addn(T,U,n);
lshiftn(C,T,6,n);
addn(T,U,n);
//
setn(C,0,n);
lshiftn(U,T,19,n);
addn(T,C,n);
lshiftn(U,T,18,n);
addn(T,C,n);
lshiftn(U,T,17,n);
addn(T,C,n);
lshiftn(U,T,16,n);
addn(T,C,n);
lshiftn(U,T,14,n);
addn(T,C,n);
lshiftn(U,T,9,n);
addn(T,C,n);
lshiftn(U,T,6,n);
addn(T,C,n);
//
lshiftn(C,T,19,n);
addn(T,A,n);
lshiftn(C,T,18,n);
addn(T,A,n);
lshiftn(C,T,17,n);
addn(T,A,n);
lshiftn(C,T,16,n);
addn(T,A,n);
lshiftn(C,T,14,n);
addn(T,A,n);
lshiftn(C,T,9,n);
addn(T,A,n);
lshiftn(C,T,6,n);
addn(T,A,n);
//
// next word
//
setn(C,0,n);
setn(B,228866,n);
lshiftn(B,T,19,n);
addn(T,C,n);
lshiftn(B,T,18,n);
addn(T,C,n);
lshiftn(B,T,17,n);
addn(T,C,n);
lshiftn(B,T,16,n);
addn(T,C,n);
lshiftn(B,T,14,n);
addn(T,C,n);
lshiftn(B,T,9,n);
addn(T,C,n);
lshiftn(B,T,6,n);
addn(T,C,n);
//
setn(U,0,n);
lshiftn(C,T,19,n);
addn(T,U,n);
lshiftn(C,T,18,n);
addn(T,U,n);
lshiftn(C,T,17,n);
addn(T,U,n);
lshiftn(C,T,16,n);
addn(T,U,n);
lshiftn(C,T,14,n);
addn(T,U,n);
lshiftn(C,T,9,n);
addn(T,U,n);
lshiftn(C,T,6,n);
addn(T,U,n);
//
setn(C,0,n);
lshiftn(U,T,19,n);
addn(T,C,n);
lshiftn(U,T,18,n);
addn(T,C,n);
lshiftn(U,T,17,n);
addn(T,C,n);
lshiftn(U,T,16,n);
addn(T,C,n);
lshiftn(U,T,14,n);
addn(T,C,n);
lshiftn(U,T,9,n);
addn(T,C,n);
lshiftn(U,T,6,n);
addn(T,C,n);
//
setn(U,0,n);
lshiftn(C,T,19,n);
addn(T,U,n);
lshiftn(C,T,18,n);
addn(T,U,n);
lshiftn(C,T,17,n);
addn(T,U,n);
lshiftn(C,T,16,n);
addn(T,U,n);
lshiftn(C,T,14,n);
addn(T,U,n);
lshiftn(C,T,9,n);
addn(T,U,n);
lshiftn(C,T,6,n);
addn(T,U,n);
//
lshiftn(U,T,19,n);
addn(T,A,n);
lshiftn(U,T,18,n);
addn(T,A,n);
lshiftn(U,T,17,n);
addn(T,A,n);
lshiftn(U,T,16,n);
addn(T,A,n);
lshiftn(U,T,14,n);
addn(T,A,n);
lshiftn(U,T,9,n);
addn(T,A,n);
lshiftn(U,T,6,n);
addn(T,A,n);
//
// next word
//
setn(C,0,n);
setn(B,109861,n);
lshiftn(B,T,19,n);
addn(T,C,n);
lshiftn(B,T,18,n);
addn(T,C,n);
lshiftn(B,T,17,n);
addn(T,C,n);
lshiftn(B,T,16,n);
addn(T,C,n);
lshiftn(B,T,14,n);
addn(T,C,n);
lshiftn(B,T,9,n);
addn(T,C,n);
lshiftn(B,T,6,n);
addn(T,C,n);
//
setn(U,0,n);
lshiftn(C,T,19,n);
addn(T,U,n);
lshiftn(C,T,18,n);
addn(T,U,n);
lshiftn(C,T,17,n);
addn(T,U,n);
lshiftn(C,T,16,n);
addn(T,U,n);
lshiftn(C,T,14,n);
addn(T,U,n);
lshiftn(C,T,9,n);
addn(T,U,n);
lshiftn(C,T,6,n);
addn(T,U,n);
//
setn(C,0,n);
lshiftn(U,T,19,n);
addn(T,C,n);
lshiftn(U,T,18,n);
addn(T,C,n);
lshiftn(U,T,17,n);
addn(T,C,n);
lshiftn(U,T,16,n);
addn(T,C,n);
lshiftn(U,T,14,n);
addn(T,C,n);
lshiftn(U,T,9,n);
addn(T,C,n);
lshiftn(U,T,6,n);
addn(T,C,n);
//
setn(U,0,n);
lshiftn(C,T,19,n);
addn(T,U,n);
lshiftn(C,T,18,n);
addn(T,U,n);
lshiftn(C,T,17,n);
addn(T,U,n);
lshiftn(C,T,16,n);
addn(T,U,n);
lshiftn(C,T,14,n);
addn(T,U,n);
lshiftn(C,T,9,n);
addn(T,U,n);
lshiftn(C,T,6,n);
addn(T,U,n);
//
setn(C,0,n);
lshiftn(U,T,19,n);
addn(T,C,n);
lshiftn(U,T,18,n);
addn(T,C,n);
lshiftn(U,T,17,n);
addn(T,C,n);
lshiftn(U,T,16,n);
addn(T,C,n);
lshiftn(U,T,14,n);
addn(T,C,n);
lshiftn(U,T,9,n);
addn(T,C,n);
lshiftn(U,T,6,n);
addn(T,C,n);
//
lshiftn(C,T,19,n);
addn(T,A,n);
lshiftn(C,T,18,n);
addn(T,A,n);
lshiftn(C,T,17,n);
addn(T,A,n);
lshiftn(C,T,16,n);
addn(T,A,n);
lshiftn(C,T,14,n);
addn(T,A,n);
lshiftn(C,T,9,n);
addn(T,A,n);
lshiftn(C,T,6,n);
addn(T,A,n);
//
lshiftn(A,T,3,n);
lshiftn(A,U,1,n);
addn(T,U,n);
lshiftn(U,C,196,n);
for (i=0; i<7; i++) {
div512_32(C[0],C[1],C[2],C[3],C[4],C[5],C[6],C[7],C[8],C[9],C[10],C[11],
C[12],C[13],C[14],C[15],A,1000000);
copyn(A,C,n);
}
lshiftn(C,A,288+25,n);
for (i=0; i<n; i++)
printf(" %#010x \n",A[i]);
fclose(Outfp);
return(0);
}