1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
|
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
typedef enum { use_2pn, use_1_2pn, use_2pm_2pn, use_1_2pm_2pn } tr_flag;
typedef struct {
double factor;
int n;
int m;
int nsign;
tr_flag flag;
} table_record;
static int g_records;
static int g_null;
static table_record g_tr[8192];
static int compare_record( const void *a, const void *b ) {
if( ((table_record*)a)->factor > ((table_record*)b)->factor )
return 1;
if( ((table_record*)a)->factor < ((table_record*)b)->factor )
return -1;
return 0;
}
static int get_record_score( int record )
{
switch( g_tr[record].flag ) {
case use_2pn: case use_1_2pn: return 3;
default: return 4;
}
}
static void add_record( double v, int n, int nsign, int m, tr_flag flag )
{
g_tr[g_records].factor = v;
g_tr[g_records].n = n;
g_tr[g_records].m = m;
g_tr[g_records].nsign = nsign;
g_tr[g_records].flag = flag;
++g_records;
}
static void dump_record( int record )
{
if( record >= g_records ) return; // should warn?
double factor = g_tr[record].factor;
int n = g_tr[record].n;
int m = g_tr[record].m;
int nsign = g_tr[record].nsign;
tr_flag flag = g_tr[record].flag;
switch( flag ) {
case use_2pn:
printf( "%0+25.21lf x = x %s %2d\n", factor, n>=0?"<<":">>", (int)abs(n) );
break;
case use_1_2pn:
printf( "%0+25.21lf x%s= x %s %2d\n", factor, nsign==-1?"-":"+", n>=0?"<<":">>", (int)abs(n) );
break;
case use_2pm_2pn:
case use_1_2pm_2pn:
printf( "%0+25.21lf x%s= ( x %s %2d ) %s ( x %s %2d )\n", factor, flag==use_2pm_2pn?" ":"-", m>=0?"<<":">>", (int)abs(m), nsign==-1?"-":"+", n>=0?"<<":">>", (int)abs(n));
break;
}
}
static void find_tail( int off, int record, int score ) {
static int trail[10240];
int i;
// printf( "%d %d %d\n", off, record, score );
trail[off] = record;
// My naive implementation found a solution with score 51
// Don't look deeper when it's worse
if( score > 54 ) return;
if( g_tr[record].factor >= -0.0001220703125 ) {
printf( "%d\n", score );
for( i=0; i<=off; ++i)
dump_record( trail[i] );
return;
}
for( i = record + 1; i < g_null; ++i )
if( ( g_tr[record].factor / g_tr[i].factor < 2.0 ) && ( g_tr[record].factor / g_tr[i].factor >= 1.6 ) )
find_tail( off + 1, i, get_record_score( record ) + score );
// Retry each entry once
if( off && ( trail[off-1] != record ) )
find_tail( off+1, record, get_record_score( record ) + score );
};
int main() {
int n, m, nsign;
// loop over all constructs in the form +-1*x^+-2n, n=-31..+31
for (n=-21; n<= 21; ++n )
for ( nsign=-1; nsign<=1; nsign+=2 )
{
// The one term only case
double v = (double)nsign * pow( 2., (double)n );
if( v > 0.0 )
add_record( log(v), n, nsign, m, use_2pn );
// Loop over second term
for (m=-13; m<=13; ++m )
{
double v = pow( 2.0, (double)m ) + (double)nsign * pow( 2., (double)n );
if( v <= 0.0 )
continue;
add_record( log(v), n, nsign, m, m ? use_2pm_2pn : use_1_2pn );
if( v < 1.0 )
add_record( log(1.0 - v), n, nsign, m, use_1_2pm_2pn );
}
}
// sort records by value
qsort( g_tr, g_records, sizeof(*g_tr), compare_record );
printf( "All potential solutions:\n" );
for (n=0; n<g_records; ++n )
dump_record(n);
printf( "Now looking for the optima of the algorithm for negative exp().\nThis can take a while... \n" );
// Remove redundant entries
for (n=1; n<g_records; ++n ) {
if (fabs(g_tr[n].factor) / fabs(g_tr[n-1].factor) > 0.999 ) {
if (g_tr[n-1].flag == use_2pn || g_tr[n-1].flag == use_1_2pn || g_tr[n].factor > 0.0 )
memmove(g_tr+n,g_tr+n+1,(g_records-n-1)*sizeof(*g_tr));
else
memmove(g_tr+n-1,g_tr+n,(g_records-n-1)*sizeof(*g_tr));
--g_records; --n;
}
}
for (g_null=0; g_tr[g_null].factor <= 0.0; ++g_null);
// Now find all paths through factors ensuring that the gap between two
// steps is not larger than factor 2.0
// Tail terminates when we reach 0.0001220703125 and prints results
for (n=0; n<g_records; ++n ) {
if (g_tr[n].factor <= -0.693 && g_tr[n].factor >= -0.694 )
find_tail( 0, n, 0);
}
return 0;
}
|
