Skip to content

Commit b382e42

Browse files
totalspectrumtotalspectrum
authored
Merge pull request #490 from Wuerfel21/W21-softcordic
Add hardware-accurate QEXP and QLOG
2 parents 429d0e9 + 1d5b241 commit b382e42

File tree

5 files changed

+181
-1
lines changed

5 files changed

+181
-1
lines changed

‎Makefile‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -129,7 +129,7 @@ LEXHEADERS = $(BUILD)/spin.tab.h $(BUILD)/basic.tab.h $(BUILD)/cgram.tab.h ast.h
129129

130130
PROGS = $(BUILD)/testlex$(EXT) $(BUILD)/spin2cpp$(EXT) $(BUILD)/flexspin$(EXT) $(BUILD)/flexcc$(EXT)
131131

132-
UTIL = dofmt.c flexbuf.c lltoa_prec.c strupr.c strrev.c strdupcat.c to_utf8.c from_utf8.c sha256.c
132+
UTIL = dofmt.c flexbuf.c lltoa_prec.c strupr.c strrev.c strdupcat.c to_utf8.c from_utf8.c sha256.c softcordic.c
133133

134134
MCPP = directive.c expand.c mbchar.c mcpp_eval.c mcpp_main.c mcpp_system.c mcpp_support.c
135135

‎expr.c‎

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,7 @@
99
#include "spinc.h"
1010
#include "backends/bytecode/bc_bedata.h"
1111
#include "becommon.h" // for PrintExpr
12+
#include "util/softcordic.h"
1213
#include <ctype.h>
1314
#include <string.h>
1415
#include <math.h>
@@ -1533,17 +1534,25 @@ EvalIntOperator(int op, ExprInt lval, ExprInt rval, int *valid, bool truncMath)
15331534
return popcount(rval);
15341535
case K_QLOG:
15351536
{
1537+
#if 1
1538+
return qlog((uint32_t)rval);
1539+
#else
15361540
double e = log2((double)(uint32_t)rval);
15371541
uint32_t i = (int)e;
15381542
return (i << 27) | (uint32_t)round((e - (double)i)*(1<<27));
1543+
#endif
15391544
}
15401545
break;
15411546
case K_QEXP:
15421547
{
1548+
#if 1
1549+
return qexp((uint32_t)rval);
1550+
#else
15431551
double e = ((uint32_t)rval) >> 27;
15441552
double f = ((double)(rval & ((1<<27)-1))) / (1<<27);
15451553
e = pow(2.0, e+f);
15461554
return (uint32_t)round(e);
1555+
#endif
15471556
}
15481557
break;
15491558
case K_POWER:

‎include/math.h‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -49,6 +49,7 @@
4949
#define exp2f(x) __builtin_exp2f(x)
5050
#define log(x) __builtin_logf(x)
5151
#define logf(x) __builtin_logf(x)
52+
#define log2(x) __builtin_log2f(x)
5253
#define log2f(x) __builtin_log2f(x)
5354
#define log10(x) __builtin_log10f(x)
5455
#define log10f(x) __builtin_log10f(x)

‎util/softcordic.c‎

Lines changed: 165 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,165 @@
1+
// P2 CORDIC Library
2+
// TODO: anything that isn't QLOG or QEXP, possibly also bit-accurate QDIV overflow
3+
// 2025 Ada Gottensträter
4+
//
5+
6+
#include "softcordic.h"
7+
8+
#ifdef __propeller2__
9+
# ifndef FORCE_SOFTCORDIC
10+
# define P2_HARDCORDIC
11+
# endif
12+
#endif
13+
14+
#ifdef P2_HARDCORDIC
15+
uint32_t qexp(uint32_t d) {
16+
uint32_t r;
17+
__asm {
18+
qexp d
19+
getqx r
20+
}
21+
return r;
22+
}
23+
uint32_t qlog(uint32_t d) {
24+
uint32_t r;
25+
__asm {
26+
qlog d
27+
getqx r
28+
}
29+
return r;
30+
}
31+
#else
32+
33+
static const int64_t zdelta_logexp[] = { // hyberbolic cordic deltas (official PNut values)
34+
0x32B803473F,0x179538DEA7,0x0B9A2C912F,0x05C73F7233,0x02E2E683F7,0x01715C285F,0x00B8AB3164,0x005C553C5C,0x002E2A92A3,0x00171547E0,0x000B8AA3C2,0x0005C551DB,0x0002E2A8ED,0x0001715476,0x0000B8AA3B,0x00005C551E,0x00002E2A8F,0x0000171547,0x00000B8AA4,0x000005C552,0x000002E2A9,0x0000017154,0x000000B8AA,0x0000005C55,0x0000002E2B,0x0000001715,0x0000000B8B,0x00000005C5,0x00000002E3,0x0000000171,0x00000000B9
35+
};
36+
37+
#define QLOGEXP_ADJ(N) {x -= x >> N;y -= y >> N;}
38+
#define QEXP_SEC(N) {xd = y>>N, yd = x>>N, zd = zdelta_logexp[N-1]; if (z>=0) {xd=-xd;yd=-yd;zd=-zd;} x-=xd;y-=yd;z+=zd;}
39+
#define QLOG_SEC(N) {xd = y>>N, yd = x>>N, zd = zdelta_logexp[N-1]; if (y<0) {xd=-xd;yd=-yd;zd=-zd;} x-=xd;y-=yd;z+=zd;}
40+
41+
uint32_t qlog(uint32_t d) {
42+
unsigned mag = d ? __builtin_clz(d) : 31;
43+
int64_t y = (int64_t)((d << mag)&0x7FFFFFFF)<<6;
44+
int64_t x = y | 2ull<<37;
45+
int64_t z = 0;
46+
int64_t xd,yd,zd;
47+
48+
QLOG_SEC(1);
49+
QLOGEXP_ADJ(2); // adjust
50+
QLOG_SEC(2);
51+
QLOGEXP_ADJ(3); // adjust
52+
QLOG_SEC(3);
53+
QLOGEXP_ADJ(4); // adjust
54+
QLOG_SEC(4);
55+
QLOG_SEC(4); // double iteration
56+
QLOG_SEC(5);
57+
QLOG_SEC(6);
58+
QLOGEXP_ADJ(7); // adjust
59+
QLOG_SEC(7);
60+
QLOGEXP_ADJ(8); // adjust
61+
QLOG_SEC(8);
62+
QLOG_SEC(9);
63+
QLOGEXP_ADJ(10); // adjust
64+
QLOG_SEC(10);
65+
QLOG_SEC(11);
66+
QLOGEXP_ADJ(12); // adjust
67+
QLOG_SEC(12);
68+
QLOG_SEC(13);
69+
QLOG_SEC(13); // double iteration
70+
QLOGEXP_ADJ(14); // adjust
71+
QLOG_SEC(14);
72+
QLOG_SEC(15);
73+
QLOGEXP_ADJ(16); // adjust
74+
QLOG_SEC(16);
75+
QLOG_SEC(17);
76+
QLOG_SEC(18);
77+
QLOGEXP_ADJ(19); // adjust
78+
QLOG_SEC(19);
79+
QLOGEXP_ADJ(20); // adjust
80+
QLOG_SEC(20);
81+
QLOG_SEC(21);
82+
QLOGEXP_ADJ(22); // adjust
83+
QLOG_SEC(22);
84+
QLOGEXP_ADJ(23); // adjust
85+
QLOG_SEC(23);
86+
QLOGEXP_ADJ(24); // adjust
87+
QLOG_SEC(24);
88+
QLOGEXP_ADJ(25); // adjust
89+
QLOG_SEC(25);
90+
QLOG_SEC(26);
91+
QLOG_SEC(27);
92+
QLOG_SEC(28);
93+
QLOG_SEC(29);
94+
QLOGEXP_ADJ(30); // adjust
95+
QLOG_SEC(30);
96+
QLOG_SEC(31);
97+
98+
int64_t tmp = ((z>>2)&~0x7Fll)+0x80+((int64_t)(mag^31)<<35);
99+
if (mag == 0 && tmp&(1ll<<39) == 0){
100+
tmp = -1;
101+
}
102+
return tmp >> 8;
103+
}
104+
105+
uint32_t qexp(uint32_t d) {
106+
unsigned mag = (d>>27)^31;
107+
int64_t x = 0x7F42E61C5A; // magic number ??
108+
int64_t y = 0;
109+
int64_t z = (int64_t)(d&0x07FFFFFF) << 11;
110+
int64_t xd,yd,zd;
111+
112+
QEXP_SEC(1);
113+
QLOGEXP_ADJ(2); // adjust
114+
QEXP_SEC(2);
115+
QLOGEXP_ADJ(3); // adjust
116+
QEXP_SEC(3);
117+
QLOGEXP_ADJ(4); // adjust
118+
QEXP_SEC(4);
119+
QEXP_SEC(4); // double iteration
120+
QEXP_SEC(5);
121+
QEXP_SEC(6);
122+
QLOGEXP_ADJ(7); // adjust
123+
QEXP_SEC(7);
124+
QLOGEXP_ADJ(8); // adjust
125+
QEXP_SEC(8);
126+
QEXP_SEC(9);
127+
QLOGEXP_ADJ(10); // adjust
128+
QEXP_SEC(10);
129+
QEXP_SEC(11);
130+
QLOGEXP_ADJ(12); // adjust
131+
QEXP_SEC(12);
132+
QEXP_SEC(13);
133+
QEXP_SEC(13); // double iteration
134+
QLOGEXP_ADJ(14); // adjust
135+
QEXP_SEC(14);
136+
QEXP_SEC(15);
137+
QLOGEXP_ADJ(16); // adjust
138+
QEXP_SEC(16);
139+
QEXP_SEC(17);
140+
QEXP_SEC(18);
141+
QLOGEXP_ADJ(19); // adjust
142+
QEXP_SEC(19);
143+
QLOGEXP_ADJ(20); // adjust
144+
QEXP_SEC(20);
145+
QEXP_SEC(21);
146+
QLOGEXP_ADJ(22); // adjust
147+
QEXP_SEC(22);
148+
QLOGEXP_ADJ(23); // adjust
149+
QEXP_SEC(23);
150+
QLOGEXP_ADJ(24); // adjust
151+
QEXP_SEC(24);
152+
QLOGEXP_ADJ(25); // adjust
153+
QEXP_SEC(25);
154+
QEXP_SEC(26);
155+
QEXP_SEC(27);
156+
QEXP_SEC(28);
157+
QEXP_SEC(29);
158+
QLOGEXP_ADJ(30); // adjust
159+
QEXP_SEC(30);
160+
QEXP_SEC(31);
161+
162+
return (((x>>mag)+(y>>mag))+0x40)>>7; // PNut's impl has buggy rounding here (0x20 instead of 0x40)
163+
}
164+
#endif
165+

‎util/softcordic.h‎

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,5 @@
1+
2+
#include <stdint.h>
3+
4+
extern uint32_t qexp(uint32_t d);
5+
extern uint32_t qlog(uint32_t d);

0 commit comments

Comments
 (0)