Golang FloatMatrix.Row方法代碼示例

// Solve secular function arising in symmetric eigenproblems. On exit 'Y' contains new
// eigenvalues and 'V' the rank-one update vector corresponding new eigenvalues.
// The matrix Qd holds for each eigenvalue then computed deltas as row vectors.
// On entry 'D' holds original eigenvalues and 'Z' is the rank-one update vector.
func TRDSecularSolveAll(y, v, Qd, d, z *cmat.FloatMatrix, rho float64, confs ...*gomas.Config) (err *gomas.Error) {
	var delta cmat.FloatMatrix
	var lmbda float64
	var e, ei int

	ei = 0
	err = nil
	if y.Len() != d.Len() || z.Len() != d.Len() || m(Qd) != n(Qd) || m(Qd) != d.Len() {
		err = gomas.NewError(gomas.ESIZE, "TRDSecularSolveAll")
		return
	}
	for i := 0; i < d.Len(); i++ {
		delta.Row(Qd, i)
		lmbda, e = trdsecRoot(d, z, &delta, i, rho)
		if e < 0 && ei == 0 {
			ei = -(i + 1)
		}
		y.SetAt(i, lmbda)
	}
	if ei == 0 {
		trdsecUpdateVecDelta(v, Qd, d, rho)
	} else {
		err = gomas.NewError(gomas.ECONVERGE, "TRDSecularSolveAll", ei)
	}
	return
}

func sortEigenVec(D, U, V, C *cmat.FloatMatrix, updown int) {
	var sD, m0, m1 cmat.FloatMatrix

	N := D.Len()
	for k := 0; k < N-1; k++ {
		sD.SubVector(D, k, N-k)
		pk := vecMinMax(&sD, -updown)
		if pk != 0 {
			t0 := D.GetAt(k)
			D.SetAt(k, D.GetAt(pk+k))
			D.SetAt(k+pk, t0)
			if U != nil {
				m0.Column(U, k)
				m1.Column(U, k+pk)
				blasd.Swap(&m1, &m0)
			}
			if V != nil {
				m0.Row(V, k)
				m1.Row(V, k+pk)
				blasd.Swap(&m1, &m0)
			}
			if C != nil {
				m0.Column(C, k)
				m1.Column(C, k+pk)
				blasd.Swap(&m1, &m0)
			}
		}
	}
}

/*
 * Compute
 *   B = B*diag(D).-1      flags & RIGHT == true
 *   B = diag(D).-1*B      flags & LEFT  == true
 *
 * If flags is LEFT (RIGHT) then element-wise divides columns (rows) of B with vector D.
 *
 * Arguments:
 *   B     M-by-N matrix if flags&RIGHT == true or N-by-M matrix if flags&LEFT == true
 *
 *   D     N element column or row vector or N-by-N matrix
 *
 *   flags Indicator bits, LEFT or RIGHT
 */
func SolveDiag(B, D *cmat.FloatMatrix, flags int, confs ...*gomas.Config) *gomas.Error {
	var c, d0 cmat.FloatMatrix
	var d *cmat.FloatMatrix

	conf := gomas.CurrentConf(confs...)
	d = D
	if !D.IsVector() {
		d0.Diag(D)
		d = &d0
	}
	dn := d0.Len()
	br, bc := B.Size()
	switch flags & (gomas.LEFT | gomas.RIGHT) {
	case gomas.LEFT:
		if br != dn {
			return gomas.NewError(gomas.ESIZE, "SolveDiag")
		}
		// scale rows;
		for k := 0; k < dn; k++ {
			c.Row(B, k)
			blasd.InvScale(&c, d.GetAt(k), conf)
		}
	case gomas.RIGHT:
		if bc != dn {
			return gomas.NewError(gomas.ESIZE, "SolveDiag")
		}
		// scale columns
		for k := 0; k < dn; k++ {
			c.Column(B, k)
			blasd.InvScale(&c, d.GetAt(k), conf)
		}
	}
	return nil
}

// Compute eigenmatrix Q for updated eigenvalues in 'dl'.
func trdsecEigenBuild(Q, z, Q2 *cmat.FloatMatrix) {
	var qi, delta cmat.FloatMatrix

	for k := 0; k < z.Len(); k++ {
		qi.Column(Q, k)
		delta.Row(Q2, k)
		trdsecEigenVecDelta(&qi, &delta, z)
	}
}

func mNormInf(A *cmat.FloatMatrix) float64 {
	var amax float64 = 0.0
	var row cmat.FloatMatrix
	arows, _ := A.Size()
	for k := 0; k < arows; k++ {
		row.Row(A, k)
		rmax := blasd.ASum(&row)
		if rmax > amax {
			amax = rmax
		}
	}
	return amax
}

/*
 * Generate one of the orthogonal matrices Q or P.T determined by BDReduce() when
 * reducing a real matrix A to bidiagonal form. Q and P.T are defined as products
 * elementary reflectors H(i) or G(i) respectively.
 *
 * Orthogonal matrix Q is generated if flag WANTQ is set. And matrix P respectively
 * if flag WANTP is set.
 */
func BDBuild(A, tau, W *cmat.FloatMatrix, K, flags int, confs ...*gomas.Config) *gomas.Error {
	var Qh, Ph, tauh, d, s cmat.FloatMatrix
	var err *gomas.Error = nil

	if m(A) == 0 || n(A) == 0 {
		return nil
	}

	if m(A) > n(A) || (m(A) == n(A) && flags&gomas.LOWER == 0) {
		switch flags & (gomas.WANTQ | gomas.WANTP) {
		case gomas.WANTQ:
			tauh.SubMatrix(tau, 0, 0, n(A), 1)
			err = QRBuild(A, &tauh, W, K, confs...)

		case gomas.WANTP:
			// Shift P matrix embedded in A down and fill first column and row
			// to unit vector
			for j := n(A) - 1; j > 0; j-- {
				s.SubMatrix(A, j-1, j, 1, n(A)-j)
				d.SubMatrix(A, j, j, 1, n(A)-j)
				blasd.Copy(&d, &s)
				A.Set(j, 0, 0.0)
			}
			// zero  first row and set first entry to one
			d.Row(A, 0)
			blasd.Scale(&d, 0.0)
			d.Set(0, 0, 1.0)

			Ph.SubMatrix(A, 1, 1, n(A)-1, n(A)-1)
			tauh.SubMatrix(tau, 0, 0, n(A)-1, 1)
			if K > n(A)-1 {
				K = n(A) - 1
			}
			err = LQBuild(&Ph, &tauh, W, K, confs...)
		}
	} else {
		switch flags & (gomas.WANTQ | gomas.WANTP) {
		case gomas.WANTQ:
			// Shift Q matrix embedded in A right and fill first column and row
			// to unit vector
			for j := m(A) - 1; j > 0; j-- {
				s.SubMatrix(A, j, j-1, m(A)-j, 1)
				d.SubMatrix(A, j, j, m(A)-j, 1)
				blasd.Copy(&d, &s)
				A.Set(0, j, 0.0)
			}
			// zero first column and set first entry to one
			d.Column(A, 0)
			blasd.Scale(&d, 0.0)
			d.Set(0, 0, 1.0)

			Qh.SubMatrix(A, 1, 1, m(A)-1, m(A)-1)
			tauh.SubMatrix(tau, 0, 0, m(A)-1, 1)
			if K > m(A)-1 {
				K = m(A) - 1
			}
			err = QRBuild(&Qh, &tauh, W, K, confs...)

		case gomas.WANTP:
			tauh.SubMatrix(tau, 0, 0, m(A), 1)
			err = LQBuild(A, &tauh, W, K, confs...)
		}
	}
	if err != nil {
		err.Update("BDBuild")
	}
	return err
}

//.........這裏部分代碼省略.........
		if bits&gomas.WANTV != 0 {
			if m(A) == 1 {
				if m(V) == m(A) {
					// V is 1-by-N
					V.Copy(A)
					LQBuild(V, tau, W, m(A), conf)
				} else {
					// V is N-by-N
					eye := cmat.FloatDiagonalSource{1.0}
					V.SetFrom(&eye, cmat.SYMM)
					if err = LQMult(V, A, tau, W, gomas.RIGHT, conf); err != nil {
						return
					}
				}
			} else {
				// V is 1-by-1
				V.Set(0, 0, -1.0)
			}
		}
		return
	}

	// Use bdSvd2x2 functions
	d0 = A.Get(0, 0)
	d1 = A.Get(1, 1)
	if m(A) >= n(A) {
		e0 = A.Get(0, 1)
	} else {
		e0 = A.Get(1, 0)
	}

	if bits&(gomas.WANTU|gomas.WANTV) == 0 {
		// no vectors
		smin, smax := bdSvd2x2(d0, e0, d1)
		S.SetAt(0, math.Abs(smax))
		S.SetAt(1, math.Abs(smin))
		return
	}

	// at least left or right eigenvector wanted
	smin, smax, cosl, sinl, cosr, sinr := bdSvd2x2Vec(d0, e0, d1)
	if bits&gomas.WANTU != 0 {
		// left eigenvectors
		if m(A) >= n(A) {
			if n(U) == n(A) {
				U.Copy(A)
				if err = QRBuild(U, tau, W, n(A), conf); err != nil {
					return
				}
			} else {
				// U is M-by-M
				eye := cmat.FloatDiagonalSource{1.0}
				U.SetFrom(&eye, cmat.SYMM)
				if err = QRMult(U, A, tau, W, gomas.RIGHT, conf); err != nil {
					return
				}
			}
			ApplyGivensRight(U, 0, 1, 0, m(A), cosl, sinl)
		} else {
			// U is 2-by-2
			eye := cmat.FloatDiagonalSource{1.0}
			U.SetFrom(&eye, cmat.SYMM)
			ApplyGivensRight(U, 0, 1, 0, m(A), cosr, sinr)
		}
	}

	if bits&gomas.WANTV != 0 {
		if n(A) > m(A) {
			if m(V) == m(A) {
				V.Copy(A)
				if err = LQBuild(V, tau, W, m(A), conf); err != nil {
					return
				}
			} else {
				eye := cmat.FloatDiagonalSource{1.0}
				V.SetFrom(&eye, cmat.SYMM)
				if err = LQMult(V, A, tau, W, gomas.RIGHT, conf); err != nil {
					return
				}
			}
			ApplyGivensLeft(V, 0, 1, 0, n(A), cosl, sinl)
		} else {
			// V is 2-by-2
			eye := cmat.FloatDiagonalSource{1.0}
			V.SetFrom(&eye, cmat.SYMM)
			ApplyGivensLeft(V, 0, 1, 0, n(A), cosr, sinr)
		}
		if smax < 0.0 {
			r.Row(V, 0)
			blasd.Scale(&r, -1.0)
		}
		if smin < 0.0 {
			r.Row(V, 1)
			blasd.Scale(&r, -1.0)
		}
	}
	S.SetAt(0, math.Abs(smax))
	S.SetAt(1, math.Abs(smin))
	return
}