- static int highest_result = 0;
- // algorithm not applicable to B(0), so just store it
- if (results.empty())
- results.push_back(cln::cl_RA(1));
-
- int n = nn.to_long();
- for (int i=highest_result; i<n/2; ++i) {
- cln::cl_RA B = 0;
- long n = 8;
- long m = 5;
- long d1 = i;
- long d2 = 2*i-1;
- for (int j=i; j>0; --j) {
- B = cln::cl_I(n*m) * (B+results[j]) / (d1*d2);
- n += 4;
- m += 2;
- d1 -= 1;
- d2 -= 2;
+ static unsigned next_r = 0;
+
+ // algorithm not applicable to B(2), so just store it
+ if (!next_r) {
+ results.push_back(cln::recip(cln::cl_RA(6)));
+ next_r = 4;
+ }
+ for (unsigned p=next_r; p<=n; p+=2) {
+ cln::cl_I c = 1;
+ cln::cl_RA b = cln::cl_RA(1-p)/2;
+ const unsigned p3 = p+3;
+ const unsigned p2 = p+2;
+ const unsigned pm = p-2;
+ unsigned i, k;
+ for (i=2, k=0; i <= pm; i += 2, k++) {
+ c = cln::exquo(c * ((p3 - i)*(p2 - i)), (i - 1)*i);
+ b = b + c * results[k];