Annotation of brogue-ce/src/brogue/Math.c, Revision 1.1.1.1
1.1 rubenllo 1: /*
2: * RogueMain.c
3: * Brogue
4: *
5: * Created by Brian Walker on 12/26/08.
6: * Copyright 2012. All rights reserved.
7: *
8: * This file is part of Brogue.
9: *
10: * This program is free software: you can redistribute it and/or modify
11: * it under the terms of the GNU Affero General Public License as
12: * published by the Free Software Foundation, either version 3 of the
13: * License, or (at your option) any later version.
14: *
15: * This program is distributed in the hope that it will be useful,
16: * but WITHOUT ANY WARRANTY; without even the implied warranty of
17: * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18: * GNU Affero General Public License for more details.
19: *
20: * You should have received a copy of the GNU Affero General Public License
21: * along with this program. If not, see <http://www.gnu.org/licenses/>.
22: */
23:
24: #include <time.h>
25: #include <limits.h>
26: #include <stdint.h> // C99
27:
28: #include "Rogue.h"
29: #include "IncludeGlobals.h"
30:
31:
32: // Random number generation
33:
34: short randClump(randomRange theRange) {
35: return randClumpedRange(theRange.lowerBound, theRange.upperBound, theRange.clumpFactor);
36: }
37:
38: // Get a random int between lowerBound and upperBound, inclusive, with probability distribution
39: // affected by clumpFactor.
40: short randClumpedRange(short lowerBound, short upperBound, short clumpFactor) {
41: if (upperBound <= lowerBound) {
42: return lowerBound;
43: }
44: if (clumpFactor <= 1) {
45: return rand_range(lowerBound, upperBound);
46: }
47:
48: short i, total = 0, numSides = (upperBound - lowerBound) / clumpFactor;
49:
50: for(i=0; i < (upperBound - lowerBound) % clumpFactor; i++) {
51: total += rand_range(0, numSides + 1);
52: }
53:
54: for(; i< clumpFactor; i++) {
55: total += rand_range(0, numSides);
56: }
57:
58: return (total + lowerBound);
59: }
60:
61: // Get a random int between lowerBound and upperBound, inclusive
62: boolean rand_percent(short percent) {
63: return (rand_range(0, 99) < clamp(percent, 0, 100));
64: }
65:
66: void shuffleList(short *list, short listLength) {
67: short i, r, buf;
68: for (i=0; i<listLength; i++) {
69: r = rand_range(0, listLength-1);
70: if (i != r) {
71: buf = list[r];
72: list[r] = list[i];
73: list[i] = buf;
74: }
75: }
76: }
77:
78: void fillSequentialList(short *list, short listLength) {
79: short i;
80: for (i=0; i<listLength; i++) {
81: list[i] = i;
82: }
83: }
84:
85: //typedef unsigned long int u4;
86: typedef uint32_t u4;
87: typedef struct ranctx { u4 a; u4 b; u4 c; u4 d; } ranctx;
88:
89: static ranctx RNGState[2];
90:
91: #define rot(x,k) (((x)<<(k))|((x)>>(32-(k))))
92: u4 ranval( ranctx *x ) {
93: u4 e = x->a - rot(x->b, 27);
94: x->a = x->b ^ rot(x->c, 17);
95: x->b = x->c + x->d;
96: x->c = x->d + e;
97: x->d = e + x->a;
98: return x->d;
99: }
100:
101: void raninit( ranctx *x, u4 seed ) {
102: u4 i;
103: x->a = 0xf1ea5eed, x->b = x->c = x->d = seed;
104: for (i=0; i<20; ++i) {
105: (void)ranval(x);
106: }
107: }
108:
109: /* ----------------------------------------------------------------------
110: range
111:
112: returns a number between 0 and N-1
113: without any bias.
114:
115: */
116:
117: #define RAND_MAX_COMBO ((unsigned long) UINT32_MAX)
118:
119: long range(long n, short RNG) {
120: unsigned long div;
121: long r;
122:
123: div = RAND_MAX_COMBO/n;
124:
125: do {
126: r = ranval(&(RNGState[RNG])) / div;
127: } while (r >= n);
128:
129: return r;
130: }
131:
132: // Get a random int between lowerBound and upperBound, inclusive, with uniform probability distribution
133:
134: #ifdef AUDIT_RNG // debug version
135: long rand_range(long lowerBound, long upperBound) {
136: int retval;
137: char RNGMessage[100];
138: if (upperBound <= lowerBound) {
139: return lowerBound;
140: }
141: long interval = upperBound - lowerBound + 1;
142: brogueAssert(interval > 1); // to verify that we didn't wrap around
143: retval = lowerBound + range(interval, rogue.RNG);
144: if (rogue.RNG == RNG_SUBSTANTIVE) {
145: randomNumbersGenerated++;
146: if (1) { //randomNumbersGenerated >= 1128397) {
147: sprintf(RNGMessage, "\n#%lu, %ld to %ld: %ld", randomNumbersGenerated, lowerBound, upperBound, retval);
148: RNGLog(RNGMessage);
149: }
150: }
151: return retval;
152: }
153: #else // normal version
154: long rand_range(long lowerBound, long upperBound) {
155: if (upperBound <= lowerBound) {
156: return lowerBound;
157: }
158: if (rogue.RNG == RNG_SUBSTANTIVE) {
159: randomNumbersGenerated++;
160: }
161: long interval = upperBound - lowerBound + 1;
162: brogueAssert(interval > 1); // to verify that we didn't wrap around
163: return lowerBound + range(interval, rogue.RNG);
164: }
165: #endif
166:
167: // seeds with the time if called with a parameter of 0; returns the seed regardless.
168: // All RNGs are seeded simultaneously and identically.
169: unsigned long seedRandomGenerator(unsigned long seed) {
170: if (seed == 0) {
171: seed = (unsigned long) time(NULL) - 1352700000;
172: }
173: raninit(&(RNGState[RNG_SUBSTANTIVE]), seed);
174: raninit(&(RNGState[RNG_COSMETIC]), seed);
175: return seed;
176: }
177:
178:
179: // Fixed-point arithmetic
180:
181: fixpt fp_round(fixpt x) {
182: long long div = x / FP_FACTOR, rem = x % FP_FACTOR;
183: int sign = (x >= 0) - (x < 0);
184:
185: if (rem >= FP_FACTOR / 2 || rem <= -FP_FACTOR / 2) {
186: return div + sign;
187: } else {
188: return div;
189: }
190: }
191:
192: // Returns the bit position of the most significant bit of x, where the unit
193: // bit has position 1. Returns 0 if x=0.
194: static int msbpos(unsigned long long x) {
195: if (x == 0) return 0;
196: int n = 0;
197: do {
198: n += 1;
199: } while (x >>= 1);
200: return n;
201: }
202:
203: static fixpt fp_exp2(int n) {
204: return (n >= 0 ? FP_FACTOR << n : FP_FACTOR >> -n);
205: }
206:
207: // Calculates sqrt(u) using the bisection method to find the root of
208: // f(x) = x^2 - u.
209: fixpt fp_sqrt(fixpt u) {
210:
211: static const fixpt SQUARE_ROOTS[] = { // values were computed by the code that follows
212: 0, 65536, 92682, 113511, 131073, 146543, 160529, 173392, 185363, 196608, 207243, 217359, 227023, 236293, 245213, 253819,
213: 262145, 270211, 278045, 285665, 293086, 300323, 307391, 314299, 321059, 327680, 334169, 340535, 346784, 352923, 358955, 364889,
214: 370727, 376475, 382137, 387717, 393216, 398640, 403991, 409273, 414487, 419635, 424721, 429749, 434717, 439629, 444487, 449293,
215: 454047, 458752, 463409, 468021, 472587, 477109, 481589, 486028, 490427, 494786, 499107, 503391, 507639, 511853, 516031, 520175,
216: 524289, 528369, 532417, 536435, 540423, 544383, 548313, 552217, 556091, 559939, 563762, 567559, 571329, 575077, 578799, 582497,
217: 586171, 589824, 593453, 597061, 600647, 604213, 607755, 611279, 614783, 618265, 621729, 625173, 628599, 632007, 635395, 638765,
218: 642119, 645455, 648773, 652075, 655360, 658629, 661881, 665117, 668339, 671545, 674735, 677909, 681071, 684215, 687347, 690465,
219: 693567, 696657, 699733, 702795, 705845, 708881, 711903, 714913, 717911, 720896, 723869, 726829, 729779, 732715, 735639, 738553
220: };
221:
222: if (u < 0) return -fp_sqrt(-u);
223:
224: if ((u & (127LL << FP_BASE)) == u) {
225: // u is an integer between 0 and 127
226: return SQUARE_ROOTS[u >> FP_BASE];
227: }
228:
229: // Find the unique k such that 2^(k-1) <= u < 2^k
230: // FP_BASE is the msbpos-1 of FP_FACTOR ("one")
231: int k = msbpos(u) - FP_BASE;
232:
233: fixpt x = 0, fx, upper, lower;
234: // Since 2^(k-1) <= u < 2^k, we have 2^(ceil(k/2)-1) <= sqrt(u) < 2^ceil(k/2).
235: // First ineq. from sqrt(u) >= 2^[(k-1)/2] = 2^[k/2 + 1/2 - 1] >= 2^(ceil(k/2) - 1)
236: // To calculate ceil(k/2), do k/2 but add 1 to k if positive.
237: upper = fp_exp2((k + (k > 0))/2);
238: lower = upper / 2;
239:
240: while (upper != lower + 1) {
241: x = (upper + lower) / 2;
242: fx = FP_MUL(x, x) - u;
243:
244: if (fx == 0) {
245: break;
246: } else if (fx > 0) {
247: upper = x;
248: } else {
249: lower = x;
250: }
251: }
252:
253: return x;
254: }
255:
256: // Returns base to the power of expn
257: fixpt fp_pow(fixpt base, int expn) {
258: if (base == 0) return 0;
259:
260: if (expn < 0) {
261: base = FP_DIV(FP_FACTOR, base);
262: expn = -expn;
263: }
264:
265: fixpt res = FP_FACTOR, err = 0;
266: while (expn--) {
267: res = res * base + (err * base) / FP_FACTOR;
268: err = res % FP_FACTOR;
269: res /= FP_FACTOR;
270: }
271:
272: return res + fp_round(err);
273: }
CVSweb
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to Math.c CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [contributed] / brogue-ce / src / brogue