ca holds coefficients from one line of table 36.A, p. 250
187536
function Ca (A, B, M0, M1) {
  this.A = A
  this.B = B
  this.M0 = M0
  this.M1 = M1
}
1207540
def _aicc(model_results, nobs):
    """"""
    aic = model_results.aic
    df_model = model_results.df_model + 1  # add one for constant term
    return aic + 2. * df_model * (nobs / (nobs - df_model - 1.) - 1.)
389159
private function computeCoefficients(): void
    {
        $samplesMatrix = $this->getSamplesMatrix();
        $targetsMatrix = $this->getTargetsMatrix();

        $ts = $samplesMatrix->transpose()->multiply($samplesMatrix)->inverse();
        $tf = $samplesMatrix->transpose()->multiply($targetsMatrix);

        $this->coefficients = $ts->multiply($tf)->getColumnValues(0);
        $this->intercept = array_shift($this->coefficients);
    }
983232
def csi(self):
        """"""
        return self.table[0, 0] / (self.table[0, 0] + self.table[0, 1] + self.table[1, 0])
1087679
private static int calcOcc(Alternative ab, boolean calcMax, AgBuffer cb) {
        int occ;
        int i;
        int max;
        State state;

        occ = 0;
        max = ab.getArgCount();
        for (i = 0; i < max; i++) {
            state = cb.lookup(ab.getArgAttribute(i));
            if (state != null) {
                occ += (calcMax)? state.maxOcc : state.minOcc;
            } else {
                occ += 1;
            }
        }
        return occ;
    }
360338
def energy_coefficients(m1, m2, s1z=0, s2z=0, phase_order=-1, spin_order=-1):
    
    implemented_phase_order = 7
    implemented_spin_order = 7
    if phase_order > implemented_phase_order:
        raise ValueError("pN coeffiecients of that order have not been implemented")
    elif phase_order == -1:
        phase_order = implemented_phase_order

    if spin_order > implemented_spin_order:
        raise ValueError("pN coeffiecients of that order have not been implemented")
    elif spin_order == -1:
        spin_order = implemented_spin_order

    qmdef1 = 1.0
    qmdef2 = 1.0

    M = m1 + m2
    dm = (m1-m2)/M
    m1M = m1 / M
    m2M = m2 / M

    s1z = s1z * m1M * m1M
    s2z = s2z * m2M * m2M

    _, eta = mass1_mass2_to_mchirp_eta(m1, m2)

    ecof = numpy.zeros(phase_order+1)
    # Orbital terms
    if phase_order >= 0:
        ecof[0] = 1.0
    if phase_order >= 1:
        ecof[1] = 0
    if phase_order >= 2:
        ecof[2] = -(1.0/12.0) * (9.0 + eta)
    if phase_order >= 3:
        ecof[3] = 0
    if phase_order >= 4:
        ecof[4] = (-81.0 + 57.0*eta - eta*eta) / 24.0
    if phase_order >= 5:
        ecof[5] = 0
    if phase_order >= 6:
        ecof[6] = - 675.0/64.0 + ( 34445.0/576.0    \
              - 205.0/96.0 * lal.PI * lal.PI ) * eta  \
              - (155.0/96.0) *eta * eta - 35.0/5184.0 * eta * eta
    # Spin terms

    ESO15s1 = 8.0/3.0 + 2.0*m2/m1
    ESO15s2 = 8.0/3.0 + 2.0*m1/m2

    ESS2 = 1.0 / eta
    EQM2s1 = qmdef1/2.0/m1M/m1M
    EQM2s1L = -qmdef1*3.0/2.0/m1M/m1M
    #EQM2s2 = qmdef2/2.0/m2M/m2M
    EQM2s2L = -qmdef2*3.0/2.0/m2M/m2M

    ESO25s1 = 11.0 - 61.0*eta/9.0 + (dm/m1M) * (-3.0 + 10.*eta/3.0)
    ESO25s2 = 11.0 - 61.0*eta/9.0 + (dm/m2M) * (3.0 - 10.*eta/3.0)

    ESO35s1 = 135.0/4.0 - 367.0*eta/4.0 + 29.0*eta*eta/12.0 + (dm/m1M) * (-27.0/4.0 + 39.0*eta - 5.0*eta*eta/4.0)
    ESO35s2 = 135.0/4.0 - 367.0*eta/4.0 + 29.0*eta*eta/12.0 - (dm/m2M) * (-27.0/4.0 + 39.0*eta - 5.0*eta*eta/4.0)

    if spin_order >=3:
        ecof[3] += ESO15s1 * s1z + ESO15s2 * s2z
    if spin_order >=4:
        ecof[4] += ESS2 * (s1z*s2z - 3.0*s1z*s2z)
        ecof[4] += EQM2s1*s1z*s1z + EQM2s1*s2z*s2z + EQM2s1L*s1z*s1z + EQM2s2L*s2z*s2z
    if spin_order >=5:
        ecof[5] = ESO25s1*s1z + ESO25s2*s2z
    if spin_order >=7:
        ecof[7] += ESO35s1*s1z + ESO35s2*s2z

    return ecof
391824
def _erfc(a):
    
    # approximation for abs(a) < 8 and abs(a) >= 1
    P = [
        2.46196981473530512524E-10,
        5.64189564831068821977E-1,
        7.46321056442269912687E0,
        4.86371970985681366614E1,
        1.96520832956077098242E2,
        5.26445194995477358631E2,
        9.34528527171957607540E2,
        1.02755188689515710272E3,
        5.57535335369399327526E2,
    ]

    Q = [
        1.32281951154744992508E1,
        8.67072140885989742329E1,
        3.54937778887819891062E2,
        9.75708501743205489753E2,
        1.82390916687909736289E3,
        2.24633760818710981792E3,
        1.65666309194161350182E3,
        5.57535340817727675546E2,
    ]

    # approximation for abs(a) >= 8
    R = [
        5.64189583547755073984E-1,
        1.27536670759978104416E0,
        5.01905042251180477414E0,
        6.16021097993053585195E0,
        7.40974269950448939160E0,
        2.97886665372100240670E0,
    ]

    S = [
        2.26052863220117276590E0,
        9.39603524938001434673E0,
        1.20489539808096656605E1,
        1.70814450747565897222E1,
        9.60896809063285878198E0,
        3.36907645100081516050E0,
    ]

    # Shortcut special cases
    if a == 0:
        return 1
    if a >= MAXVAL:
        return 0
    if a <= -MAXVAL:
        return 2

    x = a
    if a < 0:
        x = -a

    # computationally cheaper to calculate erf for small values, I guess.
    if x < 1:
        return 1 - erf(a)

    z = -a * a

    z = math.exp(z)

    if x < 8:
        p = _polevl(x, P, 8)
        q = _p1evl(x, Q, 8)
    else:
        p = _polevl(x, R, 5)
        q = _p1evl(x, S, 6)

    y = (z * p) / q

    if a < 0:
        y = 2 - y

    return y
14935
def CA(self):
        '''
        

        '''
        try:
            return self._CA
        except AttributeError:
            pass
        self._CA = [self.C, self.A]
        return self._CA
774343
def is_Linear(self):
        """"""
        coeff_dict = self.expression.as_coefficients_dict()
        for key in coeff_dict.keys():
            if len(key.free_symbols) < 2 and (key.is_Add or key.is_Mul or key.is_Atom):
                pass
            else:
                return False
            if key.is_Pow and key.args[1] != 1:
                return False
        else:
            return True
581745
public void getEachCellContent(int YNum, int cc, TableCandidate tc) {
    	int footnoteLineIndex = tc.getFootnoteBeginRow();
    	String[] contentsInRows = new String[footnoteLineIndex];
    	float[] leftX_tableColumns = tc.getLeftX_tableColumns();
    	float[] rightX_tableColumns = tc.getRightX_tableColumns();
    	String[][] cells = new String[YNum][cc];
		boolean[][] crossCells = new boolean[YNum][cc];
		for (int tt=0; tt<tc.getRows().size(); tt++) {
			for (int qq=0; qq<cc; qq++) {
				cells[tt][qq] = "";
				crossCells[tt][qq]=false; 
			}
		}
    	for (int tt=0; tt<footnoteLineIndex; tt++) {
			String cellsThisLine = "";
			TableRow row = tc.getRows().get(tt);
			int cellThisLine = row.getCells().size(); //  pmoj.objectNum_eachY[tt];
			for (int zz=0; zz<cellThisLine; zz++) {
				for (int qq=0; qq<cc; qq++) {
					if ( (row.getCells().get(zz).getX()>=leftX_tableColumns[qq]) 
					&& (row.getCells().get(zz).getEndX()<=rightX_tableColumns[qq]) ) 
							cells[tt][qq] =cells[tt][qq] + row.getCells().get(zz).getText() + " ";
					if ((row.getCells().get(zz).getX()<leftX_tableColumns[qq]) 
					&&(leftX_tableColumns[qq]<=row.getCells().get(zz).getEndX()) )
							cells[tt][qq] = cells[tt][qq] +  row.getCells().get(zz).getText() + " ";
					if ((leftX_tableColumns[qq]>=row.getCells().get(zz).getX()) 
					&& (row.getCells().get(zz).getX()>=rightX_tableColumns[qq]) 
					&&(rightX_tableColumns[qq]<row.getCells().get(zz).getEndX()) )
						cells[tt][qq] = cells[tt][qq] +  row.getCells().get(zz).getText() + " ";
				}
			}
			for (int zz=0; zz<cellThisLine; zz++) {
				int crossCellBegin = 1000;
				int crossCellEnd=0;
				for (int qq=0; qq<cc; qq++) {
					if ( (qq>0) && (cells[tt][qq-1].length()==0) 
					&& (row.getCells().get(zz).getX()<=rightX_tableColumns[qq-1]) 
					&& (row.getCells().get(zz).getX()<leftX_tableColumns[qq]) 
					&& (leftX_tableColumns[qq]<=row.getCells().get(zz).getEndX()) ) {
						crossCells[tt][qq]=true;
						crossCells[tt][qq-1]=true;
						cells[tt][qq-1] = row.getCells().get(zz).getText();
						crossCellBegin=qq-1;
					}
					if ((qq<(cc-1)) && (cells[tt][qq+1].length()==0) 
					&& (row.getCells().get(zz).getEndX()>=leftX_tableColumns[qq+1]) 
					&& (row.getCells().get(zz).getX()<=rightX_tableColumns[qq]) 
					&& (rightX_tableColumns[qq]<row.getCells().get(zz).getEndX()) ) {
						crossCells[tt][qq]=true;
						crossCells[tt][qq+1]=true;
						cells[tt][qq+1] = row.getCells().get(zz).getText();
						crossCellEnd=qq+1;
					}
				}
			}
			for (int qq=0; qq<cc; qq++) {
				cellsThisLine = cellsThisLine + cells[tt][qq] + "; ";
			}
			contentsInRows[tt]=cellsThisLine;
		}
    	tc.setCells(cells);
    	tc.setCrossCells(crossCells);
    }
201410
public static void multAddOuter( double alpha , DMatrix4x4 A , double beta , DMatrix4 u , DMatrix4 v , DMatrix4x4 C ) {
        C.a11 = alpha*A.a11 + beta*u.a1*v.a1;
        C.a12 = alpha*A.a12 + beta*u.a1*v.a2;
        C.a13 = alpha*A.a13 + beta*u.a1*v.a3;
        C.a14 = alpha*A.a14 + beta*u.a1*v.a4;
        C.a21 = alpha*A.a21 + beta*u.a2*v.a1;
        C.a22 = alpha*A.a22 + beta*u.a2*v.a2;
        C.a23 = alpha*A.a23 + beta*u.a2*v.a3;
        C.a24 = alpha*A.a24 + beta*u.a2*v.a4;
        C.a31 = alpha*A.a31 + beta*u.a3*v.a1;
        C.a32 = alpha*A.a32 + beta*u.a3*v.a2;
        C.a33 = alpha*A.a33 + beta*u.a3*v.a3;
        C.a34 = alpha*A.a34 + beta*u.a3*v.a4;
        C.a41 = alpha*A.a41 + beta*u.a4*v.a1;
        C.a42 = alpha*A.a42 + beta*u.a4*v.a2;
        C.a43 = alpha*A.a43 + beta*u.a4*v.a3;
        C.a44 = alpha*A.a44 + beta*u.a4*v.a4;
    }
1113646
def balance_ss_model(F,L,Qc,H,Pinf,P0,dF=None,dQc=None,dPinf=None,dP0=None):
    

    bF,T,T_inv = balance_matrix(F)

    bL = np.dot( T_inv, L)
    bQc = Qc # not affected

    bH = np.dot(H, T)

    bPinf = np.dot(T_inv, np.dot(Pinf, T_inv.T))

    #import pdb; pdb.set_trace()
#    LL,islower = linalg.cho_factor(Pinf)
#    inds = np.triu_indices(Pinf.shape[0],k=1)
#    LL[inds] = 0.0
#    bLL = np.dot(T_inv, LL)
#    bPinf = np.dot( bLL, bLL.T)

    bP0 = np.dot(T_inv, np.dot(P0, T_inv.T))

    if dF is not None:
        bdF = np.zeros(dF.shape)
        for i in range(dF.shape[2]):
            bdF[:,:,i] = np.dot( T_inv, np.dot( dF[:,:,i], T))

    else:
        bdF = None

    if dPinf is not None:
        bdPinf = np.zeros(dPinf.shape)
        for i in range(dPinf.shape[2]):
            bdPinf[:,:,i] = np.dot( T_inv, np.dot( dPinf[:,:,i], T_inv.T))

#            LL,islower = linalg.cho_factor(dPinf[:,:,i])
#            inds = np.triu_indices(dPinf[:,:,i].shape[0],k=1)
#            LL[inds] = 0.0
#            bLL = np.dot(T_inv, LL)
#            bdPinf[:,:,i] = np.dot( bLL, bLL.T)


    else:
        bdPinf = None

    if dP0 is not None:
        bdP0 = np.zeros(dP0.shape)
        for i in range(dP0.shape[2]):
            bdP0[:,:,i] = np.dot( T_inv, np.dot( dP0[:,:,i], T_inv.T))
    else:
        bdP0 = None


    bdQc = dQc # not affected

    # (F,L,Qc,H,Pinf,P0,dF,dQc,dPinf,dP0)

    return bF, bL, bQc, bH, bPinf, bP0, bdF, bdQc, bdPinf, bdP0
943403
public static void fe_mul121666(int[] h,int[] f)
{
  int f0 = f[0];
  int f1 = f[1];
  int f2 = f[2];
  int f3 = f[3];
  int f4 = f[4];
  int f5 = f[5];
  int f6 = f[6];
  int f7 = f[7];
  int f8 = f[8];
  int f9 = f[9];
  long h0 = f0 * (long) 121666;
  long h1 = f1 * (long) 121666;
  long h2 = f2 * (long) 121666;
  long h3 = f3 * (long) 121666;
  long h4 = f4 * (long) 121666;
  long h5 = f5 * (long) 121666;
  long h6 = f6 * (long) 121666;
  long h7 = f7 * (long) 121666;
  long h8 = f8 * (long) 121666;
  long h9 = f9 * (long) 121666;
  long carry0;
  long carry1;
  long carry2;
  long carry3;
  long carry4;
  long carry5;
  long carry6;
  long carry7;
  long carry8;
  long carry9;

  carry9 = (h9 + (long) (1<<24)) >> 25; h0 += carry9 * 19; h9 -= carry9 << 25;
  carry1 = (h1 + (long) (1<<24)) >> 25; h2 += carry1; h1 -= carry1 << 25;
  carry3 = (h3 + (long) (1<<24)) >> 25; h4 += carry3; h3 -= carry3 << 25;
  carry5 = (h5 + (long) (1<<24)) >> 25; h6 += carry5; h5 -= carry5 << 25;
  carry7 = (h7 + (long) (1<<24)) >> 25; h8 += carry7; h7 -= carry7 << 25;

  carry0 = (h0 + (long) (1<<25)) >> 26; h1 += carry0; h0 -= carry0 << 26;
  carry2 = (h2 + (long) (1<<25)) >> 26; h3 += carry2; h2 -= carry2 << 26;
  carry4 = (h4 + (long) (1<<25)) >> 26; h5 += carry4; h4 -= carry4 << 26;
  carry6 = (h6 + (long) (1<<25)) >> 26; h7 += carry6; h6 -= carry6 << 26;
  carry8 = (h8 + (long) (1<<25)) >> 26; h9 += carry8; h8 -= carry8 << 26;

  h[0] = (int)h0;
  h[1] = (int)h1;
  h[2] = (int)h2;
  h[3] = (int)h3;
  h[4] = (int)h4;
  h[5] = (int)h5;
  h[6] = (int)h6;
  h[7] = (int)h7;
  h[8] = (int)h8;
  h[9] = (int)h9;
}
117216
public static void multAddOuter( double alpha , DMatrix3x3 A , double beta , DMatrix3 u , DMatrix3 v , DMatrix3x3 C ) {
        C.a11 = alpha*A.a11 + beta*u.a1*v.a1;
        C.a12 = alpha*A.a12 + beta*u.a1*v.a2;
        C.a13 = alpha*A.a13 + beta*u.a1*v.a3;
        C.a21 = alpha*A.a21 + beta*u.a2*v.a1;
        C.a22 = alpha*A.a22 + beta*u.a2*v.a2;
        C.a23 = alpha*A.a23 + beta*u.a2*v.a3;
        C.a31 = alpha*A.a31 + beta*u.a3*v.a1;
        C.a32 = alpha*A.a32 + beta*u.a3*v.a2;
        C.a33 = alpha*A.a33 + beta*u.a3*v.a3;
    }
272806
function(p) {
            var l = this.coeffs.length;
            for(var i=0; i<l; i++) {
                var c = this.coeffs[i];
                if(c < 0) { //go borrow
                    var b; //a coefficient > 0
                    for(var j=i; j<l; j++) {//starting from where we left off
                        if(this.coeffs[j] > 0) {
                            b = this.coeffs[j];
                            break;
                        }
                    }

                    if(b) { //if such a coefficient exists
                        for(j; j>i; j--) { //go down the line and adjust using p
                            this.coeffs[j] = this.coeffs[j].subtract(new Frac(1));
                            this.coeffs[j-1] = this.coeffs[j-1].add(new Frac(p));
                        }
                        c = this.coeffs[i]; //reset c
                    }
                }

                var d = c.mod(p);
                var w = c.subtract(d).divide(p);
                if(!w.equals(0)) {
                    var up_one = i+1;
                    var next = this.coeffs[up_one] || new Frac(0);
                    next = next.add(w);
                    this.coeffs[up_one] = new Frac(next);
                    this.coeffs[i] = new Frac(d);
                }
            }

            return this;
        }
343163
def _g(self, x, tau):
        
        return tau ** 2 * (tau ** 2 + 1) ** -2 * (
                (tau ** 2 + 1 + 2 * (x ** 2 - 1)) * self.F(x) + tau * np.pi + (tau ** 2 - 1) * np.log(tau) +
                np.sqrt(tau ** 2 + x ** 2) * (-np.pi + self.L(x, tau) * (tau ** 2 - 1) * tau ** -1))
806680
def ycbcr2rgb(y__, cb_, cr_):
    

    kb_ = 0.114
    kr_ = 0.299

    r__ = 2 * cr_ / (1 - kr_) + y__
    b__ = 2 * cb_ / (1 - kb_) + y__
    g__ = (y__ - kr_ * r__ - kb_ * b__) / (1 - kr_ - kb_)

    return r__, g__, b__