Golang math.Sqrt函数代码示例

/*
   Copy x to y using packed storage.

   The vector x is an element of S, with the 's' components stored in
   unpacked storage.  On return, x is copied to y with the 's' components
   stored in packed storage and the off-diagonal entries scaled by
   sqrt(2).
*/
func pack(x, y *matrix.FloatMatrix, dims *DimensionSet, opts ...la_.Option) (err error) {
	/*DEBUGGED*/
	err = nil
	mnl := la_.GetIntOpt("mnl", 0, opts...)
	offsetx := la_.GetIntOpt("offsetx", 0, opts...)
	offsety := la_.GetIntOpt("offsety", 0, opts...)

	nlq := mnl + dims.At("l")[0] + dims.Sum("q")
	blas.Copy(x, y, &la_.IOpt{"n", nlq}, &la_.IOpt{"offsetx", offsetx},
		&la_.IOpt{"offsety", offsety})

	iu, ip := offsetx+nlq, offsety+nlq
	for _, n := range dims.At("s") {
		for k := 0; k < n; k++ {
			blas.Copy(x, y, &la_.IOpt{"n", n - k}, &la_.IOpt{"offsetx", iu + k*(n+1)},
				&la_.IOpt{"offsety", ip})
			y.SetIndex(ip, (y.GetIndex(ip) / math.Sqrt(2.0)))
			ip += n - k
		}
		iu += n * n
	}
	np := dims.SumPacked("s")
	blas.ScalFloat(y, math.Sqrt(2.0), &la_.IOpt{"n", np}, &la_.IOpt{"offset", offsety + nlq})
	return
}

// Computes the error estimate based on fNew:
func (p *peer) computeErrorModel(in *integration) (errorEstimate float64) {
	var i_n, j_stg uint

	// Loop: EmFactors
	for i_n = 0; i_n < in.n; i_n++ {
		var factor float64 = 0.0
		for j_stg = 0; j_stg < p.Stages; j_stg++ {
			factor += p.errorModelWeights[j_stg] * in.fNew[j_stg][i_n]
		}
		in.errorFactors[i_n] = math.Pow(factor/(in.AbsoluteTolerance+in.RelativeTolerance*math.Abs(in.yOld[p.Stages-1][i_n])), 2.0)
	}

	// compute error quotient/20070803
	// step ratio from error model ((1+a)^p-a^p)/est+a^p)^(1/p)-a, p=order/2:
	errorRelative := 0.0
	for i_n = 0; i_n < in.n; i_n++ {
		errorRelative += in.errorFactors[i_n]
	}

	errorEstimate = in.stepEstimate*math.Sqrt(errorRelative/float64(in.n)) + 1e-8
	errorModelDenom := math.Pow(math.Pow(in.stepRatio, 2.0)+p.errorModelA, float64(p.Order)/2.0) - p.errorModelA0
	errorStepRatio := math.Pow(errorModelDenom/errorEstimate+p.errorModelA0, 2.0/float64(p.Order)) - p.errorModelA
	in.stepEstimate = in.stepPrevious * math.Max(in.stepRatioMin, math.Min(0.95*math.Sqrt(errorStepRatio), p.stepRatioMax)) // safety interval

	return
}

// Convert the given coordinate from X,Y to polar in the given PolarSystem
func (coord Coordinate) ToPolar(system PolarSystem) (polarCoord PolarCoordinate) {

	coord.X += system.XOffset
	coord.Y += system.YOffset

	// clip coordinates to system's area
	if coord.X < system.XMin {
		fmt.Println("WARNING: X value was outside left bounds, clipping", coord.X, "to", system.XMin)
		coord.X = system.XMin
	}
	if coord.X > system.XMax {
		fmt.Println("WARNING: X value was outside right bounds, clipping", coord.X, "to", system.XMax)
		coord.X = system.XMax
	}
	if coord.Y < system.YMin {
		fmt.Println("WARNING: Y value was outside top bounds, clipping", coord.Y, "to", system.YMin)
		coord.Y = system.YMin
	}
	if coord.Y > system.YMax {
		fmt.Println("WARNING: Y value was outside bottom bounds, clipping", coord.Y, "to", system.YMax)
		coord.Y = system.YMax
	}

	polarCoord.LeftDist = math.Sqrt(coord.X*coord.X + coord.Y*coord.Y)
	xDiff := system.RightMotorDist - coord.X
	polarCoord.RightDist = math.Sqrt(xDiff*xDiff + coord.Y*coord.Y)
	polarCoord.PenUp = coord.PenUp
	return
}

func ellipse(dst *sdl.Surface, org Point, a, b int) {
	colour := image.RGBAColor{0, 255, 0, 255}

	// set the poles
	dst.Set(int(org.X), int(org.Y)+b, colour)
	dst.Set(int(org.X), int(org.Y)-b, colour)
	dst.Set(int(org.X)+a, int(org.Y), colour)
	dst.Set(int(org.X)-a, int(org.Y), colour)

	// add a bias to a and b for smoother poles
	ba := float64(a) + 0.5
	bb := float64(b) + 0.5

	// draw and rotate a quadrant
	for x := 1; x < a; x++ {
		y1 := (-bb * float64(x)) / (ba * math.Sqrt(ba*ba-float64(x*x)))
		y := int(math.Sqrt((bb * bb) * (1.0 - float64(x*x)/(ba*ba))))
		if y1 > -1.0 {
			dst.Set(int(org.X)+x, int(org.Y)+y, colour)
			dst.Set(int(org.X)-x, int(org.Y)+y, colour)
			dst.Set(int(org.X)+x, int(org.Y)-y, colour)
			dst.Set(int(org.X)-x, int(org.Y)-y, colour)
		} else {
			for dy := 1; dy <= y; dy++ {
				dx := int(math.Sqrt((ba * ba) * (1.0 - float64(dy*dy)/(bb*bb))))
				dst.Set(int(org.X)+dx, int(org.Y)+dy, colour)
				dst.Set(int(org.X)-dx, int(org.Y)+dy, colour)
				dst.Set(int(org.X)+dx, int(org.Y)-dy, colour)
				dst.Set(int(org.X)-dx, int(org.Y)-dy, colour)
			}
			break
		}
	}
}

/*
 * Compute Givens rotation such that
 *
 *   G(s,c)*v = (r)   ==  (  c s ).T ( a ) = ( r )
 *              (0)       ( -s c )   ( b )   ( 0 )
 *
 * and
 *
 *   v*G(s,c) = (r 0 ) == (a b ) (  c  s ) = ( r 0 )
 *                               ( -s  c )
 *
 */
func ComputeGivens(a, b float64) (c float64, s float64, r float64) {

	if b == 0.0 {
		c = 1.0
		s = 0.0
		r = a
	} else if a == 0.0 {
		c = 0.0
		s = 1.0
		r = b
	} else if math.Abs(b) > math.Abs(a) {
		t := a / b
		u := math.Sqrt(1.0 + t*t)
		if math.Signbit(b) {
			u = -u
		}
		s = 1.0 / u
		c = s * t
		r = b * u
	} else {
		t := b / a
		u := math.Sqrt(1.0 + t*t)
		r = a * u
		c = 1.0 / u
		s = c * t
	}
	return
}

// Correlation describes the degree of relationship between two sets of data
func Correlation(data1, data2 Float64Data) (float64, error) {

	l1 := data1.Len()
	l2 := data2.Len()

	if l1 == 0 || l2 == 0 {
		return 0, errors.New("Input data must not be empty")
	}

	if l1 != l2 {
		return 0, errors.New("Input data must be same length")
	}

	var sum_xsq, sum_ysq, sum_cross float64

	mean_x := data1.Get(0)
	mean_y := data2.Get(0)

	for i := 1; i < l1; i++ {
		ratio := float64(i) / float64(i+1)
		delta_x := data1.Get(i) - mean_x
		delta_y := data2.Get(i) - mean_y
		sum_xsq += delta_x * delta_x * ratio
		sum_ysq += delta_y * delta_y * ratio
		sum_cross += delta_x * delta_y * ratio
		mean_x += delta_x / float64(i+1)
		mean_y += delta_y / float64(i+1)
	}

	return sum_cross / (math.Sqrt(sum_xsq) * math.Sqrt(sum_ysq)), nil
}

func renderPixel(x int, y int, cx Vec, cy Vec) {
	var r1, r2 float64
	var dx, dy float64
	var radiance Vec
	var direction Vec

	for sy, i := 0, (h-y-1)*w+x; sy < 2; sy++ {
		for sx := 0; sx < 2; sx++ {
			radiance.x = 0
			radiance.y = 0
			radiance.z = 0
			for s := 0; s < samps; s++ {
				r1, r2 = 2*rand.Float64(), 2*rand.Float64()
				if r1 < 1 {
					dx = math.Sqrt(r1) - 1
				} else {
					dx = 1 - math.Sqrt(2-r1)
				}
				if r2 < 1 {
					dy = math.Sqrt(r2) - 1
				} else {
					dy = 1 - math.Sqrt(2-r2)
				}
				direction = Add(Add(SMul(cx, ((float64(sx)*.5+dx)/2+float64(x))/float64(w)-.5),
					SMul(cy, ((float64(sy)+.5+dy)/2+float64(y))/float64(h)-.5)), cam.Direction)
				radiance = Add(radiance, SMul(Radiance(&Ray{Add(cam.Origin, SMul(direction, 140.0)), Norm(direction)}, 0), 1.0/float64(samps)))
			}
			colors[i] = Add(colors[i], SMul(radiance, 0.25))
		}
	}
}

func TestGenSparseVector(t *testing.T) {

	si1 := StringIndex{"a", "b", "c"}
	v1 := []Value{1, 2, 3}

	sv1 := NewGenSparseVector(si1, v1)

	if sv1.Mag() != Value(math.Sqrt(1+4+9)) {
		t.Fatalf("Mag wrong - have %f", sv1.Mag())
	}

	si2 := StringIndex{"b", "d", "a"}
	v2 := []Value{1, 2, 7}
	sv2 := NewGenSparseVector(si2, v2)

	dot := sv1.Dot(sv2)
	if dot != 7+2 {
		t.Fatalf("Dot not as expected. Have %f", dot)
	}

	cos := sv1.Cos(sv2)
	if cos != Value((7+2)/(math.Sqrt(1+4+9)*math.Sqrt(1+4+49))) {
		t.Fatalf("cos wrong. Have %f", cos)
	}

	mean := sv1.Mean()
	if mean != 2 {
		t.Fatalf("Mean not as expected, have %f", mean)
	}

	sv1.AddConst(2)
	if sv1.Mag() != 7.071068 {
		t.Fatalf("Mag (2) not as expected, have %f", sv1.Mag())
	}
	mean = sv1.Mean()
	if mean != 4 {
		t.Fatalf("Mean (2) not as expected. Have %f", mean)
	}

	sv1.SubConst(2)
	mean = sv1.Mean()
	if mean != 2 {
		t.Fatalf("Mean (3) not as expected. Have %f", mean)
	}

	vals := make([]Value, 0, 3)
	sv1.Iter(func(index interface{}, value Value) {
		vals = append(vals, value)
	})
	if !reflect.DeepEqual(vals, v1) {
		t.Fatalf("iter values not as expected. Have %v", vals)
	}

	sv1.IterUpdate(func(index interface{}, value Value) Value {
		return 3
	})
	if sv1.Mean() != 3 {
		t.Fatalf("that didn't work")
	}
}

// Change the aspect ratio of a rectangle.
// The mode can be "area", "width", "height", "fit", "fill" or "stretch".
// Panics if mode is empty or unrecognized.
func SetAspect(w, h, aspect float64, mode string) (float64, float64) {
	switch mode {
	case "area":
		// aspect = width / height
		// width = height * aspect
		// width^2 = width * height * aspect
		// height = width / aspect
		// height^2 = width * height / aspect
		w, h = math.Sqrt(w*h*aspect), math.Sqrt(w*h/aspect)
	case "width":
		// Set height from width.
		h = w / aspect
	case "height":
		// Set width from height.
		w = h * aspect
	case "fit":
		// Shrink one dimension.
		w, h = math.Min(w, h*aspect), math.Min(h, w/aspect)
	case "fill":
		// Grow one dimension.
		w, h = math.Max(w, h*aspect), math.Max(h, w/aspect)
	case "stretch":
		// Do nothing.
	case "":
		panic("no mode specified")
	default:
		panic("unknown mode: " + mode)
	}
	return w, h
}

func testReadUniformity(t *testing.T, n int, seed int64) {
	r := New(NewSource(seed))
	buf := make([]byte, n)
	nRead, err := r.Read(buf)
	if err != nil {
		t.Errorf("Read err %v", err)
	}
	if nRead != n {
		t.Errorf("Read returned unexpected n; %d != %d", nRead, n)
	}

	// Expect a uniform distribution of byte values, which lie in [0, 255].
	var (
		mean       = 255.0 / 2
		stddev     = math.Sqrt(255.0 * 255.0 / 12.0)
		errorScale = stddev / math.Sqrt(float64(n))
	)

	expected := &statsResults{mean, stddev, 0.10 * errorScale, 0.08 * errorScale}

	// Cast bytes as floats to use the common distribution-validity checks.
	samples := make([]float64, n)
	for i, val := range buf {
		samples[i] = float64(val)
	}
	// Make sure that the entire set matches the expected distribution.
	checkSampleDistribution(t, samples, expected)
}

func calculateConfidenceInterval2(nt, nc, mt, mc, sdt, sdc float64) confInterval2 {

	var ci confInterval2
	//var z = 1.96 // http://www.dummies.com/how-to/content/creating-a-confidence-interval-for-the-difference-.html
	//var z = 2.58 // http://www.dummies.com/how-to/content/creating-a-confidence-interval-for-the-difference-.html

	ci.t1min = mt - zScore*(sdt/math.Sqrt(nt))
	ci.t0min = mc - zScore*(sdt/math.Sqrt(nc))

	ci.t1max = mt + zScore*(sdt/math.Sqrt(nt))
	ci.t0max = mt + zScore*(sdt/math.Sqrt(nc))

	ci.overlap = false

	if ci.t1max <= ci.t0max && ci.t1max >= ci.t0min {
		ci.overlap = true
	}

	if ci.t1min <= ci.t0max && ci.t1min >= ci.t0min {
		ci.overlap = true
	}

	// check for encirclement
	if ci.t0min <= ci.t1max && ci.t0min >= ci.t1min {
		ci.overlap = true
	}

	return ci
}

func Correlation(data1 interface{}, stride1 int, data2 interface{}, stride2 int, n int) (r float64) {
	var ratio, delta_x, delta_y, sum_xsq, sum_ysq, sum_cross, mean_x, mean_y float64

	arry1 := reflect.ValueOf(data1)
	arry2 := reflect.ValueOf(data2)

	/*
	 Compute:
	 sum_xsq = Sum [ (x_i - mu_x)^2 ],
	 sum_ysq = Sum [ (y_i - mu_y)^2 ] and
	 sum_cross = Sum [ (x_i - mu_x) * (y_i - mu_y) ]
	 using the above relation from Welford's paper
	*/

	mean_x = number.Float(arry1.Index(0 * stride1))
	mean_y = number.Float(arry2.Index(0 * stride2))

	for i := 0; i < n; i++ {
		v1 := number.Float(arry1.Index(i * stride1))
		v2 := number.Float(arry2.Index(i * stride2))
		ratio = float64(i) / float64(i+1)
		delta_x = v1 - mean_x
		delta_y = v2 - mean_y
		sum_xsq += delta_x * delta_x * ratio
		sum_ysq += delta_y * delta_y * ratio
		sum_cross += delta_x * delta_y * ratio
		mean_x += delta_x / float64(i+1)
		mean_y += delta_y / float64(i+1)
	}

	r = sum_cross / (math.Sqrt(sum_xsq) * math.Sqrt(sum_ysq))

	return
}

// Randomize randomizes the weights and biases
// such that the sum of the weights has a mean
// of 0 and a variance of 1.
//
// This will create d.Weights and d.Biases if
// they are nil.
func (d *DenseLayer) Randomize() {
	if d.Biases == nil {
		d.Biases = &autofunc.LinAdd{
			Var: &autofunc.Variable{
				Vector: make(linalg.Vector, d.OutputCount),
			},
		}
	}
	if d.Weights == nil {
		d.Weights = &autofunc.LinTran{
			Rows: d.OutputCount,
			Cols: d.InputCount,
			Data: &autofunc.Variable{
				Vector: make(linalg.Vector, d.OutputCount*d.InputCount),
			},
		}
	}

	sqrt3 := math.Sqrt(3)
	for i := 0; i < d.OutputCount; i++ {
		d.Biases.Var.Vector[i] = sqrt3 * ((rand.Float64() * 2) - 1)
	}

	weightCoeff := math.Sqrt(3.0 / float64(d.InputCount))
	for i := range d.Weights.Data.Vector {
		d.Weights.Data.Vector[i] = weightCoeff * ((rand.Float64() * 2) - 1)
	}
}

func C206() (interface{}, error) {
	min := int(math.Sqrt(10203040506070809))
	max := int(math.Sqrt(19293949596979899)) + 1

	checks := map[int]int{
		9: 1,
		8: 100,
		7: 10000,
		6: 1000000,
		5: 100000000,
		4: 10000000000,
		3: 1000000000000,
		2: 100000000000000,
		1: 10000000000000000,
	}

	for i := min; i < max; i++ {
		pow := i * i
		matched := true
		for target, divisor := range checks {
			if pow/divisor%10 != target {
				matched = false
				break
			}
		}
		if matched {
			return i * 10, nil
		}
	}

	return nil, nil
}

// check that the integration function works
func TestIntegrateMid(t *testing.T) {
	tests := []struct {
		fn     smoothFn
		x1, x2 float64
		Tot    float64
	}{
		// linear
		{func(x float64) float64 { return 0.5 * x }, 0.0, 1.0, 0.25},
		// normal distribution
		{func(x float64) float64 { return 1 / math.Sqrt(2*math.Pi) * math.Exp(-(x*x)/2) }, -100, 100, 1.0},
		// normal distribution half
		{func(x float64) float64 { return 1 / math.Sqrt(2*math.Pi) * math.Exp(-(x*x)/2) }, -100, 0, 0.5},
		// normal distribution segment
		{func(x float64) float64 { return 1 / math.Sqrt(2*math.Pi) * math.Exp(-(x*x)/2) }, -2, -1, .1359051219835},
		// scaled gamma distribution (similar to my dissertation experiment 3)
		{func(x float64) float64 {
			k, theta, a := 1.5, 2.0, 1.0/600
			return a / (math.Gamma(k) * math.Pow(theta, k)) * math.Sqrt(x*a) * math.Exp(-x*a/2)
		}, 0, 2400, 0.73853606463},
	}

	for i, test := range tests {
		got := integrateMid(test.fn, test.x1, test.x2, 10000)
		if diff := math.Abs(got - test.Tot); diff > 1e-10 {
			t.Errorf("case %v (integral from %v to %v): got %v, want %v", i+1, test.x1, test.x2, got, test.Tot)
		}
	}
}