You are an expert in writing code in many different languages. Your goal is to perform code completion for the following code, keeping in mind the rest of the code and the file meta data.
    Metadata details:
    filename: dlaqps.go
    path of file: ./repos/hoverfly/vendor/gonum.org/v1/gonum/lapack/gonum
    
    the code of the file until where you have to start completion: // Copyright ©2017 The Gonum Authors. All rights reserved.
// Use of this source code is governed by a BSD-style
// license that can be found in the LICENSE file.

package gonum

import (
	"math"

	"gonum.org/v1/gonum/blas"
	"gonum.org/v1/gonum/blas/blas64"
)

// Dlaqps computes a step of QR factorization with column pivoting
// of an m×n matrix A by using Blas-3. It tries to factorize nb
// columns from A starting from the row offset, and updates all
// of the matrix with Dgemm.
//
// In some cases, due to catastrophic cancellations, it cannot
// factorize nb columns. Hence, the actual number of factorized
// columns is returned in kb.
//
// Dlaqps computes a QR factorization with column pivoting of the
// block A[offset:m, 0:nb] of the m×n matrix A. The block
// A[0:offset, 0:n] is accordingly pivoted, but not factorized.
//
// On exit, the upper triangle of block A[offset:m, 0:kb] is the
// triangular factor obtained. The elements in block A[offset:m, 0:n]
// below the diagonal, together with tau, represent the orthogonal
// matrix Q as a product of elementary reflectors.
//
// offset is number of rows of the matrix A that must be pivoted but
// not factorized. offset must not be negative otherwise Dlaqps will panic.
//
// On exit, jpvt holds the permutation that was applied; the jth column
// of A*P was the jpvt[j] column of A. jpvt must have length n,
// otherwise Dlapqs will panic.
//
// On exit tau holds the scalar factors of the elementary reflectors.
// It must have length nb, otherwise Dlapqs will panic.
//
// vn1 and vn2 hold the partial and complete column norms respectively.
// They must have length n, otherwise Dlapqs will panic.
//
// auxv must have length nb, otherwise Dlaqps will panic.
//
// f and ldf represent an n×nb matrix F that is overwritten during the
// call.
//
// Dlaqps is an internal routine. It is exported for testing purposes.
func (impl Implementation) Dlaqps(m, n, offset, nb int, a []float64, lda int, jpvt []int, tau, vn1, vn2, auxv, f []float64, ldf int) (kb int) {
	switch {
	case m < 0:
		panic(mLT0)
	case n < 0:
		panic(nLT0)
	case offset < 0:
		panic(offsetLT0)
	case offset > m:
		panic(offsetGTM)
	case nb < 0:
		panic(nbLT0)
	case nb > n:
		panic(nbGTN)
	case lda < max(1, n):
		panic(badLdA)
	case ldf < max(1, nb):
		panic(badLdF)
	}

	if m == 0 || n == 0 {
		return 0
	}

	switch {
	case len(a) < (m-1)*lda+n:
		panic(shortA)
	case len(jpvt) != n:
		panic(badLenJpvt)
	case len(vn1) < n:
		panic(shortVn1)
	case len(vn2) < n:
		panic(shortVn2)
	}

	if nb == 0 {
		return 0
	}

	switch {
	case len(tau) < nb:
		panic(shortTau)
	case len(auxv) < nb:
		panic(shortAuxv)
	case len(f) < (n-1)*ldf+nb:
		
    
func (impl Implementation) Dlaqps(m, n, offset, nb int, a []float64, lda int, jpvt []int, tau, vn1, vn2, auxv, f []float64, ldf int) (kb int) {
	switch {
	case m < 0:
		panic(mLT0)
	case n < 0:
		panic(nLT0)
	case offset < 0:
		panic(offsetLT0)
	case offset > m:
		panic(offsetGTM)
	case nb < 0:
		panic(nbLT0)
	case nb > n:
		panic(nbGTN)
	case lda < max(1, n):
		panic(badLdA)
	case ldf < max(1, nb):
		panic(badLdF)
	}

	if m == 0 || n == 0 {
		return 0
	}

	switch {
	case len(a) < (m-1)*lda+n:
		panic(shortA)
	case len(jpvt) != n:
		panic(badLenJpvt)
	case len(vn1) < n:
		panic(shortVn1)
	case len(vn2) < n:
		panic(shortVn2)
	}

	if nb == 0 {
		return 0
	}

	switch {
	case len(tau) < nb:
		panic(shortTau)
	case len(auxv) < nb:
		panic(shortAuxv)
	case len(f) < (n-1)*ldf+nb:
		panic(shortF)
	}

	if offset == m {
		return 0
	}

	lastrk := min(m, n+offset)
	lsticc := -1
	tol3z := math.Sqrt(dlamchE)

	bi := blas64.Implementation()

	var k, rk int
	for ; k < nb && lsticc == -1; k++ {
		rk = offset + k

		// Determine kth pivot column and swap if necessary.
		p := k + bi.Idamax(n-k, vn1[k:], 1)
		if p != k {
			bi.Dswap(m, a[p:], lda, a[k:], lda)
			bi.Dswap(k, f[p*ldf:], 1, f[k*ldf:], 1)
			jpvt[p], jpvt[k] = jpvt[k], jpvt[p]
			vn1[p] = vn1[k]
			vn2[p] = vn2[k]
		}

		// Apply previous Householder reflectors to column K:
		//
		// A[rk:m, k] = A[rk:m, k] - A[rk:m, 0:k-1]*F[k, 0:k-1]ᵀ.
		if k > 0 {
			bi.Dgemv(blas.NoTrans, m-rk, k, -1,
				a[rk*lda:], lda,
				f[k*ldf:], 1,
				1,
				a[rk*lda+k:], lda)
		}

		// Generate elementary reflector H_k.
		if rk < m-1 {
			a[rk*lda+k], tau[k] = impl.Dlarfg(m-rk, a[rk*lda+k], a[(rk+1)*lda+k:], lda)
		} else {
			tau[k] = 0
		}

		akk := a[rk*lda+k]
		a[rk*lda+k] = 1

		// Compute kth column of F:
		//
		// Compute F[k+1:n, k] = tau[k]*A[rk:m, k+1:n]ᵀ*A[rk:m, k].
		if k < n-1 {
			bi.Dgemv(blas.Trans, m-rk, n-k-1, tau[k],
				a[rk*lda+k+1:], lda,
				a[rk*lda+k:], lda,
				0,
				f[(k+1)*ldf+k:], ldf)
		}

		// Padding F[0:k, k] with zeros.
		for j := 0; j < k; j++ {
			f[j*ldf+k] = 0
		}

		// Incremental updating of F:
		//
		// F[0:n, k] := F[0:n, k] - tau[k]*F[0:n, 0:k-1]*A[rk:m, 0:k-1]ᵀ*A[rk:m,k].
		if k > 0 {
			bi.Dgemv(blas.Trans, m-rk, k, -tau[k],
				a[rk*lda:], lda,
				a[rk*lda+k:], lda,
				0,
				auxv, 1)
			bi.Dgemv(blas.NoTrans, n, k, 1,
				f, ldf,
				auxv, 1,
				1,
				f[k:], ldf)
		}

		// Update the current row of A:
		//
		// A[rk, k+1:n] = A[rk, k+1:n] - A[rk, 0:k]*F[k+1:n, 0:k]ᵀ.
		if k < n-1 {
			bi.Dgemv(blas.NoTrans, n-k-1, k+1, -1,
				f[(k+1)*ldf:], ldf,
				a[rk*lda:], 1,
				1,
				a[rk*lda+k+1:], 1)
		}

		// Update partial column norms.
		if rk < lastrk-1 {
			for j := k + 1; j < n; j++ {
				if vn1[j] == 0 {
					continue
				}

				// The following marked lines follow from the
				// analysis in Lapack Working Note 176.
				r := math.Abs(a[rk*lda+j]) / vn1[j] // *
				temp := math.Max(0, 1-r*r)          // *
				r = vn1[j] / vn2[j]                 // *
				temp2 := temp * r * r               // *
				if temp2 < tol3z {
					// vn2 is used here as a collection of
					// indices into vn2 and also a collection
					// of column norms.
					vn2[j] = float64(lsticc)
					lsticc = j
				} else {
					vn1[j] *= math.Sqrt(temp) // *
				}
			}
		}

		a[rk*lda+k] = akk
	}
	kb = k
	rk = offset + kb

	// Apply the block reflector to the rest of the matrix:
	//
	// A[offset+kb+1:m, kb+1:n] := A[offset+kb+1:m, kb+1:n] - A[offset+kb+1:m, 1:kb]*F[kb+1:n, 1:kb]ᵀ.
	if kb < min(n, m-offset) {
		bi.Dgemm(blas.NoTrans, blas.Trans,
			m-rk, n-kb, kb, -1,
			a[rk*lda:], lda,
			f[kb*ldf:], ldf,
			1,
			a[rk*lda+kb:], lda)
	}

	// Recomputation of difficult columns.
	for lsticc >= 0 {
		itemp := int(vn2[lsticc])

		// NOTE: The computation of vn1[lsticc] relies on the fact that
		// Dnrm2 does not fail on vectors with norm below the value of
		// sqrt(dlamchS)
		v := bi.Dnrm2(m-rk, a[rk*lda+lsticc:], lda)
		vn1[lsticc] = v
		vn2[lsticc] = v

		lsticc = itemp
	}

	return kb
}