Cauchy Integral Experiment

[mike3]

Well I improved the code with this technique:

Code:

/* Do an iteration over a vector v representing values of tetration from -iA to +iA */
/* Odd number of nodes required */
OneIteration(v, A) = {
            local(N=matsize(v)[2]);
            local(newv = v);
            local(expv = exp(v));
            local(logv = log(v));
            local(zi=-I*A);
            local(dz=GetDelta(N,-I*A,I*A));
            local(I1=0);
            local(I2=0);
            local(I3=0);
            local(I4=0);
            local(node1=0);
            local(node2=0);
            local(ctrn=0);
            local(nidx1=0);
            local(nidx2=0);

            /* Update every other node from center to ends */
            ctrn = (N-1)/2+1;
            nidx1 = ctrn;
            nidx2 = ctrn;
            while(nidx1 >= 1,
                  zi = -I*A + (nidx1-1)*dz;
                  I1 = CauchyIntegrateArr(1-I*A, 1+I*A, expv, zi);
                  I2 = UpperSemicircle(L, A, zi);
                  I3 = -CauchyIntegrateArr(-1-I*A, -1+I*A, logv, zi);
                  I4 = LowerSemicircle(conj(L), A, zi);
                  node1 = I1+I2+I3+I4;

                  zi = -I*A + (nidx2-1)*dz;
                  I1 = CauchyIntegrateArr(1-I*A, 1+I*A, expv, zi);
                  I2 = UpperSemicircle(L, A, zi);
                  I3 = -CauchyIntegrateArr(-1-I*A, -1+I*A, logv, zi);
                  I4 = LowerSemicircle(conj(L), A, zi);
                  node2 = I1+I2+I3+I4;

                  newv[nidx1] = node1; /* simultaneously update both */
                  newv[nidx2] = node2;

                  expv[nidx1] = exp(node1);
                  expv[nidx2] = exp(node2);

                  logv[nidx1] = log(node1);
                  logv[nidx2] = log(node2);

                  nidx1 = nidx1 - 2;
                  nidx2 = nidx2 + 2;
               );

            /* Update remaining nodes from ends to center */
            if(ctrn%2 == 0,
               nidx1 = 1;
               nidx2 = N;,
               nidx1 = 2;
               nidx2 = N-1;
              );
            while(nidx1 <= nidx2,
                  zi = -I*A + (nidx1-1)*dz;
                  I1 = CauchyIntegrateArr(1-I*A, 1+I*A, expv, zi);
                  I2 = UpperSemicircle(L, A, zi);
                  I3 = -CauchyIntegrateArr(-1-I*A, -1+I*A, logv, zi);
                  I4 = LowerSemicircle(conj(L), A, zi);
                  node1 = I1+I2+I3+I4;

                  zi = -I*A + (nidx2-1)*dz;
                  I1 = CauchyIntegrateArr(1-I*A, 1+I*A, expv, zi);
                  I2 = UpperSemicircle(L, A, zi);
                  I3 = -CauchyIntegrateArr(-1-I*A, -1+I*A, logv, zi);
                  I4 = LowerSemicircle(conj(L), A, zi);
                  node2 = I1+I2+I3+I4;

                  newv[nidx1] = node1; /* simultaneously update both */
                  newv[nidx2] = node2;

                  expv[nidx1] = exp(node1);
                  expv[nidx2] = exp(node2);

                  logv[nidx1] = log(node1);
                  logv[nidx2] = log(node2);

                  nidx1 = nidx1 + 2;
                  nidx2 = nidx2 - 2;
               );

            return(newv);
}

I update the nodes starting from the center one then every other one out, the I go from the end nodes inward (or the next closest nodes to those end nodes if they were already done in the first pass). It seems to work good for node numbers like 211, 215, 219, etc. but not 209, 213, 217, etc. where it now diverges badly. Why is that?