lapack/testlapack: add implementation comments to Dlartg test

lapack/testlapack/dlartg.go

@@ -20,24 +20,32 @@ type Dlartger interface {
 func DlartgTest(t *testing.T, impl Dlartger) {
 	const tol = 1e-14
 	// safmn2 and safmx2 are copied from native.Dlartg.
+	// safmn2 ~ 2*10^{-146}
 	safmn2 := math.Pow(dlamchB, math.Trunc(math.Log(dlamchS/dlamchE)/math.Log(dlamchB)/2))
+	// safmx2 ~ 5*10^145
 	safmx2 := 1 / safmn2
 	rnd := rand.New(rand.NewSource(1))
 	for i := 0; i < 1000; i++ {
+		// Generate randomly huge, tiny, and "normal" input arguments to Dlartg.
 		var f float64
 		var fHuge bool
 		switch rnd.Intn(3) {
 		case 0:
 			// Huge f.
 			fHuge = true
-			f = math.Pow(10, 10-20*rnd.Float64()) * safmx2
+			// scale is in the range (10^{-10}, 10^10].
+			scale := math.Pow(10, 10-20*rnd.Float64())
+			// f is in the range (5*10^135, 5*10^155].
+			f = scale * safmx2
 		case 1:
 			// Tiny f.
+			// f is in the range (2*10^{-156}, 2*10^{-136}].
 			f = math.Pow(10, 10-20*rnd.Float64()) * safmn2
 		default:
 			f = rnd.NormFloat64()
 		}
 		if rnd.Intn(2) == 0 {
+			// Sometimes change the sign of f.
 			f *= -1
 		}
 
@@ -58,26 +66,38 @@ func DlartgTest(t *testing.T, impl Dlartger) {
 			g *= -1
 		}
 
+		// Generate a plane rotation so that
+		//  [ cs sn] * [f] = [r]
+		//  [-sn cs]   [g] = [0]
 		cs, sn, r := impl.Dlartg(f, g)
 
+		// Check that the first equation holds.
 		rWant := cs*f + sn*g
 		if !floats.EqualWithinAbsOrRel(math.Abs(rWant), math.Abs(r), tol, tol) {
 			t.Errorf("Case f=%v,g=%v: unexpected r. Want %v, got %v", f, g, rWant, r)
 		}
+		// Check that cs and sn define a plane rotation. The 2×2 matrix
+		// has orthogonal columns by construction, so only check that
+		// the columns/rows have unit norm.
 		oneTest := cs*cs + sn*sn
 		if math.Abs(oneTest-1) > tol {
 			t.Errorf("Case f=%v,g=%v: expected cs^2+sn^2==1, got %v", f, g, oneTest)
 		}
 		if !fHuge && !gHuge {
+			// Check that the second equation holds.
+			// If both numbers are huge, cancellation errors make
+			// this test unreliable.
 			zeroTest := -sn*f + cs*g
 			if math.Abs(zeroTest) > tol {
 				t.Errorf("Case f=%v,g=%v: expected zero, got %v", f, g, zeroTest)
 			}
 		}
+		// Check that cs is positive as documented.
 		if math.Abs(f) > math.Abs(g) && cs < 0 {
 			t.Errorf("Case f=%v,g=%v: unexpected negative cs %v", f, g, cs)
 		}
 	}
+	// Check other documented special cases.
 	for i := 0; i < 100; i++ {
 		cs, sn, _ := impl.Dlartg(rnd.NormFloat64(), 0)
 		if cs != 1 {