74 lines
1.9 KiB
Go
74 lines
1.9 KiB
Go
package thermal
|
|
|
|
import (
|
|
"math"
|
|
"testing"
|
|
)
|
|
|
|
// TestConeIsCutToRepose is the property the whole hillslope story rests on: after enough passes nothing may
|
|
// stand steeper than the angle of repose. A cone is the worst case, because every cell on it is steep.
|
|
func TestConeIsCutToRepose(t *testing.T) {
|
|
const (
|
|
w, h = 101, 101
|
|
cellM = 10.0
|
|
deg = 33.0
|
|
)
|
|
talus := TalusFromDegrees(deg)
|
|
field := make([]float32, w*h)
|
|
for y := 0; y < h; y++ {
|
|
for x := 0; x < w; x++ {
|
|
// A steep cone: 60 m of fall per 10 m cell at the flanks, far beyond repose.
|
|
d := math.Hypot(float64(x-w/2), float64(y-h/2))
|
|
field[y*w+x] = float32(math.Max(0, 600-6*d*cellM))
|
|
}
|
|
}
|
|
before := maxSlope(field, w, h, cellM)
|
|
scratch := make([]float32, w*h)
|
|
|
|
volBefore := sum(field)
|
|
Apply(field, w, h, cellM, talus, 400, nil, scratch)
|
|
volAfter := sum(field)
|
|
|
|
after := maxSlope(field, w, h, cellM)
|
|
t.Logf("max slope %.1f deg -> %.1f deg (repose %.1f)", degOf(before), degOf(after), deg)
|
|
if degOf(after) > deg+2 {
|
|
t.Errorf("after 400 passes the steepest slope is %.1f deg, want no more than %.1f", degOf(after), deg+2)
|
|
}
|
|
if rel := math.Abs(volAfter-volBefore) / volBefore; rel > 1e-4 {
|
|
t.Errorf("mass changed by %.4f%%; thermal weathering must conserve it", rel*100)
|
|
}
|
|
}
|
|
|
|
func sum(a []float32) float64 {
|
|
var s float64
|
|
for _, v := range a {
|
|
s += float64(v)
|
|
}
|
|
return s
|
|
}
|
|
|
|
func degOf(slope float64) float64 { return math.Atan(slope) * 180 / math.Pi }
|
|
|
|
func maxSlope(a []float32, w, h int, cellM float64) float64 {
|
|
worst := 0.0
|
|
for y := 0; y < h; y++ {
|
|
for x := 0; x < w; x++ {
|
|
i := y*w + x
|
|
for k := 0; k < 8; k++ {
|
|
nx, ny := x+dx8[k], y+dy8[k]
|
|
if nx < 0 || ny < 0 || nx >= w || ny >= h {
|
|
continue
|
|
}
|
|
d := cellM
|
|
if dx8[k] != 0 && dy8[k] != 0 {
|
|
d = cellM * math.Sqrt2
|
|
}
|
|
if s := float64(a[i]-a[ny*w+nx]) / d; s > worst {
|
|
worst = s
|
|
}
|
|
}
|
|
}
|
|
}
|
|
return worst
|
|
}
|