- // do acceleration transformation (hoelder convolution [BBB])
- if (need_hoelder) {
-
- ex result;
- const int size = x.nops();
- lst newx;
- for (lst::const_iterator it = x.begin(); it != x.end(); ++it) {
- newx.append(*it / y);
- }
-
- for (int r=0; r<=size; ++r) {
- ex buffer = pow(-1, r);
- ex p = 2;
- bool adjustp;
- do {
- adjustp = false;
- for (lst::const_iterator it = newx.begin(); it != newx.end(); ++it) {
- if (*it == 1/p) {
- p += (3-p)/2;
- adjustp = true;
- continue;
- }
- }
- } while (adjustp);
- ex q = p / (p-1);
- lst qlstx;
- lst qlsts;
- for (int j=r; j>=1; --j) {
- qlstx.append(1-newx.op(j-1));
- if (newx.op(j-1).info(info_flags::real) && newx.op(j-1) > 1 && newx.op(j-1) <= 2) {
- qlsts.append( s.op(j-1));
- } else {
- qlsts.append( -s.op(j-1));
- }
- }
- if (qlstx.nops() > 0) {
- buffer *= G_numeric(qlstx, qlsts, 1/q);
- }
- lst plstx;
- lst plsts;
- for (int j=r+1; j<=size; ++j) {
- plstx.append(newx.op(j-1));
- plsts.append(s.op(j-1));
- }
- if (plstx.nops() > 0) {
- buffer *= G_numeric(plstx, plsts, 1/p);
- }
- result += buffer;
+ // do transformation
+ Gparameter pendint;
+ ex result = G_transform(pendint, a, scale, gsyms);
+ // replace dummy symbols with their values
+ result = result.eval().expand();
+ result = result.subs(subslst).evalf();
+ if (!is_a<numeric>(result))
+ throw std::logic_error("G_do_trafo: G_transform returned non-numeric result");
+
+ cln::cl_N ret = ex_to<numeric>(result).to_cl_N();
+ return ret;
+}
+
+// handles the transformations and the numerical evaluation of G
+// the parameter x, s and y must only contain numerics
+static cln::cl_N
+G_numeric(const std::vector<cln::cl_N>& x, const std::vector<int>& s,
+ const cln::cl_N& y)
+{
+ // check for convergence and necessary accelerations
+ bool need_trafo = false;
+ bool need_hoelder = false;
+ std::size_t depth = 0;
+ for (std::size_t i = 0; i < x.size(); ++i) {
+ if (!zerop(x[i])) {
+ ++depth;
+ const cln::cl_N x_y = abs(x[i]) - y;
+ if (instanceof(x_y, cln::cl_R_ring) &&
+ realpart(x_y) < cln::least_negative_float(cln::float_format(Digits - 2)))
+ need_trafo = true;
+
+ if (abs(abs(x[i]/y) - 1) < 0.01)
+ need_hoelder = true;