Q = I - gamma*u*u<sup>H</sup>
439432
public static ZMatrixRMaj householder(ZMatrixRMaj u , double gamma ) {
        int N = u.getDataLength()/2;
        // u*u^H
        ZMatrixRMaj uut = new ZMatrixRMaj(N,N);
        VectorVectorMult_ZDRM.outerProdH(u, u, uut);
        // foo = -gamma*u*u^H
        CommonOps_ZDRM.elementMultiply(uut,-gamma,0,uut);

        // I + foo
        for (int i = 0; i < N; i++) {
            int index = (i*uut.numCols+i)*2;
            uut.data[index] = 1 + uut.data[index];
        }

        return uut;
    }
46134
@Override
    public DMatrixRMaj getQ(DMatrixRMaj Q , boolean compact ) {
        if( compact ) {
            Q = UtilDecompositons_DDRM.checkIdentity(Q,numRows,minLength);
        } else {
            Q = UtilDecompositons_DDRM.checkIdentity(Q,numRows,numRows);
        }

        for( int j = minLength-1; j >= 0; j-- ) {
            double u[] = dataQR[j];

            double vv = u[j];
            u[j] = 1;
            QrHelperFunctions_DDRM.rank1UpdateMultR(Q, u, gammas[j], j, j, numRows, v);
            u[j] = vv;
        }

        return Q;
    }
1186889
def real_rotated_Gammas(m,l,phi1,phi2,theta1,theta2,gamma_ml):
    
    
    if m>0:
        ans=(1./sqrt(2))*(rotated_Gamma_ml(m,l,phi1,phi2,theta1,theta2,gamma_ml) + \
                          (-1)**m*rotated_Gamma_ml(-m,l,phi1,phi2,theta1,theta2,gamma_ml))
        return ans.real
    if m==0:
        return rotated_Gamma_ml(0,l,phi1,phi2,theta1,theta2,gamma_ml).real
    if m<0:
        ans=(1./sqrt(2)/complex(0.,1))*(rotated_Gamma_ml(-m,l,phi1,phi2,theta1,theta2,gamma_ml) - \
                                        (-1)**m*rotated_Gamma_ml(m,l,phi1,phi2,theta1,theta2,gamma_ml))
        return ans.real
1113853
def _c0(self):
        ""
        h_df = self.df / 2
        p, S = self._p, self.S

        return h_df * (logdet(S) + p * logtwo) + lpgamma(p, h_df)
579159
def crt(self, mp, mq):
        
        u = (mq - mp) * self.p_inverse % self.q
        return mp + (u * self.p)
500190
def _J(self):
        
        pd = self.particle_distribution(self._Ep * u.GeV)
        return pd.to("1/GeV").value
574169
private double[][] estimateGamma(double[][][] xi) {
        double[][] gamma = new double[xi.length + 1][numStates];

        for (int t = 0; t < xi.length; t++) {
            for (int i = 0; i < numStates; i++) {
                for (int j = 0; j < numStates; j++) {
                    gamma[t][i] += xi[t][i][j];
                }
            }
        }

        int n = xi.length - 1;
        for (int j = 0; j < numStates; j++) {
            for (int i = 0; i < numStates; i++) {
                gamma[xi.length][j] += xi[n][i][j];
            }
        }

        return gamma;
    }
394547
private static function hsl2rgb($h, $s, $l, &$r, &$g, &$b)
    {
        if ($s == 0) {
            $r = round($l * 255);
            $g = $r;
            $b = $r;

            return;
        }

        $temp1;
        if ($l < 0.5) {
            $temp1 = $l * (1.0 + $s);
        } else {
            $temp1 = $l + $s - $l * $s;
        }

        $temp2 = 2 * $l - $temp1;

        $h = $h / 360;

        $tempR = $h + (1 / 3);
        $tempG = $h;
        $tempB = $h - (1 / 3);

        if ($tempR < 0) {
            $tempR++;
        } elseif ($tempR > 1) {
            $tempR--;
        }

        if ($tempG < 0) {
            $tempG++;
        } elseif ($tempG > 1) {
            $tempG--;
        }

        if ($tempB < 0) {
            $tempB++;
        } elseif ($tempB > 1) {
            $tempB--;
        }

        if (6 * $tempR < 1) {
            $r = $temp2 + ($temp1 - $temp2) * 6 * $tempR;
        } elseif (2 * $tempR < 1) {
            $r = $temp1;
        } elseif (3 * $tempR < 2) {
            $r = $temp2 + ($temp1 - $temp2) * ((2 / 3) - $tempR) * 6;
        } else {
            $r = $temp2;
        }

        if (6 * $tempG < 1) {
            $g = $temp2 + ($temp1 - $temp2) * 6 * $tempG;
        } elseif (2 * $tempG < 1) {
            $g = $temp1;
        } elseif (3 * $tempG < 2) {
            $g = $temp2 + ($temp1 - $temp2) * ((2 / 3) - $tempG) * 6;
        } else {
            $g = $temp2;
        }

        if (6 * $tempB < 1) {
            $b = $temp2 + ($temp1 - $temp2) * 6 * $tempB;
        } elseif (2 * $tempB < 1) {
            $b = $temp1;
        } elseif (3 * $tempB < 2) {
            $b = $temp2 + ($temp1 - $temp2) * ((2 / 3) - $tempB) * 6;
        } else {
            $b = $temp2;
        }

        $r = round($r * 255);
        $g = round($g * 255);
        $b = round($b * 255);
    }
373789
def get_Q(self):
        
        Z = self.Z.detach().clone().numpy()
        O = self.O.numpy()
        I_k = np.eye(self.k)
        return O @ Z @ np.linalg.inv(I_k + Z.T @ O @ Z) @ Z.T @ O
81036
def F(self, **kwargs):
        '''
        
        '''
        Eg = self.Eg_Gamma_0(**kwargs)
        Delta_SO = self.Delta_SO(**kwargs)
        Ep = self.Ep(**kwargs)
        meff = self.meff_e_Gamma_0(**kwargs)
        return (1./meff-1-(Ep*(Eg+2.*Delta_SO/3.))/(Eg*(Eg+Delta_SO)))/2
228652
def Gamma(k, theta, tag=None):
    
    assert (
        k > 0 and theta > 0
    ), 'Gamma "k" and "theta" parameters must be greater than zero'
    return uv(ss.gamma(k, scale=theta), tag=tag)
1370914
def lpgamma(p, a):
    
    res = p * (p - 1.0) / 4.0 * log(pi)
    for ii in range(1, p+1):
        res += gammaln(a + (1.0 - ii) / 2.0)

    return res
817536
def sX(qubit: Qubit, coefficient: complex = 1.0) -> Pauli:
    """"""
    return Pauli.sigma(qubit, 'X', coefficient)
1183589
def sigma_r2(self, r, a, gamma, rho0_r0_gamma, r_ani):
        
        # first term
        prefac1 = 4*np.pi * const.G * a**(-gamma) * rho0_r0_gamma / (3-gamma)
        prefac2 = r * (r + a)**3/(r**2 + r_ani**2)
        hyp1 = vel_util.hyp_2F1(a=2+gamma, b=gamma, c=3+gamma, z=1./(1+r/a))
        hyp2 = vel_util.hyp_2F1(a=3, b=gamma, c=1+gamma, z=-a/r)
        fac = r_ani**2/a**2 * hyp1 / ((2+gamma) * (r/a + 1)**(2+gamma)) + hyp2 / (gamma*(r/a)**gamma)
        return prefac1 * prefac2 * fac * (self._cosmo.arcsec2phys_lens(1.) * const.Mpc / 1000) ** 2
1035691
public static double gammToSigma(double gamma)
    {
        if(gamma <= 0 || Double.isNaN(gamma) || Double.isInfinite(gamma))
            throw new IllegalArgumentException("gamma must be positive, not " + gamma);
        return 1/Math.sqrt(2*gamma);
    }
950571
def calculate_u(dt, dx, dy, u, um, k):
    up = numpy.zeros((u.shape[0], u.shape[1]))
    ""
    for i in xrange(1, u.shape[0]-1):
        for j in xrange(1, u.shape[1]-1):
            up[i,j] = 2*u[i,j] - um[i,j] + \
                (dt/dx)**2*(
                        (0.5*(k[i+1,j] + k[i,j])*(u[i+1,j] - u[i,j]) -
                            0.5*(k[i,j] + k[i-1,j])*(u[i,j] - u[i-1,j]))) + \
                        (dt/dy)**2*(
                                (0.5*(k[i,j+1] + k[i,j])*(u[i,j+1] - u[i,j]) -
                                    0.5*(k[i,j] + k[i,j-1])*(u[i,j] - u[i,j-1])))
    return up
855963
function gamma( x ) {
	var w;
	var y;
	var v;

	w = 1 / x;
	w = 1 + w * polyval( w );
	y = exp( x );

	// Check `x` to avoid `pow()` overflow...
	if ( x > MAX_STIRLING ) {
		v = pow( x, 0.5*x - 0.25 );
		y = v * (v/y);
	} else {
		y = pow( x, x-0.5 ) / y;
	}
	return SQRT2PI * y * w;
}