Golang FloatMatrix.Rows方法代码示例

func runTest(A *matrix.FloatMatrix, ntest, LB int) time.Duration {

	var flags matops.Flags
	var mintime time.Duration

	N := A.Rows()
	ipiv := make([]int, N, N)
	flags = matops.LOWER
	if testUpper {
		flags = matops.UPPER
	}

	W := matrix.FloatZeros(A.Rows(), LB+2)
	fnc := func() {
		_, ERRmatops = matops.DecomposeBK(A, W, ipiv, flags, LB)
	}

	A0 := A.Copy()
	for n := 0; n < ntest; n++ {
		if n > 0 {
			// restore original A
			A0.CopyTo(A)
		}
		mperf.FlushCache()
		time0 := mperf.Timeit(fnc)
		if n == 0 || time0 < mintime {
			mintime = time0
		}
		if verbose {
			fmt.Printf("%.4f ms\n", time0.Seconds()*1000.0)
		}
	}
	return mintime
}

/*
 *  ( A11  a12 )   ( U11 u12 )( D1  0  )( U11.t 0 )
 *  ( a21  a22 )   (  0   1  )(  0  d2 )( u12.t 1 )
 *
 *   a22  =   d2
 *   a01  =   u12*d2                       => u12 = a12/d2
 *   A11  =   u12*d2*u12.t + U11*D1*U11.t  => U11 = A11 - u12*d2*u12.t
 */
func unblkUpperLDLnoPiv(A *matrix.FloatMatrix) (err error) {
	var ATL, ATR, ABL, ABR matrix.FloatMatrix
	var A00, a01, A02, a11, a12, A22 matrix.FloatMatrix

	err = nil
	partition2x2(
		&ATL, &ATR,
		&ABL, &ABR, A, 0, 0, pBOTTOMRIGHT)

	for ATL.Rows() > 0 {
		repartition2x2to3x3(&ATL,
			&A00, &a01, &A02,
			nil, &a11, &a12,
			nil, nil, &A22, A, 1, pTOPLEFT)

		// --------------------------------------------------------

		// A00 = A00 - u01*d11*u01.T = A00 - a01*a01.T/a11; triangular update
		err = MVUpdateTrm(&A00, &a01, &a01, -1.0/a11.Float(), UPPER)

		// u01 = a01/a11
		InvScale(&a01, a11.Float())
		// ---------------------------------------------------------

		continue3x3to2x2(
			&ATL, &ATR,
			&ABL, &ABR, &A00, &a11, &A22, A, pTOPLEFT)
	}
	return
}

func runRefTest(A *matrix.FloatMatrix, ntest, LB int) time.Duration {

	var mintime time.Duration

	N := A.Rows()
	ipiv := make([]int32, N, N)
	lopt := linalg.OptLower
	if testUpper {
		lopt = linalg.OptUpper
	}

	fnc := func() {
		ERRlapack = lapack.Sytrf(A, ipiv, lopt)
	}

	A0 := A.Copy()
	for n := 0; n < ntest; n++ {
		if n > 0 {
			// restore original A
			A0.CopyTo(A)
		}
		mperf.FlushCache()
		time0 := mperf.Timeit(fnc)
		if n == 0 || time0 < mintime {
			mintime = time0
		}
	}
	return mintime
}

/*
 *  ( a11  a12 )   ( 1   0   )( d1  0   )( l  l21.t )
 *  ( a21  A22 )   ( l21 L22 )(  0  A22 )( 0  L22.t )
 *
 *   a11  =   d1
 *   a21  =   l21*d1                       => l21 = a21/d1
 *   A22  =   l21*d1*l21.t + L22*D2*L22.t  => L22 = A22 - l21*d1*l21t
 */
func unblkLowerLDLnoPiv(A *matrix.FloatMatrix) (err error) {
	var ATL, ATR, ABL, ABR matrix.FloatMatrix
	var A00, a10, a11, A20, a21, A22 matrix.FloatMatrix

	err = nil
	partition2x2(
		&ATL, &ATR,
		&ABL, &ABR, A, 0, 0, pTOPLEFT)

	for ATL.Rows() < A.Rows() {
		repartition2x2to3x3(&ATL,
			&A00, nil, nil,
			&a10, &a11, nil,
			&A20, &a21, &A22, A, 1, pBOTTOMRIGHT)

		// --------------------------------------------------------

		// d11 = a11; no-op

		// A22 = A22 - l21*d11*l21.T = A22 - a21*a21.T/a11; triangular update
		err = MVUpdateTrm(&A22, &a21, &a21, -1.0/a11.Float(), LOWER)

		// l21 = a21/a11
		InvScale(&a21, a11.Float())
		// ---------------------------------------------------------

		continue3x3to2x2(
			&ATL, &ATR,
			&ABL, &ABR, &A00, &a11, &A22, A, pBOTTOMRIGHT)
	}
	return
}

// Calculate C = alpha*A*B.T + beta*C, C is M*N, A is M*P and B is N*P
func MMMultTransB(C, A, B *matrix.FloatMatrix, alpha, beta float64) error {
	psize := int64(C.NumElements() * A.Cols())
	Ar := A.FloatArray()
	ldA := A.LeadingIndex()
	Br := B.FloatArray()
	ldB := B.LeadingIndex()
	Cr := C.FloatArray()
	ldC := C.LeadingIndex()
	if nWorker <= 1 || psize <= limitOne {
		calgo.DMult(Cr, Ar, Br, alpha, beta, calgo.TRANSB, ldC, ldA, ldB,
			B.Rows(), 0, C.Cols(), 0, C.Rows(), vpLen, nB, mB)
		return nil
	}

	// here we have more than one worker available
	worker := func(cstart, cend, rstart, rend int, ready chan int) {
		calgo.DMult(Cr, Ar, Br, alpha, beta, calgo.TRANSB, ldC, ldA, ldB, B.Rows(),
			cstart, cend, rstart, rend, vpLen, nB, mB)
		ready <- 1
	}
	colworks, rowworks := divideWork(C.Rows(), C.Cols(), nWorker)
	scheduleWork(colworks, rowworks, C.Cols(), C.Rows(), worker)
	//scheduleWork(colworks, rowworks, worker)
	return nil
}

func main() {
	blas.PanicOnError(true)
	matrix.PanicOnError(true)

	var data *matrix.FloatMatrix = nil
	flag.Parse()
	if len(dataVal) > 0 {
		data, _ = matrix.FloatParse(dataVal)
		if data == nil {
			fmt.Printf("could not parse:\n%s\n", dataVal)
			return
		}
	} else {
		data = matrix.FloatNormal(20, 20)
	}

	if len(spPath) > 0 {
		checkpnt.Reset(spPath)
		checkpnt.Activate()
		checkpnt.Verbose(spVerbose)
		checkpnt.Format("%.7f")
	}

	sol, err := mcsdp(data)
	if sol != nil && sol.Status == cvx.Optimal {
		x := sol.Result.At("x")[0]
		z := sol.Result.At("z")[0]
		matrix.Reshape(z, data.Rows(), data.Rows())
		fmt.Printf("x=\n%v\n", x.ToString("%.9f"))
		//fmt.Printf("z=\n%v\n", z.ToString("%.9f"))
		check(x, z)
	} else {
		fmt.Printf("status: %v\n", err)
	}
}

func updateBlas(t *testing.T, Y1, Y2, C1, C2, T, W *matrix.FloatMatrix) {
	if W.Rows() != C1.Cols() {
		panic("W.Rows != C1.Cols")
	}
	// W = C1.T
	ScalePlus(W, C1, 0.0, 1.0, TRANSB)
	//fmt.Printf("W = C1.T:\n%v\n", W)
	// W = C1.T*Y1
	blas.TrmmFloat(Y1, W, 1.0, linalg.OptLower, linalg.OptUnit, linalg.OptRight)
	t.Logf("W = C1.T*Y1:\n%v\n", W)
	// W = W + C2.T*Y2
	blas.GemmFloat(C2, Y2, W, 1.0, 1.0, linalg.OptTransA)
	t.Logf("W = W + C2.T*Y2:\n%v\n", W)

	// --- here: W == C.T*Y ---
	// W = W*T
	blas.TrmmFloat(T, W, 1.0, linalg.OptUpper, linalg.OptRight)
	t.Logf("W = C.T*Y*T:\n%v\n", W)

	// --- here: W == C.T*Y*T ---
	// C2 = C2 - Y2*W.T
	blas.GemmFloat(Y2, W, C2, -1, 1.0, linalg.OptTransB)
	t.Logf("C2 = C2 - Y2*W.T:\n%v\n", C2)
	//  W = Y1*W.T ==> W.T = W*Y1.T
	blas.TrmmFloat(Y1, W, 1.0, linalg.OptLower,
		linalg.OptUnit, linalg.OptRight, linalg.OptTrans)
	t.Logf("W.T = W*Y1.T:\n%v\n", W)

	// C1 = C1 - W.T
	ScalePlus(C1, W, 1.0, -1.0, TRANSB)
	//fmt.Printf("C1 = C1 - W.T:\n%v\n", C1)

	// --- here: C = (I - Y*T*Y.T).T * C ---
}

func runRefTest(A *matrix.FloatMatrix, ntest, LB int) time.Duration {

	var flags matops.Flags
	var mintime time.Duration

	N := A.Rows()
	ipiv := make([]int, N, N)
	flags = matops.LOWER
	if testUpper {
		flags = matops.UPPER
	}

	W := matrix.FloatZeros(A.Rows(), LB+2)
	fnc := func() {
		_, ERRref = matops.DecomposeLDL(A, W, ipiv, flags, 0)
	}

	A0 := A.Copy()
	for n := 0; n < ntest; n++ {
		if n > 0 {
			// restore original A
			A0.CopyTo(A)
		}
		mperf.FlushCache()
		time0 := mperf.Timeit(fnc)
		if n == 0 || time0 < mintime {
			mintime = time0
		}
	}
	return mintime
}

// Inverse UNIT diagonal tridiagonal matrix
func unblockedInverseUnitLower(A *matrix.FloatMatrix) (err error) {
	var ATL, ATR, ABL, ABR matrix.FloatMatrix
	var A00, a10t, a11, A20, a21, A22 matrix.FloatMatrix

	err = nil
	partition2x2(
		&ATL, &ATR,
		&ABL, &ABR, A, 0, 0, pTOPLEFT)

	for ATL.Rows() < A.Rows() {
		repartition2x2to3x3(&ATL,
			&A00, nil, nil,
			&a10t, &a11, nil,
			&A20, &a21, &A22, A, 1, pBOTTOMRIGHT)
		// -------------------------------------------------

		// a21 = -a21
		Scale(&a21, -1.0)
		// A20 = A20 + a21*a10.t
		MVRankUpdate(&A20, &a21, &a10t, 1.0)

		// -------------------------------------------------
		continue3x3to2x2(
			&ATL, &ATR,
			&ABL, &ABR, &A00, &a11, &A22, A, pBOTTOMRIGHT)
	}
	return
}

// Inverse NON-UNIT diagonal tridiagonal matrix
func unblockedInverseUpper(A *matrix.FloatMatrix) (err error) {
	var ATL, ATR, ABL, ABR matrix.FloatMatrix
	var A00, a01, A02, a11, a12t, A22 matrix.FloatMatrix

	err = nil
	partition2x2(
		&ATL, &ATR,
		&ABL, &ABR, A, 0, 0, pTOPLEFT)

	for ATL.Rows() < A.Rows() {
		repartition2x2to3x3(&ATL,
			&A00, &a01, &A02,
			nil, &a11, &a12t,
			nil, nil, &A22, A, 1, pBOTTOMRIGHT)
		// -------------------------------------------------
		aval := a11.Float()

		// a12 = -a12/a11
		InvScale(&a12t, -aval)
		// A02 = A02 + a01*a12
		MVRankUpdate(&A02, &a01, &a12t, 1.0)
		// a01 = a01/a11
		InvScale(&a01, aval)
		// a11 = 1.0/a11
		a11.SetAt(0, 0, 1.0/aval)

		// -------------------------------------------------
		continue3x3to2x2(
			&ATL, &ATR,
			&ABL, &ABR, &A00, &a11, &A22, A, pBOTTOMRIGHT)
	}
	return
}

// unblocked LU decomposition w/o pivots, FLAME LU nopivots variant 5
func unblockedLUnoPiv(A *matrix.FloatMatrix) (err error) {
	var ATL, ATR, ABL, ABR matrix.FloatMatrix
	var A00, a01, A02, a10, a11, a12, A20, a21, A22 matrix.FloatMatrix

	err = nil
	partition2x2(
		&ATL, &ATR,
		&ABL, &ABR, A, 0, 0, pTOPLEFT)

	for ATL.Rows() < A.Rows() {
		repartition2x2to3x3(&ATL,
			&A00, &a01, &A02,
			&a10, &a11, &a12,
			&A20, &a21, &A22, A, 1, pBOTTOMRIGHT)

		// a21 = a21/a11
		//a21.Scale(1.0/a11.Float())
		InvScale(&a21, a11.Float())
		// A22 = A22 - a21*a12
		err = MVRankUpdate(&A22, &a21, &a12, -1.0)

		continue3x3to2x2(
			&ATL, &ATR,
			&ABL, &ABR, &A00, &a11, &A22, A, pBOTTOMRIGHT)
	}
	return
}

func runTest(A *matrix.FloatMatrix, ntest, LB int) time.Duration {
	var W *matrix.FloatMatrix = nil
	var mintime time.Duration

	N := A.Cols()
	tau := matrix.FloatZeros(N, 1)
	if LB > 0 {
		W = matrix.FloatZeros(A.Rows(), LB)
	}
	fnc := func() {
		_, ERRmatops = matops.DecomposeQR(A, tau, W, LB)
	}

	A0 := A.Copy()
	for n := 0; n < ntest; n++ {
		if n > 0 {
			// restore original A
			A0.CopyTo(A)
			tau.Scale(0.0)
		}
		mperf.FlushCache()
		time0 := mperf.Timeit(fnc)
		if n == 0 || time0 < mintime {
			mintime = time0
		}
		if verbose {
			fmt.Printf("%.4f ms\n", time0.Seconds()*1000.0)
		}
	}
	return mintime
}

// blocked LU decomposition w/o pivots, FLAME LU nopivots variant 5
func blockedLUnoPiv(A *matrix.FloatMatrix, nb int) (err error) {
	var ATL, ATR, ABL, ABR matrix.FloatMatrix
	var A00, A01, A02, A10, A11, A12, A20, A21, A22 matrix.FloatMatrix

	err = nil
	partition2x2(
		&ATL, &ATR,
		&ABL, &ABR, A, 0, 0, pTOPLEFT)

	for ATL.Rows() < A.Rows() {
		repartition2x2to3x3(&ATL,
			&A00, &A01, &A02,
			&A10, &A11, &A12,
			&A20, &A21, &A22, A, nb, pBOTTOMRIGHT)

		// A00 = LU(A00)
		unblockedLUnoPiv(&A11)
		// A12 = trilu(A00)*A12.-1  (TRSM)
		SolveTrm(&A12, &A11, 1.0, LEFT|LOWER|UNIT)
		// A21 = A21.-1*triu(A00) (TRSM)
		SolveTrm(&A21, &A11, 1.0, RIGHT|UPPER)
		// A22 = A22 - A21*A12
		Mult(&A22, &A21, &A12, -1.0, 1.0, NOTRANS)

		continue3x3to2x2(
			&ATL, &ATR,
			&ABL, &ABR, &A00, &A11, &A22, A, pBOTTOMRIGHT)
	}
	return
}

func runTest(A *matrix.FloatMatrix, ntest, LB int) time.Duration {

	var mintime time.Duration

	M := A.Rows()
	N := A.Cols()
	nN := N
	if M < N {
		nN = M
	}
	ipiv := make([]int, nN, nN)

	fnc := func() {
		_, ERRmatops = matops.DecomposeLU(A, ipiv, LB)
	}

	A0 := A.Copy()
	for n := 0; n < ntest; n++ {
		if n > 0 {
			// restore original A
			A0.CopyTo(A)
		}
		mperf.FlushCache()
		time0 := mperf.Timeit(fnc)
		if n == 0 || time0 < mintime {
			mintime = time0
		}
		if verbose {
			fmt.Printf("%.4f ms\n", time0.Seconds()*1000.0)
		}
	}
	return mintime
}

func runRefTest(A *matrix.FloatMatrix, ntest, LB int) time.Duration {

	var mintime time.Duration

	M := A.Rows()
	N := A.Cols()
	nN := N
	if M < N {
		nN = M
	}
	ipiv := make([]int32, nN, nN)

	fnc := func() {
		ERRlapack = lapack.Getrf(A, ipiv)
	}

	A0 := A.Copy()
	for n := 0; n < ntest; n++ {
		if n > 0 {
			// restore original A
			A0.CopyTo(A)
		}
		mperf.FlushCache()
		time0 := mperf.Timeit(fnc)
		if n == 0 || time0 < mintime {
			mintime = time0
		}
	}
	return mintime
}