Successors function for solving the cross.
435055
def cross_successors(state, last_action=None):
        
        centres, edges = state
        acts = sum([
            [s, s.inverse(), s * 2] for s in
            map(Step, "RUFDRB".replace(last_action.face if last_action else "", "", 1))
            ], [])
        for step in acts:
            yield step, (centres, CrossSolver._rotate(edges, step))
1288571
def cross(v1, v2):
    '''
    '''
    return (
        v1[1] * v2[2] - v1[2] * v2[1],
        v1[2] * v2[0] - v1[0] * v2[2],
        v1[0] * v2[1] - v1[1] * v2[0],
    )
368902
@Override
    public void solve(ZMatrixRMaj B, ZMatrixRMaj X) {
        if( X.numRows != numCols )
            throw new IllegalArgumentException("Unexpected dimensions for X: X rows = "+X.numRows+" expected = "+numCols);
        else if( B.numRows != numRows || B.numCols != X.numCols )
            throw new IllegalArgumentException("Unexpected dimensions for B");

        U = decomposer.getR(U,true);
        final double gammas[] = decomposer.getGammas();
        final double dataQR[] = QR.data;

        final int BnumCols = B.numCols;

        // solve each column one by one
        for( int colB = 0; colB < BnumCols; colB++ ) {

            // make a copy of this column in the vector
            for( int i = 0; i < numRows; i++ ) {
                int indexB = (i*BnumCols + colB)*2;
                a[i*2]   = B.data[indexB];
                a[i*2+1] = B.data[indexB+1];
            }

            // Solve Qa=b
            // a = Q'b
            // a = Q_{n-1}...Q_2*Q_1*b
            //
            // Q_n*b = (I-gamma*u*u^H)*b = b - u*(gamma*U^H*b)
            for( int n = 0; n < numCols; n++ ) {
                int indexU = (n*numRows + n + 1)*2;

                double realUb = a[n*2];
                double imagUb = a[n*2+1];

                // U^H*b
                for( int i = n+1; i < numRows; i++ ) {
                    double realU = dataQR[indexU++];
                    double imagU = -dataQR[indexU++];

                    double realB = a[i*2];
                    double imagB = a[i*2+1];

                    realUb += realU*realB - imagU*imagB;
                    imagUb += realU*imagB + imagU*realB;
                }

                // gamma*U^T*b
                realUb *= gammas[n];
                imagUb *= gammas[n];

                a[n*2]   -= realUb;
                a[n*2+1] -= imagUb;

                indexU = (n*numRows + n + 1)*2;
                for( int i = n+1; i < numRows; i++) {
                    double realU = dataQR[indexU++];
                    double imagU = dataQR[indexU++];

                    a[i*2]   -= realU*realUb - imagU*imagUb;
                    a[i*2+1] -= realU*imagUb + imagU*realUb;
                }
            }

            // solve for Rx = b using the standard upper triangular solver
            TriangularSolver_ZDRM.solveU(U.data, a, numCols);

            // save the results

            for( int i = 0; i < numCols; i++ ) {
                int indexX = (i*X.numCols+colB)*2;

                X.data[indexX]   = a[i*2];
                X.data[indexX+1] = a[i*2+1];
            }
        }
    }
346688
def solve(m,c):
    
    G={ (m,c,1):[] }
    frontier=[ (m,c,1) ]  # 1 as boat starts on left bank 
    while len(frontier) > 0:
        hold=list(frontier)
        for node in hold:
            newnode=[]
            frontier.remove(node)
            newnode.extend(pick_next_boat_trip(node, m,c, frontier))
            for neighbor in newnode:
                if neighbor not in G:
                    G[node].append(neighbor)
                    G[neighbor]=[node]
                    frontier.append(neighbor)
    return mod_plan.find_path_BFS(G,(m,c,1),(0,0,0))
1172491
function bnModInverse(m)
{
  var ac = m.isEven();
  if ((this.isEven() && ac) || m.signum() == 0) return Int128.ZERO;
  var u = m.clone(),
    v = this.clone();
  var a = nbv(1),
    b = nbv(0),
    c = nbv(0),
    d = nbv(1);
  while (u.signum() != 0)
  {
    while (u.isEven())
    {
      u.rShiftTo(1, u);
      if (ac)
      {
        if (!a.isEven() || !b.isEven())
        {
          a.addTo(this, a);
          b.subTo(m, b);
        }
        a.rShiftTo(1, a);
      }
      else if (!b.isEven()) b.subTo(m, b);
      b.rShiftTo(1, b);
    }
    while (v.isEven())
    {
      v.rShiftTo(1, v);
      if (ac)
      {
        if (!c.isEven() || !d.isEven())
        {
          c.addTo(this, c);
          d.subTo(m, d);
        }
        c.rShiftTo(1, c);
      }
      else if (!d.isEven()) d.subTo(m, d);
      d.rShiftTo(1, d);
    }
    if (u.compareTo(v) >= 0)
    {
      u.subTo(v, u);
      if (ac) a.subTo(c, a);
      b.subTo(d, b);
    }
    else
    {
      v.subTo(u, v);
      if (ac) c.subTo(a, c);
      d.subTo(b, d);
    }
  }
  if (v.compareTo(Int128.ONE) != 0) return Int128.ZERO;
  if (d.compareTo(m) >= 0) return d.subtract(m);
  if (d.signum() < 0) d.addTo(m, d);
  else return d;
  if (d.signum() < 0) return d.add(m);
  else return d;
}
509895
def tqxn(mt, x, n, t):
    
    return tpx(mt, x, t) * qx(mt, x + n)
699930
def _multicomplex2(f, fx, x, h):
        """"""
        n = len(x)
        ee = np.diag(h)
        hess = np.outer(h, h)
        cmplx_wrap = Bicomplex.__array_wrap__
        for i in range(n):
            for j in range(i, n):
                zph = Bicomplex(x + 1j * ee[i, :], ee[j, :])
                hess[i, j] = cmplx_wrap(f(zph)).imag12 / hess[j, i]
                hess[j, i] = hess[i, j]
        return hess
659027
private static Triangle findnext2(Vector3 p, Triangle v) {
		if (v.abnext != null && !v.abnext.halfplane)
			return v.abnext;
		if (v.bcnext != null && !v.bcnext.halfplane)
			return v.bcnext;
		if (v.canext != null && !v.canext.halfplane)
			return v.canext;
		return null;
	}
189750
def dot(v1, v2):
    '''
    '''
    x1, y1, z1 = v1
    x2, y2, z2 = v2
    return x1 * x2 + y1 * y2 + z1 * z2
212602
def xstep(self):
        
        """

        self.YU[:] = self.Y - self.U

        b = self.DSf + self.rho*sl.rfftn(self.YU, None, self.cri.axisN)
        if self.cri.Cd == 1:
            self.Xf[:] = sl.solvedbi_sm(self.Df, self.mu + self.rho,
                                        b, self.c, self.cri.axisM)
        else:
            self.Xf[:] = sl.solvemdbi_ism(self.Df, self.mu + self.rho, b,
                                          self.cri.axisM, self.cri.axisC)

        self.X = sl.irfftn(self.Xf, self.cri.Nv, self.cri.axisN)

        if self.opt['LinSolveCheck']:
            Dop = lambda x: sl.inner(self.Df, x, axis=self.cri.axisM)
            if self.cri.Cd == 1:
                DHop = lambda x: np.conj(self.Df) * x
            else:
                DHop = lambda x: sl.inner(np.conj(self.Df), x,
                                          axis=self.cri.axisC)
            ax = DHop(Dop(self.Xf)) + (self.mu + self.rho)*self.Xf
            self.xrrs = sl.rrs(ax, b)
        else:
            self.xrrs = None
760169
def factorize_and_solve(self, A, b):
        

        self.factorize(A)

        return self.solve(b)
929977
def cross(self,p_i,p_j, max_depth = 2):
        """"""
        # only choose crossover points for out_types available in both programs
        # pdb.set_trace()
        # determine possible outttypes
        types_p_i = [t for t in [p.out_type for p in p_i]]
        types_p_j = [t for t in [p.out_type for p in p_j]]
        types = set(types_p_i).intersection(types_p_j)

        # grab subtree of p_i
        p_i_sub = [i for i,n in enumerate(p_i) if n.out_type in types]
        x_i_end = self.random_state.choice(p_i_sub)
        x_i_begin = x_i_end
        arity_sum = p_i[x_i_end].arity[p_i[x_i_end].in_type]
        # print("x_i_end:",x_i_end)
        # i = 0
        while (arity_sum > 0): #and i < 1000:
            if x_i_begin == 0:
                print("arity_sum:",arity_sum,"x_i_begin:",x_i_begin,"x_i_end:",x_i_end)
            x_i_begin -= 1
            arity_sum += p_i[x_i_begin].arity[p_i[x_i_begin].in_type]-1
        #     i += 1
        # if i == 1000:
        #     print("in variation")
        #     pdb.set_trace()

        # grab subtree of p_j with matching out_type to p_i[x_i_end]
        p_j_sub = [i for i,n in enumerate(p_j) if n.out_type == p_i[x_i_end].out_type]
        x_j_end = self.random_state.choice(p_j_sub)
        x_j_begin = x_j_end
        arity_sum = p_j[x_j_end].arity[p_j[x_j_end].in_type]
        # i = 0
        while (arity_sum > 0): #and i < 1000:
            if x_j_begin == 0:
                print("arity_sum:",arity_sum,"x_j_begin:",x_j_begin,"x_j_end:",x_j_end)
                print("p_j:",p_j)
            x_j_begin -= 1
            arity_sum += p_j[x_j_begin].arity[p_j[x_j_begin].in_type]-1
        #     i += 1
        # if i == 1000:
        #     print("in variation")
        #     pdb.set_trace()
        #swap subtrees
        tmpi = p_i[:]
        tmpj = p_j[:]
        tmpi[x_i_begin:x_i_end+1:],tmpj[x_j_begin:x_j_end+1:] = \
                        tmpj[x_j_begin:x_j_end+1:],tmpi[x_i_begin:x_i_end+1:]

        if not self.is_valid_program(p_i) or not self.is_valid_program(p_j):
            # pdb.set_trace()
            print("parent 1:",p_i,"x_i_begin:",x_i_begin,"x_i_end:",x_i_end)
            print("parent 2:",p_j,"x_j_begin:",x_j_begin,"x_j_end:",x_j_end)
            print("child 1:",tmpi)
            print("child 2:",tmpj)
            raise ValueError('Crossover produced an invalid program.')

        # size check, then assignment
        if len(tmpi) <= 2**max_depth-1:
            p_i[:] = tmpi
        if len(tmpj) <= 2**max_depth-1:
            p_j[:] = tmpj
1361820
function cross(A, B) {
    if (A.shape.length !== 1 || B.shape.length !== 1) {
        throw new Error('A or B is not 1D');
    }
    if (A.length < 3 || B.length < 3) {
        throw new Error('A or B is less than 3 in length');
    }
    var a1 = A.getN(0);
    var a2 = A.getN(1);
    var a3 = A.getN(2);
    var b1 = B.getN(0);
    var b2 = B.getN(1);
    var b3 = B.getN(2);
    return new ndarray_1.NDArray([
        a2 * b3 - a3 * b2,
        a3 * b1 - a1 * b3,
        a1 * b2 - a2 * b1
    ]);
}
559212
function solveIntersection (a, b, c, d) {
  var ba = rvSub(b, a)
  var dc = rvSub(d, c)

  var baXdc = ratPerp(ba, dc)

  if (ratSign(baXdc) === 0) {
    return null
  }

  var ac = rvSub(a, c)
  var dcXac = ratPerp(dc, ac)

  var t = ratDiv(dcXac, baXdc)
  var s = rvMuls(ba, t)
  var r = rvAdd(a, s)

  return r
}
700286
def cross(b,c)
        unless b.size == 3 and c.size == 3
          raise "Vectors must be of length 3"
        end
        Vector[b[1]*c[2] - b[2]*c[1], b[2]*c[0] - b[0]*c[2], b[0]*c[1] - b[1]*c[0]]
      end
164521
def step(self, y, u, t, h):
		
		k1 = h * self.func(t, y, u)
		k2 = h * self.func(t + .5*h, y + .5*h*k1, u)
		k3 = h * self.func(t + .5*h, y + .5*h*k2, u)
		k4 = h * self.func(t + h, y + h*k3, u)
		return y + (k1 + 2*k2 + 2*k3 + k4) / 6.0
117438
private void ssMergeBackward(int PA, int first, int middle, int last, int buf,
                                 int depth) {
        // PA, first, middle, last, buf are pointers in SA
        int p1, p2;// pointers in SA
        int a, b, c, bufend;// pointers in SA
        int t, r, x;

        bufend = buf + (last - middle) - 1;
        ssBlockSwap(buf, middle, last - middle);

        x = 0;
        if (SA[bufend] < 0) {
            p1 = PA + ~SA[bufend];
            x |= 1;
        } else {
            p1 = PA + SA[bufend];
        }
        if (SA[middle - 1] < 0) {
            p2 = PA + ~SA[middle - 1];
            x |= 2;
        } else {
            p2 = PA + SA[middle - 1];
        }
        for (t = SA[a = last - 1], b = bufend, c = middle - 1; ; ) {
            r = ssCompare(p1, p2, depth);
            if (0 < r) {
                if ((x & 1) != 0) {
                    do {
                        SA[a--] = SA[b];
                        SA[b--] = SA[a];
                    }
                    while (SA[b] < 0);
                    x ^= 1;
                }
                SA[a--] = SA[b];
                if (b <= buf) {
                    SA[buf] = t;
                    break;
                }
                SA[b--] = SA[a];
                if (SA[b] < 0) {
                    p1 = PA + ~SA[b];
                    x |= 1;
                } else {
                    p1 = PA + SA[b];
                }
            } else if (r < 0) {
                if ((x & 2) != 0) {
                    do {
                        SA[a--] = SA[c];
                        SA[c--] = SA[a];
                    }
                    while (SA[c] < 0);
                    x ^= 2;
                }
                SA[a--] = SA[c];
                SA[c--] = SA[a];
                if (c < first) {
                    while (buf < b) {
                        SA[a--] = SA[b];
                        SA[b--] = SA[a];
                    }
                    SA[a] = SA[b];
                    SA[b] = t;
                    break;
                }
                if (SA[c] < 0) {
                    p2 = PA + ~SA[c];
                    x |= 2;
                } else {
                    p2 = PA + SA[c];
                }
            } else {
                if ((x & 1) != 0) {
                    do {
                        SA[a--] = SA[b];
                        SA[b--] = SA[a];
                    }
                    while (SA[b] < 0);
                    x ^= 1;
                }
                SA[a--] = ~SA[b];
                if (b <= buf) {
                    SA[buf] = t;
                    break;
                }
                SA[b--] = SA[a];
                if ((x & 2) != 0) {
                    do {
                        SA[a--] = SA[c];
                        SA[c--] = SA[a];
                    }
                    while (SA[c] < 0);
                    x ^= 2;
                }
                SA[a--] = SA[c];
                SA[c--] = SA[a];
                if (c < first) {
                    while (buf < b) {
                        SA[a--] = SA[b];
                        SA[b--] = SA[a];
                    }
                    SA[a] = SA[b];
                    SA[b] = t;
                    break;
                }
                if (SA[b] < 0) {
                    p1 = PA + ~SA[b];
                    x |= 1;
                } else {
                    p1 = PA + SA[b];
                }
                if (SA[c] < 0) {
                    p2 = PA + ~SA[c];
                    x |= 2;
                } else {
                    p2 = PA + SA[c];
                }
            }
        }
    }