Ensure random weight matrix is valid.

    TODO:  The diagonally dominant tuning currently doesn't make sense.
           Our weight matrix has zeros along the diagonal, so multiplying by
           a diagonal matrix results in a zero-matrix.
431904
def _fix_weights(weight_fun, *args):
    
    weights = weight_fun(*args)

    # TODO: fix this
    # disable checks for now
    return weights

    # if positive semidefinite, then we're good as is
    if _check_psd(weights):
        return weights

    # make diagonally dominant
    off_diag_sums = np.sum(weights, axis=1)  # NOTE: assumes diag is zero
    mod_mat = np.linalg.inv(np.sqrt(np.diag(off_diag_sums)))
    return np.dot(mod_mat, weights, mod_mat)
661185
public static boolean cholL( DMatrix2x2 A ) {

        A.a11 = Math.sqrt(A.a11);
        A.a12 = 0;
        A.a21 = (A.a21)/A.a11;
        A.a22 = Math.sqrt(A.a22-A.a21*A.a21);
        return !UtilEjml.isUncountable(A.a22);
    }
96453
def _eig_complex_symmetric(M: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    
    if not np.allclose(M, M.transpose()):
        raise np.linalg.LinAlgError('Not a symmetric matrix')

    # The matrix of eigenvectors should be orthogonal.
    # But the standard 'eig' method will fail to return an orthogonal
    # eigenvector matrix when the eigenvalues are degenerate. However,
    # both the real and
    # imaginary part of M must be symmetric with the same orthogonal
    # matrix of eigenvectors. But either the real or imaginary part could
    # vanish. So we use a randomized algorithm where we diagonalize a
    # random linear combination of real and imaginary parts to find the
    # eigenvectors, taking advantage of the 'eigh' subroutine for
    # diagonalizing symmetric matrices.
    # This can fail if we're very unlucky with our random coefficient, so we
    # give the algorithm a few chances to succeed.

    # Empirically, never seems to fail on randomly sampled complex
    # symmetric 4x4 matrices.
    # If failure rate is less than 1 in a million, then 16 rounds
    # will have overall failure rate less than 1 in a googol.
    # However, cannot (yet) guarantee that there aren't special cases
    # which have much higher failure rates.

    # GEC 2018

    max_attempts = 16
    for _ in range(max_attempts):
        c = np.random.uniform(0, 1)
        matrix = c * M.real + (1-c) * M.imag
        _, eigvecs = np.linalg.eigh(matrix)
        eigvecs = np.array(eigvecs, dtype=complex)
        eigvals = np.diag(eigvecs.transpose() @ M @ eigvecs)

        # Finish if we got a correct answer.
        reconstructed = eigvecs @ np.diag(eigvals) @ eigvecs.transpose()
        if np.allclose(M, reconstructed):
            return eigvals, eigvecs

    # Should never happen. Hopefully.
    raise np.linalg.LinAlgError(
        'Cannot diagonalize complex symmetric matrix.')
481934
def initialize_weights(self):
        """"""
        n = self._outputSize
        m = self._inputSize
        self._Q = self._random.sample((n,m))

        # Normalize the weights of each units
        for i in range(n):
            self._Q[i] /=  np.sqrt( np.dot(self._Q[i], self._Q[i]) )
880892
public static DMatrixRMaj symmetricPosDef(int width, Random rand) {
        // This is not formally proven to work.  It just seems to work.
        DMatrixRMaj a = new DMatrixRMaj(width,1);
        DMatrixRMaj b = new DMatrixRMaj(width,width);

        for( int i = 0; i < width; i++ ) {
            a.set(i,0,rand.nextDouble());
        }

        CommonOps_DDRM.multTransB(a,a,b);

        for( int i = 0; i < width; i++ ) {
            b.add(i,i,1);
        }

        return b;
    }
1173424
function skewness( out, x ) {
	var len = x.length,
		i;
	if ( out.length !== len ) {
		throw new Error( 'skewness()::invalid input arguments. Input and output matrices must be the same length.' );
	}
	for ( i = 0; i < len; i++ ) {
		out.data[ i ] = SKEWNESS( x.data[ i ] );
	}
	return out;
}
272591
@Override
    public void setWeights(Map<String, INDArray> weights) throws InvalidKerasConfigurationException {
        this.weights = new HashMap<>();
        if (center) {
            if (weights.containsKey(PARAM_NAME_BETA))
                this.weights.put(BatchNormalizationParamInitializer.BETA, weights.get(PARAM_NAME_BETA));
            else
                throw new InvalidKerasConfigurationException("Parameter " + PARAM_NAME_BETA + " does not exist in weights");
        } else {
            INDArray dummyBeta = Nd4j.zerosLike(weights.get(PARAM_NAME_BETA));
            this.weights.put(BatchNormalizationParamInitializer.BETA, dummyBeta);
        }
        if (scale) {
            if (weights.containsKey(PARAM_NAME_GAMMA))
                this.weights.put(BatchNormalizationParamInitializer.GAMMA, weights.get(PARAM_NAME_GAMMA));
            else
                throw new InvalidKerasConfigurationException(
                        "Parameter " + PARAM_NAME_GAMMA + " does not exist in weights");
        } else {
            INDArray dummyGamma = weights.containsKey(PARAM_NAME_GAMMA)
                    ? Nd4j.onesLike(weights.get(PARAM_NAME_GAMMA))
                    : Nd4j.onesLike(weights.get(PARAM_NAME_BETA));
            this.weights.put(BatchNormalizationParamInitializer.GAMMA, dummyGamma);
        }
        if (weights.containsKey(conf.getLAYER_FIELD_BATCHNORMALIZATION_MOVING_MEAN()))
            this.weights.put(BatchNormalizationParamInitializer.GLOBAL_MEAN, weights.get(conf.getLAYER_FIELD_BATCHNORMALIZATION_MOVING_MEAN()));
        else
            throw new InvalidKerasConfigurationException(
                    "Parameter " + conf.getLAYER_FIELD_BATCHNORMALIZATION_MOVING_MEAN() + " does not exist in weights");
        if (weights.containsKey(conf.getLAYER_FIELD_BATCHNORMALIZATION_MOVING_VARIANCE()))
            this.weights.put(BatchNormalizationParamInitializer.GLOBAL_VAR, weights.get(conf.getLAYER_FIELD_BATCHNORMALIZATION_MOVING_VARIANCE()));
        else
            throw new InvalidKerasConfigurationException(
                    "Parameter " + conf.getLAYER_FIELD_BATCHNORMALIZATION_MOVING_VARIANCE() + " does not exist in weights");
        if (weights.size() > 4) {
            Set<String> paramNames = weights.keySet();
            paramNames.remove(PARAM_NAME_BETA);
            paramNames.remove(PARAM_NAME_GAMMA);
            paramNames.remove(conf.getLAYER_FIELD_BATCHNORMALIZATION_MOVING_MEAN());
            paramNames.remove(conf.getLAYER_FIELD_BATCHNORMALIZATION_MOVING_VARIANCE());
            String unknownParamNames = paramNames.toString();
            log.warn("Attempting to set weights for unknown parameters: "
                    + unknownParamNames.substring(1, unknownParamNames.length() - 1));
        }
    }
559630
public double[] buildWeightMatrix(double[][] nodeMatrix, double[][] bondMatrix) {
        for (int i = 0; i < nodeNumber; i++) {
            weight[i + 1] = nodeMatrix[i][0];
            for (int j = 0; j < layerNumber; j++) {
                weight[i + 1] += nodeMatrix[i][j + 1] * bondMatrix[i][j] * Math.pow(10.0, (double) -(j + 1));
            }
        }
        weight[0] = 0.0;
        return weight;
    }
448620
def initialize(self):
        
        self.randomize()
        self.dweight = Numeric.zeros((self.fromLayer.size, \
                                      self.toLayer.size), 'f')
        self.wed = Numeric.zeros((self.fromLayer.size, \
                                  self.toLayer.size), 'f')
        self.wedLast = Numeric.zeros((self.fromLayer.size, \
                                      self.toLayer.size), 'f')
145357
function selectWithoutValidation(choices) {
	// Generate a list of options with positive non-zero chances
	let choicesWithNonZeroChances = [];
	let totalWeight = 0.0;
	for(let i = 0; i < choices.length; i++) {
		if(choices[i].chance > 0.0) {
			choicesWithNonZeroChances.push(choices[i]);
			totalWeight += choices[i].chance;
		}
	}

	// Pick a random value from [0,1)
	let value = Math.random() * totalWeight;

	// Iterate over possibilities until we find the selected one
	let chanceCovered = 0.0;
	for(let i = 0; i < choicesWithNonZeroChances.length; i++) {
		chanceCovered += choicesWithNonZeroChances[i].chance;
		if(value < chanceCovered) {
			return choicesWithNonZeroChances[i].result;
		}
	}
}
1218898
public static boolean cholL( DMatrix4x4 A ) {

        A.a11 = Math.sqrt(A.a11);
        A.a12 = 0;
        A.a13 = 0;
        A.a14 = 0;
        A.a21 = (A.a21)/A.a11;
        A.a22 = Math.sqrt(A.a22-A.a21*A.a21);
        A.a23 = 0;
        A.a24 = 0;
        A.a31 = (A.a31)/A.a11;
        A.a32 = (A.a32-A.a31*A.a21)/A.a22;
        A.a33 = Math.sqrt(A.a33-A.a31*A.a31-A.a32*A.a32);
        A.a34 = 0;
        A.a41 = (A.a41)/A.a11;
        A.a42 = (A.a42-A.a41*A.a21)/A.a22;
        A.a43 = (A.a43-A.a41*A.a31-A.a42*A.a32)/A.a33;
        A.a44 = Math.sqrt(A.a44-A.a41*A.a41-A.a42*A.a42-A.a43*A.a43);
        return !UtilEjml.isUncountable(A.a44);
    }
123708
@Override
    public void learn(double[] x, double y) {
        if (x.length != p) {
            throw new IllegalArgumentException(String.format("Invalid input vector size: %d, expected: %d", x.length, p));
        }

        System.arraycopy(x, 0, x1, 0, p);
        double v = 1 + V.xax(x1);
        // If 1/v is NaN, then the update to V will no longer be invertible.
        // See https://en.wikipedia.org/wiki/Sherman%E2%80%93Morrison_formula#Statement
        if (Double.isNaN(1/v)){
            throw new IllegalStateException("The updated V matrix is no longer invertible.");
        }

        V.ax(x1, Vx);
        for (int j = 0; j <= p; j++) {
            for (int i = 0; i <= p; i++) {
                double tmp = V.get(i, j) - ((Vx[i] * Vx[j])/v);
                V.set(i, j, tmp/lambda);
            }
        }

        // V has been updated. Compute Vx again.
        V.ax(x1, Vx);
        
        double err = y - predict(x);
        for (int i = 0; i <= p; i++){
            w[i] += Vx[i] * err;
        }
    }
232250
public final void setFromMatrix(Matrix4f m1) {
		double ww = 0.25f*(m1.m00 + m1.m11 + m1.m22 + m1.m33);

		if (ww >= 0) {
			if (ww >= EPS2) {
				this.w =  Math.sqrt(ww);
				ww =  0.25f/this.w;
				this.x = (m1.m21 - m1.m12)*ww;
				this.y = (m1.m02 - m1.m20)*ww;
				this.z = (m1.m10 - m1.m01)*ww;
				return;
			} 
		}
		else {
			this.w = 0;
			this.x = 0;
			this.y = 0;
			this.z = 1;
			return;
		}

		this.w = 0;
		ww = -0.5f*(m1.m11 + m1.m22);

		if (ww >= 0) {
			if (ww >= EPS2) {
				this.x =  Math.sqrt(ww);
				ww = 1.0f/(2.0f*this.x);
				this.y = m1.m10*ww;
				this.z = m1.m20*ww;
				return;
			}
		} else {
			this.x = 0;
			this.y = 0;
			this.z = 1;
			return;
		}

		this.x = 0;
		ww = 0.5f*(1.0f - m1.m22);

		if (ww >= EPS2) {
			this.y =  Math.sqrt(ww);
			this.z = m1.m21/(2.0f*this.y);
			return;
		}

		this.y = 0;
		this.z = 1;
	}
234903
private void initMatrix(int quality) {
        double[] AANscaleFactor = {1.0, 1.387039845, 1.306562965, 1.175875602,
                1.0, 0.785694958, 0.541196100, 0.275899379};
        int i;
        int j;
        int index;
        int Quality;
        int temp;

// converting quality setting to that specified in the jpeg_quality_scaling
// method in the IJG Jpeg-6a C libraries

        Quality = quality;
        if (Quality <= 0)
            Quality = 1;
        if (Quality > 100)
            Quality = 100;
        if (Quality < 50)
            Quality = 5000 / Quality;
        else
            Quality = 200 - Quality * 2;

// Creating the luminance matrix

        quantum_luminance[0] = 16;
        quantum_luminance[1] = 11;
        quantum_luminance[2] = 10;
        quantum_luminance[3] = 16;
        quantum_luminance[4] = 24;
        quantum_luminance[5] = 40;
        quantum_luminance[6] = 51;
        quantum_luminance[7] = 61;
        quantum_luminance[8] = 12;
        quantum_luminance[9] = 12;
        quantum_luminance[10] = 14;
        quantum_luminance[11] = 19;
        quantum_luminance[12] = 26;
        quantum_luminance[13] = 58;
        quantum_luminance[14] = 60;
        quantum_luminance[15] = 55;
        quantum_luminance[16] = 14;
        quantum_luminance[17] = 13;
        quantum_luminance[18] = 16;
        quantum_luminance[19] = 24;
        quantum_luminance[20] = 40;
        quantum_luminance[21] = 57;
        quantum_luminance[22] = 69;
        quantum_luminance[23] = 56;
        quantum_luminance[24] = 14;
        quantum_luminance[25] = 17;
        quantum_luminance[26] = 22;
        quantum_luminance[27] = 29;
        quantum_luminance[28] = 51;
        quantum_luminance[29] = 87;
        quantum_luminance[30] = 80;
        quantum_luminance[31] = 62;
        quantum_luminance[32] = 18;
        quantum_luminance[33] = 22;
        quantum_luminance[34] = 37;
        quantum_luminance[35] = 56;
        quantum_luminance[36] = 68;
        quantum_luminance[37] = 109;
        quantum_luminance[38] = 103;
        quantum_luminance[39] = 77;
        quantum_luminance[40] = 24;
        quantum_luminance[41] = 35;
        quantum_luminance[42] = 55;
        quantum_luminance[43] = 64;
        quantum_luminance[44] = 81;
        quantum_luminance[45] = 104;
        quantum_luminance[46] = 113;
        quantum_luminance[47] = 92;
        quantum_luminance[48] = 49;
        quantum_luminance[49] = 64;
        quantum_luminance[50] = 78;
        quantum_luminance[51] = 87;
        quantum_luminance[52] = 103;
        quantum_luminance[53] = 121;
        quantum_luminance[54] = 120;
        quantum_luminance[55] = 101;
        quantum_luminance[56] = 72;
        quantum_luminance[57] = 92;
        quantum_luminance[58] = 95;
        quantum_luminance[59] = 98;
        quantum_luminance[60] = 112;
        quantum_luminance[61] = 100;
        quantum_luminance[62] = 103;
        quantum_luminance[63] = 99;

        for (j = 0; j < 64; j++) {
            temp = (quantum_luminance[j] * Quality + 50) / 100;
            if (temp <= 0) temp = 1;
            if (temp > 255) temp = 255;
            quantum_luminance[j] = temp;
        }
        index = 0;
        for (i = 0; i < 8; i++) {
            for (j = 0; j < 8; j++) {
// The divisors for the LL&M method (the slow integer method used in
// jpeg 6a library).  This method is currently (04/04/98) incompletely
// implemented.
//                        DivisorsLuminance[index] = ((double) quantum_luminance[index]) << 3;
// The divisors for the AAN method (the float method used in jpeg 6a library.
                DivisorsLuminance[index] = (1.0 / ((double) quantum_luminance[index] * AANscaleFactor[i] * AANscaleFactor[j] * 8.0));
                index++;
            }
        }


// Creating the chrominance matrix

        quantum_chrominance[0] = 17;
        quantum_chrominance[1] = 18;
        quantum_chrominance[2] = 24;
        quantum_chrominance[3] = 47;
        quantum_chrominance[4] = 99;
        quantum_chrominance[5] = 99;
        quantum_chrominance[6] = 99;
        quantum_chrominance[7] = 99;
        quantum_chrominance[8] = 18;
        quantum_chrominance[9] = 21;
        quantum_chrominance[10] = 26;
        quantum_chrominance[11] = 66;
        quantum_chrominance[12] = 99;
        quantum_chrominance[13] = 99;
        quantum_chrominance[14] = 99;
        quantum_chrominance[15] = 99;
        quantum_chrominance[16] = 24;
        quantum_chrominance[17] = 26;
        quantum_chrominance[18] = 56;
        quantum_chrominance[19] = 99;
        quantum_chrominance[20] = 99;
        quantum_chrominance[21] = 99;
        quantum_chrominance[22] = 99;
        quantum_chrominance[23] = 99;
        quantum_chrominance[24] = 47;
        quantum_chrominance[25] = 66;
        quantum_chrominance[26] = 99;
        quantum_chrominance[27] = 99;
        quantum_chrominance[28] = 99;
        quantum_chrominance[29] = 99;
        quantum_chrominance[30] = 99;
        quantum_chrominance[31] = 99;
        quantum_chrominance[32] = 99;
        quantum_chrominance[33] = 99;
        quantum_chrominance[34] = 99;
        quantum_chrominance[35] = 99;
        quantum_chrominance[36] = 99;
        quantum_chrominance[37] = 99;
        quantum_chrominance[38] = 99;
        quantum_chrominance[39] = 99;
        quantum_chrominance[40] = 99;
        quantum_chrominance[41] = 99;
        quantum_chrominance[42] = 99;
        quantum_chrominance[43] = 99;
        quantum_chrominance[44] = 99;
        quantum_chrominance[45] = 99;
        quantum_chrominance[46] = 99;
        quantum_chrominance[47] = 99;
        quantum_chrominance[48] = 99;
        quantum_chrominance[49] = 99;
        quantum_chrominance[50] = 99;
        quantum_chrominance[51] = 99;
        quantum_chrominance[52] = 99;
        quantum_chrominance[53] = 99;
        quantum_chrominance[54] = 99;
        quantum_chrominance[55] = 99;
        quantum_chrominance[56] = 99;
        quantum_chrominance[57] = 99;
        quantum_chrominance[58] = 99;
        quantum_chrominance[59] = 99;
        quantum_chrominance[60] = 99;
        quantum_chrominance[61] = 99;
        quantum_chrominance[62] = 99;
        quantum_chrominance[63] = 99;

        for (j = 0; j < 64; j++) {
            temp = (quantum_chrominance[j] * Quality + 50) / 100;
            if (temp <= 0) temp = 1;
            if (temp >= 255) temp = 255;
            quantum_chrominance[j] = temp;
        }
        index = 0;
        for (i = 0; i < 8; i++) {
            for (j = 0; j < 8; j++) {
// The divisors for the LL&M method (the slow integer method used in
// jpeg 6a library).  This method is currently (04/04/98) incompletely
// implemented.
//                        DivisorsChrominance[index] = ((double) quantum_chrominance[index]) << 3;
// The divisors for the AAN method (the float method used in jpeg 6a library.
                DivisorsChrominance[index] = (1.0 / ((double) quantum_chrominance[index] * AANscaleFactor[i] * AANscaleFactor[j] * 8.0));
                index++;
            }
        }

// quantum and Divisors are objects used to hold the appropriate matices

        quantum[0] = quantum_luminance;
        Divisors[0] = DivisorsLuminance;
        quantum[1] = quantum_chrominance;
        Divisors[1] = DivisorsChrominance;


    }
623391
public static boolean hasNaN( DMatrixD1 m )
    {
        int length = m.getNumElements();

        for( int i = 0; i < length; i++ ) {
            if( Double.isNaN(m.get(i)))
                return true;
        }
        return false;
    }
1321172
private void fillKsubMatrix( Coordinate[] controlPoints, GeneralMatrix L ) {
        double alfa = 0;
        int controlPointsNum = controlPoints.length;
        for( int i = 0; i < controlPointsNum; i++ ) {
            for( int j = i + 1; j < controlPointsNum; j++ ) {
                double u = calculateFunctionU(controlPoints[i], controlPoints[j]);
                L.setElement(i, j, u);
                L.setElement(j, i, u);
                alfa = alfa + (u * 2); // same for upper and lower part
            }
        }

        alfa = alfa / (controlPointsNum * controlPointsNum);
    }
59872
private double computeExplainedVariance(double[] eigenValues, int filteredEigenPairs) {
    double strongsum = 0., weaksum = 0.;
    for(int i = 0; i < filteredEigenPairs; i++) {
      strongsum += eigenValues[i];
    }
    for(int i = filteredEigenPairs; i < eigenValues.length; i++) {
      weaksum += eigenValues[i];
    }
    return strongsum / (strongsum + weaksum);
  }