Calculates the pfs based on stem with in G.  Uses volume Index as guide
179615
function (stemG, wsVI, laVI, SLA) {
    if (stemG < 10) {
      stemG=10;
    }
    var VI = Math.pow( (stemG / wsVI.stems_per_stump / wsVI.constant),(1 / wsVI.power) );

    // Add up for all stems
    var la = laVI.constant * Math.pow(VI,laVI.power) * wsVI.stems_per_stump;
    var wf = 1000 * (la / SLA);  // Foilage Weight in g;
    var pfs = wf/stemG;
    return pfs;
  }
1236915
def money_flow(close_data, high_data, low_data, volume):
    
    catch_errors.check_for_input_len_diff(
        close_data, high_data, low_data, volume
        )
    mf = volume * tp(close_data, high_data, low_data)
    return mf
1160712
def psf_general(vx_size=(1,1,1), fwhm=(5, 5, 6), hradius=8, scale=2):
    '''
    
    '''
    xSig = (scale*fwhm[0]/vx_size[0]) / (2*(2*np.log(2))**.5)
    ySig = (scale*fwhm[1]/vx_size[1]) / (2*(2*np.log(2))**.5)
    zSig = (scale*fwhm[2]/vx_size[2]) / (2*(2*np.log(2))**.5)

    # get the separated kernels
    x = np.arange(-hradius,hradius+1)
    xKrnl = np.exp(-(x**2/(2*xSig**2)))
    yKrnl = np.exp(-(x**2/(2*ySig**2)))
    zKrnl = np.exp(-(x**2/(2*zSig**2)))

    # normalise kernels
    xKrnl /= np.sum(xKrnl)
    yKrnl /= np.sum(yKrnl)
    zKrnl /= np.sum(zKrnl)

    krnl = np.array([zKrnl, yKrnl, xKrnl], dtype=np.float32)

    # for checking if the normalisation worked
    # np.prod( np.sum(krnl,axis=1) )

    # return all kernels together
    return krnl
989122
def G(w, Xs):
    ''''''
    n = len(Xs)
    P = projection_matrix(w)
    Ys = [np.dot(P, X) for X in Xs]
    A = calc_A(Ys)
    A_hat = calc_A_hat(A, skew_matrix(w))

    
    u = sum(np.dot(Y, Y) for Y in Ys) / n
    v = np.dot(A_hat, sum(np.dot(Y, Y) * Y for Y in Ys)) / np.trace(np.dot(A_hat, A))

    return sum((np.dot(Y, Y) - u - 2 * np.dot(Y, v)) ** 2 for Y in Ys)
893184
def stem_frequencies(self):
        """"""
        stemfreq = self.array.sum(axis=1).astype(np.float)
        stemfreq /= stemfreq.sum()
        return stemfreq
773514
private void dfs(Visitor visitor, int v, int[] cc, int id) {
        visitor.visit(v);
        cc[v] = id;
        for (int t = 0; t < n; t++) {
            if (graph[v][t] != 0.0) {
                if (cc[t] == -1) {
                    dfs(visitor, t, cc, id);
                }
            }
        }
    }
741166
function bgrat(a, b, x) {
      var T, fTs, u, logH, logM, lnx2s, lnx2sn, ps, cs, j, ser;

      T = a + 0.5 * (b - 1);
      fTs = 1 / (4 * T * T);
      u = -T * Math.log(x);
      logH = Math.log(b) + log1p(gam1(b)) + b * Math.log(u) - u;
      logM = logH - (algdiv(b, a) + b * Math.log(T));
      cs = [1];   // cs[n] = 1 / (2n + 1)!
      ps = [1];   // ps[n] = p_n in (9.3)
      lnx2s = Math.log(x) * 0.5;   // (ln(x) / 2) ^ 2
      lnx2s *= lnx2s;
      lnx2sn = 1;   // lnx2s ^ n
      j = gratR(b, u, logH);
      ser =
      Math.log(
         series(function(n) {
            var mbmn; // holds m*b - n

            mbmn = -n;
            if (n > 0) {
               j = (b + 2 * n - 2) * (b + 2 * n - 2 + 1) * j +
                     (u + b + 2 * n - 2 + 1) * lnx2sn;
               j *= fTs;
               cs[n] = cs[n - 1] / (2 * n) / (2 * n + 1);
               ps[n] = (b - 1) * cs[n] + 1 / n * series(function(m) {
                  if (m === 0) { return 0; }
                  mbmn += b;
                  return mbmn * cs[m] * ps[n - m];
               }, n);
               lnx2sn *= lnx2s;
            }
            return ps[n] * j;
         })
      );
      return logM + ser;
   }
126598
def _gpinv(probs, kappa, sigma):
    """"""
    x = np.full_like(probs, np.nan)
    if sigma <= 0:
        return x
    ok = (probs > 0) & (probs < 1)
    if np.all(ok):
        if np.abs(kappa) < np.finfo(float).eps:
            x = -np.log1p(-probs)  # pylint: disable=invalid-unary-operand-type
        else:
            x = np.expm1(-kappa * np.log1p(-probs)) / kappa
        x *= sigma
    else:
        if np.abs(kappa) < np.finfo(float).eps:
            x[ok] = -np.log1p(-probs[ok])  # pylint: disable=unsupported-assignment-operation, E1130
        else:
            x[ok] = (  # pylint: disable=unsupported-assignment-operation
                np.expm1(-kappa * np.log1p(-probs[ok])) / kappa
            )
        x *= sigma
        x[probs == 0] = 0
        if kappa >= 0:
            x[probs == 1] = np.inf  # pylint: disable=unsupported-assignment-operation
        else:
            x[probs == 1] = -sigma / kappa  # pylint: disable=unsupported-assignment-operation

    return x
778783
def money_flow_index(close_data, high_data, low_data, volume, period):
    
    catch_errors.check_for_input_len_diff(
        close_data, high_data, low_data, volume
        )
    catch_errors.check_for_period_error(close_data, period)

    mf = money_flow(close_data, high_data, low_data, volume)
    tp = typical_price(close_data, high_data, low_data)

    flow = [tp[idx] > tp[idx-1] for idx in range(1, len(tp))]
    pf = [mf[idx] if flow[idx] else 0 for idx in range(0, len(flow))]
    nf = [mf[idx] if not flow[idx] else 0 for idx in range(0, len(flow))]

    pmf = [sum(pf[idx+1-period:idx+1]) for idx in range(period-1, len(pf))]
    nmf = [sum(nf[idx+1-period:idx+1]) for idx in range(period-1, len(nf))]

    # Dividing by 0 is not an issue, it turns the value into NaN which we would
    # want in that case
    with warnings.catch_warnings():
        warnings.simplefilter("ignore", category=RuntimeWarning)
        money_ratio = np.array(pmf) / np.array(nmf)

    mfi = 100 - (100 / (1 + money_ratio))

    mfi = fill_for_noncomputable_vals(close_data, mfi)

    return mfi
20642
def control_volumes(self):
        
        if self._control_volumes is None:
            #   1/3. * (0.5 * edge_length) * covolume
            # = 1/6 * edge_length**2 * ce_ratio_edge_ratio
            v = self.ei_dot_ei * self.ce_ratios / 6.0
            # Explicitly sum up contributions per cell first. Makes
            # numpy.add.at faster.
            # For every node k (range(4)), check for which edges k appears in
            # local_idx, and sum() up the v's from there.
            idx = self.local_idx
            vals = numpy.array(
                [
                    sum([v[i, j] for i, j in zip(*numpy.where(idx == k)[1:])])
                    for k in range(4)
                ]
            ).T
            #
            self._control_volumes = numpy.zeros(len(self.node_coords))
            numpy.add.at(self._control_volumes, self.cells["nodes"], vals)
        return self._control_volumes
1069640
private int dfs(int[][] merge, int index, int[] order, int i) {
        int n = merge.length + 1;

        if (merge[index][0] > merge.length) {
            i = dfs(merge, merge[index][0] - n, order, i);
        } else {
            order[i++] = merge[index][0];
        }

        if (merge[index][1] > merge.length) {
            i = dfs(merge, merge[index][1] - n, order, i);
        } else {
            order[i++] = merge[index][1];
        }

        return i;
    }
1009588
private void dfs(int v, int[] cc, int id) {
        cc[v] = id;
        for (int t = 0; t < n; t++) {
            if (graph[v][t] != 0.0) {
                if (cc[t] == -1) {
                    dfs(t, cc, id);
                }
            }
        }
    }
929281
def _sigma_pi_hiE(self, Tp, a):
        
        m_p = self._m_p
        csip = (Tp - 3.0) / m_p
        m1 = a[0] * csip ** a[3] * (1 + np.exp(-a[1] * csip ** a[4]))
        m2 = 1 - np.exp(-a[2] * csip ** 0.25)
        multip = m1 * m2
        return self._sigma_inel(Tp) * multip
1328822
def value(self):
        
        # This slightly convoluted calculation is equivalent to df's "Use%";
        # it calculates the percentage of FS usage as a proportion of the
        # space available to *non-root users*. Technically this means it can
        # exceed 100% (when FS is filled to the point that only root can write
        # to it), hence the clamp.
        vfs = os.statvfs(self.filesystem)
        used = vfs.f_blocks - vfs.f_bfree
        total = used + vfs.f_bavail
        return min(1.0, used / total)
1331279
def _cmp_fstruct(self, s1, s2, frac_tol, mask):
        
        if len(s2) > len(s1):
            raise ValueError("s1 must be larger than s2")
        if mask.shape != (len(s2), len(s1)):
            raise ValueError("mask has incorrect shape")

        return is_coord_subset_pbc(s2, s1, frac_tol, mask)
18855
function(T, SA, C, B, n, k) {
        var b, i, j;
        var c0, c1;
        /* compute SAl */
        if (C === B) { getCounts(T, C, n, k); }
        getBuckets(C, B, k, false); /* find starts of buckets */
        j = n - 1;
        b = B[c1 = T[j]];
        j--;
        SA[b++] = (T[j] < c1) ? ~j : j;
        for (i = 0; i < n; i++) {
            if ((j = SA[i]) > 0) {
                ASSERT(T[j] >= T[j+1]);
                if ((c0 = T[j]) !== c1) { B[c1] = b; b = B[c1 = c0]; }
                ASSERT(i < b);
                j--;
                SA[b++] = (T[j] < c1) ? ~j : j;
                SA[i] = 0;
            } else if (j < 0) {
                SA[i] = ~j;
            }
        }
        /* compute SAs */
        if (C === B) { getCounts(T, C, n, k); }
        getBuckets(C, B, k, 1); /* find ends of buckets */
        for (i = n-1, b = B[c1 = 0]; i >= 0; i--) {
            if ((j = SA[i]) > 0) {
                ASSERT(T[j] <= T[j+1]);
                if ((c0 = T[j]) !== c1) { B[c1] = b; b = B[c1 = c0]; }
                ASSERT(b <= i);
                j--;
                SA[--b] = (T[j] > c1) ? ~(j+1) : j;
                SA[i] = 0;
            }
        }
    }
159468
function erf( x ) {
	var sign;
	var ax;
	var z;
	var r;
	var s;
	var y;
	var p;
	var q;

	// Special case: NaN
	if ( x !== x ) {
		return NaN;
	}
	// Special case: +infinity
	if ( x === PINF ) {
		return 1;
	}
	// Special case: -infinity
	if ( x === NINF ) {
		return -1;
	}
	// Special case: +-0
	if ( x === 0 ) {
		return x;
	}
	if ( x < 0 ) {
		sign = true;
		ax = -x;
	} else {
		sign = false;
		ax = x;
	}
	// |x| < 0.84375
	if ( ax < 0.84375 ) {
		if ( ax < SMALL ) {
			if ( ax < VERY_TINY ) {
				// Avoid underflow:
				return 0.125 * (8.0*x + EFX8*x);
			}
			return x + EFX*x;
		}
		z = x * x;
		r = PPC + z*polyvalPP( z );
		s = QQC + z*polyvalQQ( z );
		y = r / s;
		return x + x*y;
	}
	// 0.84375 <= |x| < 1.25
	if ( ax < 1.25 ) {
		s = ax - 1;
		p = PAC + s*polyvalPA( s );
		q = QAC + s*polyvalQA( s );
		if ( sign ) {
			return -ERX - p/q;
		}
		return ERX + p/q;
	}
	// +inf > |x| >= 6
	if ( ax >= 6 ) {
		if ( sign ) {
			return TINY - 1.0; // raise inexact
		}
		return 1.0 - TINY; // raise inexact
	}
	s = 1.0 / (ax*ax);

	// |x| < 1/0.35 ~ 2.857143
	if ( ax < 2.857142857142857 ) {
		r = RAC + s*polyvalRA( s );
		s = SAC + s*polyvalSA( s );
	}
	// |x| >= 1/0.35 ~ 2.857143
	else {
		r = RBC + s*polyvalRB( s );
		s = SBC + s*polyvalSB( s );
	}
	z = setLowWord( ax, 0 ); // pseudo-single (20-bit) precision x
	r = exp( -z*z - 0.5625 ) * exp( (z-ax)*(z+ax) + r/s );
	if ( sign ) {
		return r/ax - 1;
	}
	return 1 - r/ax;
}