GameCult/Aetheria
Watch
9
0
Fork
You've already forked Aetheria
 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Merge eureka/aetheria-gravity-well-law (cut well-law; Soul verdicts s1-s4)

Browse source 
Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
This commit is contained in:
Meta Vens 2026-10-17 03:49:59 +02:00
commit 09a7f18985
22 changed files with 1616 additions and 437 deletions
Download patch file Download diff file Expand all files Collapse all files

BIN
Assets/Resources/Settings.asset (Stored with Git LFS)
View file

Binary file not shown.

113
Assets/Scripts/Editor/GravityFieldParityVerify.cs
View file

@ -1,4 +1,5 @@
using System;
using System.Collections.Generic;
using System.IO;
using System.Linq;
using GameCult.Caching;
@ -27,14 +28,16 @@ namespace Aetheria.EditorTools
const string Tag = "[GravityFieldParityVerify]";
const int Size = 256;
// Measured on this tree (the "max error" line this run prints): the largest |texel + GetHeight| over the windows of the three
// frames and windows is 0.0065 u, in the last frame's planetoid window (float32 on both sides, the GPU's log against the CPU's, the
// quad's interpolated world position). The tolerance is 4.6x that. A 0.1% error in the tail constant (.5 -> .501)
// moves the drawn height by 0.053 u where all three tails reach (the planetoid and zone windows), 1.8x the
// tolerance; a reference of zero moves it by about 150 u.
// frames is 0.0103 u, in the moved-again frame's sun window (float32 on both sides, the GPU's log and the closed form's
// single-precision floor of 1.5e-5 per unit strength against the CPU's, the quad's interpolated world position). The
// tolerance is 2.9x that. A band whose terms were drawn from the wrong edge (reference, inner or outer) moves the
// drawn height by whole units.
const float Tolerance = 0.03f;
// The worst angle between the texture's central difference and GetForce (B28: 1.0 degrees median zone, 2.9 worst).
const float GradientDegrees = 3f;
static float _maxError, _maxAngle;
// Below this force (u/s^2) the drawn height differs between neighbouring texels by less than float32 resolves.
const float FlatForce = 1e-3f;
static float _maxError, _maxAngle, _maxTau;
static bool _ok = true;
static void OnLog(string message, string stack, LogType type)
@ -70,22 +73,29 @@ namespace Aetheria.EditorTools
var material = AssetDatabase.LoadAssetAtPath<Material>("Assets/Materials/Brushes/Gravity Field.mat")
?? throw new InvalidOperationException("No Gravity Field.mat.");
var zone = BuildZone(settings, work, out var orbits);
var zone = BuildZone(settings, work, out var orbits, multiRoot: true);
GravityLowering.Lower(zone, 0f, ref wells);
field = GravityLowering.CreateField(material);
// Three frames through the one buffer, as ZoneRenderer.Update runs it: the bodies as built, then two moves
// (an orbit's Position is what Zone.GetHeight reads). A lowering that mirrors an axis or stops refreshing
// the upload agrees with the sim in the first frame only where it is symmetric, and not after a move.
// Three frames through the one buffer, as ZoneRenderer.Update runs it: the bodies as built, then two moves (an
// orbit's Position is what Zone.GetHeight reads). A move keeps the tree's geometry, since a rail never leaves
// its parent's distance and the walk skips a subtree outside its reach: the sun is the root, the planet is
// its Distance from the sun, the planetoid its Distance from the planet. A lowering that mirrors an axis or
// stops refreshing the upload agrees with the sim in the first frame only where it is symmetric.
CheckFrame("built", zone, orbits);
orbits[1].Position = new float2(-420f, 560f);
orbits[2].Position = new float2(900f, -700f);
orbits[1].Position = Around(orbits[0], orbits[1], .31f);
orbits[2].Position = Around(orbits[1], orbits[2], .77f);
orbits[5].Position = Around(orbits[1], orbits[5], .41f);
GravityLowering.Lower(zone, 0f, ref wells);
CheckFrame("moved", zone, orbits);
orbits[0].Position = new float2(210f, -95f);
orbits[2].Position = new float2(-1800f, 1300f);
orbits[1].Position = Around(orbits[0], orbits[1], .52f);
orbits[2].Position = Around(orbits[1], orbits[2], .05f);
orbits[5].Position = Around(orbits[1], orbits[5], .93f);
GravityLowering.Lower(zone, 0f, ref wells);
CheckFrame("moved-again", zone, orbits);
if (_maxTau <= 1f) Fail("series-branch", $"the largest tau over the drawn bounded terms is {_maxTau:G4}; no window draws the series branch");
Debug.Log($"{Tag} largest drawn tau {_maxTau:G4}");
Debug.Log($"{Tag} max error {_maxError:G4} u (tolerance {Tolerance:G4}); max gradient angle {_maxAngle:F2} degrees (limit {GradientDegrees:F0})");
CheckWiring(settings);
@ -108,35 +118,69 @@ namespace Aetheria.EditorTools
static float Core(Zone zone, float mass) => zone.Settings.RocheLimit.Evaluate(mass);
// Where an orbit lies when it stands `turns` round its parent at its Distance.
static float2 Around(Orbit parent, Orbit child, float turns) => parent.Position + child.Data.Distance * OrbitData.Evaluate(turns);
static Vector2 At(Orbit orbit) => new Vector2(orbit.Position.x, orbit.Position.y);
// Orbits 0-2 are the sun, the planet and the planetoid, each window centred where the body is now; orbit 3 is the
// body of no strength. The zone window holds all four.
// body of no strength, orbit 4 the second sun and orbit 5 the planet's moon inside its Roche limit (whose band has
// tau above 1: the series branch of the shader). The zone window holds all of them; the far window is centred on the
// zone node, so its band (nothing inside A, the zone's mass beyond Beta A) is drawn. A window draws a term when
// the term's disc reaches it; the largest tau over the drawn terms must exceed 1 or the series branch went unchecked.
static void CheckFrame(string frame, Zone zone, Orbit[] orbits)
{
CheckTexels(frame + "-planetoid", zone, At(orbits[2]), 4f * Core(zone, 29.2f));
CheckTexels(frame + "-planet", zone, At(orbits[1]), 4f * Core(zone, 273.9f));
CheckTexels(frame + "-sun", zone, At(orbits[0]), 4f * Core(zone, 14152f));
CheckTexels(frame + "-zone", zone, Vector2.zero, 1.05f * zone.Pack.Radius);
var centre = new Vector2(zone.Wells.Centre.x, zone.Wells.Centre.y);
var span = 0f;
foreach (var root in new[] { orbits[0], orbits[4] })
span = Math.Max(span, (At(root) - centre).magnitude + zone.Wells.Reach(zone.Orbits.First(o => o.Value == root).Key));
// The wells have their own windows: out here a texel is 40 u wide, so only the ring from 0.9 A outward, where the
// zone node's band rises and the zone becomes one body, is measured.
CheckTexels(frame + "-far", zone, centre, 1.6f * span, .9f * span);
CheckGradient(zone, orbits.Select(At).ToArray());
}
// A sun (14152 at (-130, 85)), a planet (273.9 at (610, -240)), a planetoid (29.2 at (670, -270)) and a body of
// mass 0, whose well has no strength (a legal body; its curves give c = 0 and R = 0), at (1700, 300): parentless,
// fixed, off both axes and asymmetric, so the zone's wells sit where the verify says before any Update and a
// mirrored or swapped axis cannot agree with the sim.
internal static Zone BuildZone(GameSettings settings, string work, out Orbit[] orbits, float radius = 4182f)
// A sun (14152 at (-130, 85)), its planet (273.9, placed at (610, -240) at the pack's time), the planet's planetoid
// (29.2, placed at (670, -270)) and a body of mass 0, whose well has no strength (a legal body; its curves give
// c = 0 and R = 0), at (1700, 300): a family, so the band terms are drawn (the planet's and the sun's cusps and
// fades), off both axes and asymmetric, so the zone's wells sit where the verify says before any Update and a
// mirrored or swapped axis cannot agree with the sim. The sun and the zero-mass body are fixed roots.
internal static Zone BuildZone(GameSettings settings, string work, out Orbit[] orbits, float radius = 4182f, bool multiRoot = false)
{
var cache = AetheriaStores.Open(Path.Combine(work, "catalog.cc"), Path.Combine(work, "run.cc"), catalogWritable: true);
cache.Upsert(new VerseGrammar { Revision = 1 });
cache.FlushAsync().Wait();
var items = new ItemManager(cache, new ProvenanceLedger(), settings.GameplaySettings, _ => { });
var positions = new[] { new Vector2(-130f, 85f), new Vector2(610f, -240f), new Vector2(670f, -270f), new Vector2(1700f, 300f) };
var masses = new[] { 14152f, 273.9f, 29.2f, 0f };
var keys = new CultRecordKey[positions.Length];
var pack = new ZonePack { Radius = radius };
for (var i = 0; i < positions.Length; i++)
var positions = new List<Vector2> { new Vector2(-130f, 85f), new Vector2(610f, -240f), new Vector2(670f, -270f), new Vector2(1700f, 300f) };
var masses = new List<float> { 14152f, 273.9f, 29.2f, 0f };
var parents = new List<int> { -1, 0, 1, -1 };
if (multiRoot)
{
var orbit = keys[i] = cache.Upsert(new OrbitData { FixedPosition = new float2(positions[i].x, positions[i].y) }).Key;
// A second sun (a second subzone, so the zone node exists) and a moon of the planet 0.6 of its Roche limit away.
positions.Add(new Vector2(-2400f, -1900f));
masses.Add(3000f);
parents.Add(-1);
positions.Add(positions[1] + new Vector2(-.6f * settings.PlanetSettings.RocheLimit.Evaluate(273.9f), 0f));
masses.Add(1f);
parents.Add(1);
}
var keys = new CultRecordKey[positions.Count];
var pack = new ZonePack { Radius = radius };
for (var i = 0; i < positions.Count; i++)
{
var data = new OrbitData { FixedPosition = new float2(positions[i].x, positions[i].y) };
if (parents[i] >= 0)
{
var offset = positions[i] - positions[parents[i]];
data.Parent = new CultRecordRef<OrbitData>(keys[parents[i]]);
data.Distance = offset.magnitude;
data.Phase = Mathf.Atan2(offset.y, offset.x) / (2f * Mathf.PI);
}
var orbit = keys[i] = cache.Upsert(data).Key;
var body = cache.Upsert(new PlanetData
{ Name = "Body " + i, Orbit = new CultRecordRef<OrbitData>(orbit), Mass = masses[i] }).Key;
pack.Orbits.Add(new CultRecordRef<OrbitData>(orbit));
@ -195,13 +239,29 @@ namespace Aetheria.EditorTools
}
}
static void CheckTexels(string window, Zone zone, Vector2 centre, float half)
// The largest tau = (a^2 + eps^2) / (b^2 - a^2) over the bounded terms whose disc reaches the window's square.
static void NoteTau(Zone zone, Vector2 centre, float half)
{
foreach (var term in zone.Wells.Terms)
{
if (float.IsPositiveInfinity(term.Outer)) continue;
var at = term.Node.IsSet() ? zone.Orbits[term.Node].Position : zone.Wells.Centre;
if (Math.Abs(at.x - centre.x) > half + term.Outer || Math.Abs(at.y - centre.y) > half + term.Outer) continue;
_maxTau = Math.Max(_maxTau, (term.Inner * term.Inner + term.Softening * term.Softening) / (term.Outer * term.Outer - term.Inner * term.Inner));
}
}
static void CheckTexels(string window, Zone zone, Vector2 centre, float half, float ringFrom = 0f)
{
Render(centre, half, out var depth, out var world);
NoteTau(zone, centre, half);
var worst = 0f;
var live = 0;
var measured = 0;
for (var i = 0; i < depth.Length; i++)
{
if ((world[i] - centre).magnitude < ringFrom) continue;
measured++;
var expected = -zone.GetHeight(new float2(world[i].x, world[i].y));
if (float.IsNaN(depth[i])) { Fail($"{window}-texels", $"texel {i}: NaN"); return; }
var error = Math.Abs(depth[i] - expected);
@ -215,7 +275,7 @@ namespace Aetheria.EditorTools
}
_maxError = Math.Max(_maxError, worst);
// A window that draws nothing agrees with a zone that is flat: the check must be live.
if (live < depth.Length / 2) Fail($"{window}-texels", $"only {live} of {depth.Length} texels are deeper than 1 u; the check is vacuous");
if (live < measured / 2) Fail($"{window}-texels", $"only {live} of {measured} texels are deeper than 1 u; the check is vacuous");
Debug.Log($"{Tag} {window}: {depth.Length} texels, max error {worst:G4} u");
}
@ -239,6 +299,7 @@ namespace Aetheria.EditorTools
var det = au.x * av.y - au.y * av.x;
var g = new Vector2((du * av.y - dv * au.y) / det, (dv * au.x - du * av.x) / det);
var f = zone.GetForce(new float2(here.x, here.y));
if (CultMath.math.length(f) < FlatForce) continue; // between the subzones the ground is flat: no direction to compare
var angle = Vector2.Angle(g, new Vector2(f.x, f.y));
if (float.IsNaN(angle) || angle > GradientDegrees)
{

18
Assets/Scripts/ServerShared/DrawAhead.cs
View file

@ -51,24 +51,6 @@ public static class DrawAhead
return float3(at.x + bearing.x * distance, at.y, at.z + bearing.y * distance);
}
// One well as the shader takes it: its body's rail at the drawn time, with the strength, core and reference the sim's
// GetHeight sums with. GravityLowering copies it into the shader's record and decides nothing.
public readonly struct DrawnWell
{
public readonly float2 Position;
public readonly float Strength, CoreRadius, Reference;
public DrawnWell(float2 position, float strength, float coreRadius, float reference)
{
Position = position;
Strength = strength;
CoreRadius = coreRadius;
Reference = reference;
}
}
public static DrawnWell Well(Zone zone, Planet body, float lead) =>
new DrawnWell(Orbit(zone, body.BodyData.Orbit.Key, lead), body.WellStrength, body.RocheLimit, zone.WellReference(body));
// A free body (a vapour cloud) drawn where the frame falls lead sim seconds after the last step, by the velocity
// that step moved it with: the launch velocity and the drift together.
public static float2 Position(KinematicBody body, float lead) => body.Position + body.TotalVelocity * lead;

158
Assets/Scripts/ServerShared/GravityField.cs
View file

@ -1,32 +1,150 @@
using System;
using CultMath;
using static CultMath.math;
using float2 = CultMath.float2;
/// <summary>
/// The law of one gravity well: a 2D point mass with a softened core. Outside the core radius the potential is
/// logarithmic and the pull is c / r; inside it the pull falls linearly to zero at the centre and meets the tail in
/// value and slope at the core radius. c is the well's strength (u^2/s^2), R its core radius (u). No threshold
/// anywhere: a well pulls at every distance. Zone.GetForce and Zone.GetHeight sum this over a zone's bodies.
/// The law of one band term U(r; eps, a, b), the only primitive of a zone's gravity (WellTree decides which terms
/// exist). A Plummer-softened log source whose pull toward the centre is exact up to a and fades to nothing over [a, b]:
/// pull r / (r^2 + eps^2) * B(s), s = r^2, B = 1 for s &lt;= a^2, 1 - S(t) on the band with t = (s - a^2) / (b^2 - a^2),
/// 0 for s &gt;= b^2; S is the order-Order smoothstep. The potential is 0 beyond b. With a = b = +inf it is the global term
/// 1/2 ln(r^2 + eps^2). The potential is a function of r^2 (smooth at the centre), and C^Order on the circles r = a and
/// r = b (C4 at Order 4: the player feels higher-order seams), so the drawn height is exactly the sim's potential.
/// On the band the closed form is 1/2 [R ln((t + tau) / (1 + tau)) - (I(1) - I(t))], tau = (a^2 + eps^2) / (b^2 - a^2),
/// (1 - S(t)) = (t + tau) Q(t) + R, I the antiderivative of Q; for tau &gt; 1 it cancels catastrophically in float32
/// (7e-3 at tau 3, 12 at tau 10), so the band is summed as a series in (1 - t) / (1 + tau), whose ratio is below 1/2.
/// Pull and Potential are per unit strength; the caller scales by the term's strength (WellStrength at a mass).
/// Gravity.cginc is generated from this file's constants and coefficients (GravityShader).
/// </summary>
public static class GravityField
{
// The tail's additive constant: it makes the tail meet the core in value at the core radius. Gravity.cginc is generated from it.
public const float TailOffset = .5f;
/// <summary>The smoothstep's order: the potential is C^Order on every band edge.</summary>
public const int Order = 4;
/// <summary>A family's charge is fully on by its support, Beta x its outermost child.</summary>
public const float Beta = 1.5f;
/// <summary>A family is felt whole out to max(support, Gamma x its reach) and fades from there to its reach.</summary>
public const float Gamma = .5f;
/// <summary>A body's reach in its parent's field, in its distance (the nominal reach before the fade's room).</summary>
public const float Reach = .5f;
/// <summary>A ring member of a bodyless node reaches RingReach x its distance: whole out to the ring's diameter.</summary>
public const float RingReach = 4f;
// The potential at distance r from the body, defined up to a constant (callers subtract their reference).
// A well of no strength (a body whose curves are unauthored, so c = 0 and R = 0) digs nothing.
public static float Potential(float c, float coreRadius, float r) =>
c == 0f ? 0f :
r < coreRadius
? c * r * r / (2 * coreRadius * coreRadius)
: c * (log(r / coreRadius) + TailOffset);
/// <summary>The band is summed as a series for tau above this. float32 error per unit strength against the double integral is at most 1.8e-6 anywhere (worst 1.7e-6 just above the seam, where the series has not reached float32 resolution); the closed form is float32-clean only to about 0.45.</summary>
public const float SeriesFrom = .4f;
// The acceleration at `offsetFromBody`, toward the body: zero at the body itself.
public static float2 Acceleration(float c, float coreRadius, float2 offsetFromBody)
/// <summary>The series branch's length: its ratio is below 1/1.4, so the truncated tail is below 1.4e-6 at the seam and falls geometrically above it.</summary>
public const int SeriesTerms = 40;
/// <summary>S(t) as polynomial coefficients in t, lowest first: t^Order sum_j C(Order-1+j, j) (1-t)^j.</summary>
public static readonly float[] S = BuildS();
/// <summary>1 - S(t) as polynomial coefficients in t, lowest first: the band's profile B.</summary>
public static readonly float[] OneMinusS = BuildOneMinusS();
private static float[] BuildS()
{
var r = length(offsetFromBody);
if (r == 0f) return float2.zero;
var pull = r < coreRadius ? c * r / (coreRadius * coreRadius) : c / r;
return -offsetFromBody / r * pull;
var c = new double[2 * Order];
for (var j = 0; j < Order; j++)
{
var cj = Binomial(Order - 1 + j, j);
for (var i = 0; i <= j; i++)
c[Order + i] += cj * Binomial(j, i) * (i % 2 == 0 ? 1 : -1);
}
var s = new float[2 * Order];
for (var i = 0; i < s.Length; i++) s[i] = (float) c[i];
return s;
}
private static float[] BuildOneMinusS()
{
var f = new float[2 * Order];
for (var i = 0; i < f.Length; i++) f[i] = -S[i];
f[0] += 1f;
return f;
}
private static double Binomial(int n, int k)
{
var r = 1.0;
for (var i = 1; i <= k; i++) r = r * (n - k + i) / i;
return r;
}
// B(s) for the band term (a, b); b = +inf is no fade at all.
private static float Profile(float inner, float outer, float s)
{
if (float.IsPositiveInfinity(outer)) return 1f;
var a2 = inner * inner;
var b2 = outer * outer;
if (s <= a2) return 1f;
if (s >= b2) return 0f;
var t = (s - a2) / (b2 - a2);
var acc = 0f;
for (var i = OneMinusS.Length - 1; i >= 0; i--) acc = acc * t + OneMinusS[i];
return acc;
}
/// <summary>The pull toward the centre at distance r, per unit strength.</summary>
public static float Pull(float softening, float inner, float outer, float r)
{
var s = r * r;
return r * Profile(inner, outer, s) / (s + softening * softening);
}
/// <summary>The potential at distance r, per unit strength: 0 at and beyond outer, the global log when outer is +inf.</summary>
public static float Potential(float softening, float inner, float outer, float r)
{
var s = r * r;
var e2 = softening * softening;
if (float.IsPositiveInfinity(outer)) return .5f * log(s + e2);
var a2 = inner * inner;
var b2 = outer * outer;
if (s >= b2) return 0f;
var tau = (a2 + e2) / (b2 - a2);
if (s <= a2) return .5f * log((s + e2) / (a2 + e2)) + Band(tau, 0f);
return Band(tau, (s - a2) / (b2 - a2));
}
// The band's potential at t in [0, 1], 0 at t = 1.
private static float Band(float tau, float t) => tau <= SeriesFrom ? BandClosed(tau, t) : BandSeries(tau, t);
private static float BandClosed(float tau, float t)
{
// (1 - S(t)) = (t + tau) Q(t) + R by synthetic division; I is Q's antiderivative.
var n = OneMinusS.Length - 1;
Span<float> q = stackalloc float[n];
q[n - 1] = OneMinusS[n];
for (var i = n - 1; i > 0; i--) q[i - 1] = OneMinusS[i] - tau * q[i];
var remainder = OneMinusS[0] - tau * q[0];
return .5f * (remainder * log((t + tau) / (1f + tau)) - (Antiderivative(q, 1f) - Antiderivative(q, t)));
}
private static float Antiderivative(ReadOnlySpan<float> q, float t)
{
var acc = 0f;
for (var i = q.Length - 1; i >= 0; i--) acc = (acc + q[i] / (i + 1)) * t;
return acc;
}
private static float BandSeries(float tau, float t)
{
// With y = 1 - t, 1 - S(1 - y) = S(y) and 1 / (t + tau) = sum_k y^k / (1 + tau)^(k+1):
// band = -1/2 sum_k (1 + tau)^-(k+1) sum_j S_j y^(j+k+1) / (j+k+1).
var y = 1f - t;
Span<float> yj = stackalloc float[S.Length];
yj[0] = 1f;
for (var j = 1; j < yj.Length; j++) yj[j] = yj[j - 1] * y;
var w = 1f / (1f + tau);
var wk = w;
var yk = y;
var acc = 0f;
for (var k = 0; k < SeriesTerms; k++)
{
var inner = 0f;
for (var j = Order; j < S.Length; j++) inner += S[j] * yj[j] / (j + k + 1);
acc += wk * yk * inner;
wk *= w;
yk *= y;
}
return -.5f * acc;
}
}

101
Assets/Scripts/ServerShared/GravityShader.cs
View file

@ -1,10 +1,12 @@
using System.Globalization;
using System.Linq;
/// <summary>
/// The single author of Assets/Shaders/Gravity.cginc, the GPU lowering of the zone's potential. The law is
/// GravityField's (its constants are read here, not respelled) and the height is Zone.GetHeight: each well's
/// potential minus its reference, summed, divided by the exaggeration. Unity does not run C# when it imports a
/// shader, so the committed cginc is this output; GravityShaderTests fails unless it equals Source() exactly.
/// GravityField's (its constants and smoothstep coefficients are read here, not respelled) and the height is
/// Zone.GetHeight: each term's potential minus its reference, summed over the terms whose edge the point is inside,
/// divided by the exaggeration. Unity does not run C# when it imports a shader, so the committed cginc is this output;
/// GravityShaderTests fails unless it equals Source() exactly.
/// To change the shader, change this file and regenerate: AETHERIA_REGENERATE_GRAVITY_CGINC=1 dotnet test --filter GravityShaderTests.
/// </summary>
public static class GravityShader
@ -12,48 +14,111 @@ public static class GravityShader
public const string Path = "Assets/Shaders/Gravity.cginc";
/// <summary>The whole file, with LF line ends.</summary>
public static string Source() => (Template
.Replace("@TAIL@", GravityField.TailOffset.ToString("R", CultureInfo.InvariantCulture)))
public static string Source() => Template
.Replace("@ORDER@", GravityField.Order.ToString(CultureInfo.InvariantCulture))
.Replace("@N@", GravityField.S.Length.ToString(CultureInfo.InvariantCulture))
.Replace("@Q@", (GravityField.S.Length - 1).ToString(CultureInfo.InvariantCulture))
.Replace("@QLAST@", (GravityField.S.Length - 2).ToString(CultureInfo.InvariantCulture))
.Replace("@SERIESFROM@", GravityField.SeriesFrom.ToString("R", CultureInfo.InvariantCulture))
.Replace("@SERIES@", GravityField.SeriesTerms.ToString(CultureInfo.InvariantCulture))
.Replace("@ONEMINUSS@", Coefficients(GravityField.OneMinusS))
.Replace("@S@", Coefficients(GravityField.S))
.Replace("\r\n", "\n");
private static string Coefficients(float[] c) =>
string.Join(", ", c.Select(x => x.ToString("R", CultureInfo.InvariantCulture)));
private const string Template = @"#ifndef AETHERIA_GRAVITY_INCLUDED
#define AETHERIA_GRAVITY_INCLUDED
// GENERATED by ServerShared/GravityShader.cs from GravityField; do not edit (GravityShaderTests fails on any difference).
// The lowering of the zone's potential (aetheria-gravity, cut grid-lowering). The law is ServerShared's GravityField
// and the height is Zone.GetHeight; this file decides nothing. GravityLowering.Lower fills the buffer and the globals
// from the zone (positions from Orbit.Position, reference from Zone.WellReference) each frame.
// from the zone's WellTree (positions from the nodes' rails at the drawn time) each frame.
// GravityFieldParityVerify renders it and compares every texel with -Zone.GetHeight.
struct GravityWell
struct GravityTerm
{
float2 position;
float strength;
float coreRadius;
float softening;
float inner;
float outer;
float reference;
};
StructuredBuffer<GravityWell> _GravityWells;
StructuredBuffer<GravityTerm> _GravityWells;
int _GravityWellCount;
float _HeightExaggeration;
// GravityField.Potential: c is the well's strength, R its core radius, r the distance.
float gravityPotential(float c, float R, float r)
// GravityField.OneMinusS and GravityField.S: the band profile B(t) = 1 - S(t) and the order-@ORDER@ smoothstep, lowest first.
static const float GravityOneMinusS[@N@] = { @ONEMINUSS@ };
static const float GravityS[@N@] = { @S@ };
// GravityField.BandClosed: the band's potential at t for tau <= @SERIESFROM@, per unit strength.
float gravityBandClosed(float tau, float t)
{
return c == 0 ? 0 :
r < R
? c * r * r / (2 * R * R)
: c * (log(r / R) + @TAIL@);
float q[@Q@];
q[@QLAST@] = GravityOneMinusS[@Q@];
[unroll] for (int i = @QLAST@; i > 0; i--) q[i - 1] = GravityOneMinusS[i] - tau * q[i];
float remainder = GravityOneMinusS[0] - tau * q[0];
float whole = 0;
float part = 0;
[unroll] for (int j = @QLAST@; j >= 0; j--)
{
whole = whole + q[j] / (j + 1);
part = (part + q[j] / (j + 1)) * t;
}
return 0.5 * (remainder * log((t + tau) / (1 + tau)) - (whole - part));
}
// Zone.GetHeight at world xz: each well's potential minus its own reference, summed, divided by the exaggeration.
// GravityField.BandSeries: the band's potential at t for tau > @SERIESFROM@, a series in (1 - t) / (1 + tau).
float gravityBandSeries(float tau, float t)
{
float y = 1 - t;
float yj[@N@];
yj[0] = 1;
[unroll] for (int j = 1; j < @N@; j++) yj[j] = yj[j - 1] * y;
float w = 1 / (1 + tau);
float wk = w;
float yk = y;
float acc = 0;
for (int k = 0; k < @SERIES@; k++)
{
float inner = 0;
[unroll] for (int m = @ORDER@; m < @N@; m++) inner += GravityS[m] * yj[m] / (m + k + 1);
acc += wk * yk * inner;
wk *= w;
yk *= y;
}
return -0.5 * acc;
}
// GravityField.Potential: one band term, per unit strength; the global term when b is infinite.
float gravityTermPotential(float softening, float a, float b, float r)
{
float s = r * r;
float e2 = softening * softening;
if (isinf(b)) return 0.5 * log(s + e2);
float a2 = a * a;
float b2 = b * b;
if (s >= b2) return 0;
float tau = (a2 + e2) / (b2 - a2);
if (s <= a2) return 0.5 * log((s + e2) / (a2 + e2)) + (tau <= @SERIESFROM@ ? gravityBandClosed(tau, 0) : gravityBandSeries(tau, 0));
float t = (s - a2) / (b2 - a2);
return tau <= @SERIESFROM@ ? gravityBandClosed(tau, t) : gravityBandSeries(tau, t);
}
// Zone.GetHeight at world xz: each term's potential minus its own reference, summed, divided by the exaggeration.
float gravityHeight(float2 p)
{
float result = 0;
for (int i = 0; i < _GravityWellCount; i++)
{
GravityWell well = _GravityWells[i];
result += gravityPotential(well.strength, well.coreRadius, length(p - well.position)) - well.reference;
GravityTerm term = _GravityWells[i];
float r = length(p - term.position);
if (r < term.outer)
result += term.strength * gravityTermPotential(term.softening, term.inner, term.outer, r) - term.reference;
}
return result / _HeightExaggeration;
}

6
Assets/Scripts/ServerShared/Scenarios/LockOn.cs
View file

@ -12,9 +12,9 @@ public sealed class LockOn : Scenario
// The hostile starts about 4000 m off, ahead and to port of the player, facing the player. At that range its passive
// sensors hold no track of the player (the test asks the sim's own PassiveInfo against the detection threshold). The
// hostile closes, sees the player some fourteen seconds in and locks. The zone's pull bends the hostile's run and moves the
// player, and the warning must still point at where the hostile was staged to within 5 degrees at the lock: this side keeps
// the bend under one degree (on the far side it is more than seven). At Duel's 900 m (and at 600 m) it sees and locks
// hostile closes, sees the player some fourteen seconds in and locks. Its flanking approach and the zone's pull carry it
// off the staged line (about 5 degrees by the lock on the merged tree), and the warning must point at where the hostile
// is at the lock, which the test checks against the live bearing. At Duel's 900 m (and at 600 m) it sees and locks
// within a second.
public static readonly CultMath.float2 HostileAt = float2(-2000, 3450);

7
Assets/Scripts/ServerShared/Settings.cs
View file

@ -16,11 +16,12 @@ public class PlanetSettings
public ExponentialCurve LightRadius;
public ExponentialCurve BodyRadius;
// The closest a child may orbit its parent, in u, by the parent's mass, and the softened core radius of the parent's
// well (the steepest drawn slope sits at the Roche limit). Authored like BodyRadius, never per body.
// The closest a child may orbit its parent, in u, by the parent's mass, and the Plummer softening of the parent's
// well (the steepest slope, 0.5 c / eps, sits at the Roche limit). Authored like BodyRadius, never per body.
public ExponentialCurve RocheLimit = new ExponentialCurve();
// A body's well strength by its mass (GravityField): c in u^2/s^2 (the pull is c / r beyond the core, which is the Roche limit).
// A body's well strength by its mass (GravityField): c in u^2/s^2. The pull is c d / (d^2 + eps^2) inside a family and
// c(family mass) beyond it (WellTree); eps is the Roche limit.
public ExponentialCurve WellStrength = new ExponentialCurve();
// The one global exaggeration of gravity height: height = potential / K, in u/s^2 per unit of height.

431
Assets/Scripts/ServerShared/WellTree.cs Normal file
View file

@ -0,0 +1,431 @@
using System;
using System.Collections.Generic;
using GameCult.Caching;
using CultMath;
using static CultMath.math;
using float2 = CultMath.float2;
/// <summary>
/// An orbit as the well tree reads it: where it hangs on the tree and the mass of the well-bearing body on it (0 for none).
/// Extent is the outermost distance from the parent of anything riding this orbit's rail (a belt orbit: its outermost
/// asteroid; any other orbit: its Distance). Anchor is a parentless orbit's stationary position, ignored with a parent.
/// </summary>
public readonly struct WellNode
{
public readonly CultRecordKey Key, Parent;
public readonly float Distance, Phase, Mass, Extent;
public readonly float2 Anchor;
public WellNode(CultRecordKey key, CultRecordKey parent, float distance, float phase, float mass, float extent, float2 anchor)
{
Key = key;
Parent = parent;
Distance = distance;
Phase = phase;
Mass = mass;
Extent = extent;
Anchor = anchor;
}
}
/// <summary>
/// The one owner of which gravity terms a zone has: the orbit tree turned into band terms (GravityField is the law of one
/// term). Built once from the pack's orbits and the authored curves, never reading time; the caller passes positions.
///
/// Every node F (a body, or a bodyless node with its ring group standing in for the body) contributes at its own position
/// c(m_F) U(eps(m_F), L_F, B_F) + c(M_F) [ U(eps(M_F), max(B_F, Gamma Lam_F), Lam_F) - U(eps(M_F), L_F, B_F) ]
/// m_F is the body's mass, M_F the family's, c and eps the WellStrength and RocheLimit curves. L_F is the outermost
/// child's extent: inside it the node pulls exactly as its body alone. B_F, the support, bounds every descendant term:
/// beyond it the family pulls exactly as one body of mass M_F, out to max(B_F, Gamma Lam_F), and is gone at its reach
/// Lam_F. A zone's sole massive root has infinite reach, so its far field is exactly c(M_zone). A zone of several massive
/// roots (subzones) has one zone node, which is not an orbit: at the roots' centre of mass it owns the only global term,
/// c(M_zone) U_global less c(M_zone) U(eps, A, Beta A), nothing inside A and the whole zone as one body beyond Beta A;
/// each root is bounded to its own subzone (SubzoneReach).
/// A body with no massive descendants is one term c(m) U(eps(m), Gamma Lam, Lam). Beyond a node's reach its whole
/// subtree contributes exactly nothing, so Force and Potential skip it: the walk is exact, not an approximation.
/// </summary>
public sealed class WellTree
{
/// <summary>Orbits whose distances differ by no more than this are on one ring (a binary, a ring of siblings).</summary>
public const float RingWindow = .1f;
/// <summary>A band term of the field: Strength x GravityField's U(Softening, Inner, Outer) at Node's position, less Reference.</summary>
public readonly struct Term
{
public readonly CultRecordKey Node;
public readonly float Strength, Softening, Inner, Outer, Reference;
public Term(CultRecordKey node, float strength, float softening, float inner, float outer, float reference)
{
Node = node;
Strength = strength;
Softening = softening;
Inner = inner;
Outer = outer;
Reference = reference;
}
}
private readonly Dictionary<CultRecordKey, int> _index = new Dictionary<CultRecordKey, int>();
private CultRecordKey[] _key;
private int[][] _kids;
private int[] _roots;
private int[][] _ring;
private int[] _termStart, _termCount;
private float[] _distance, _phase, _reach, _support, _extent;
private float2[] _anchor;
private int[] _massiveRoots;
private float _packRadius;
private int _zoneStart, _zoneCount;
private Term[] _terms;
/// <summary>Where the zone node's terms sit (the family-mass-weighted mean of the massive roots' anchors); meaningful when it has terms.</summary>
public float2 Centre { get; private set; }
/// <summary>Every term of the zone, grouped by node, the zone node's last (its Node is unset: it sits at Centre): what the lowering uploads and the shader sums.</summary>
public IReadOnlyList<Term> Terms => _terms;
private WellTree() { }
/// <summary>A node's reach in its parent's field (+inf for a zone's sole massive root; a subzone's radius for a root of several); 0 for an orbit that carries no mass.</summary>
public float Reach(CultRecordKey node) => _reach[_index[node]];
/// <summary>A node's support B: every term of its descendants lies inside it, so beyond it the family is one body.</summary>
public float Support(CultRecordKey node) => _support[_index[node]];
/// <summary>
/// A well's strength and softening at a mass, from the authored curves (the one place they are read and checked).
/// A mass of 0 is no well.
/// </summary>
public static (float Strength, float Softening) Curves(PlanetSettings settings, float mass)
{
var strength = settings.WellStrength.Evaluate(mass);
var softening = settings.RocheLimit.Evaluate(mass);
if (!(strength >= 0f) || float.IsInfinity(strength))
throw new InvalidOperationException("A body's well strength must be finite and not negative (PlanetSettings.WellStrength).");
if (strength > 0 && !(softening > 0f && !float.IsInfinity(softening)))
throw new InvalidOperationException("A body's Roche limit, which is its softening, must be positive and finite (PlanetSettings.RocheLimit).");
return (strength, softening);
}
/// <summary>
/// The tree of a zone's orbits. `packRadius` is where the global terms are zero: each carries its own potential there
/// as its Reference (a bounded term is already 0 beyond its outer edge, so it needs none).
/// </summary>
public static WellTree Build(IReadOnlyList<WellNode> nodes, PlanetSettings settings, float packRadius)
{
var tree = new WellTree();
var n = nodes.Count;
tree._key = new CultRecordKey[n];
tree._distance = new float[n];
tree._phase = new float[n];
tree._extent = new float[n];
tree._anchor = new float2[n];
tree._packRadius = packRadius;
for (var i = 0; i < n; i++)
{
tree._extent[i] = nodes[i].Extent;
tree._anchor[i] = nodes[i].Anchor;
tree._index[nodes[i].Key] = i;
tree._key[i] = nodes[i].Key;
tree._distance[i] = nodes[i].Distance;
tree._phase[i] = nodes[i].Phase;
}
var parent = new int[n];
var kids = new List<int>[n];
var depth = new int[n];
for (var i = 0; i < n; i++) kids[i] = new List<int>();
for (var i = 0; i < n; i++)
{
parent[i] = nodes[i].Parent.IsSet() && tree._index.TryGetValue(nodes[i].Parent, out var p) ? p : -1;
if (parent[i] >= 0) kids[parent[i]].Add(i);
}
for (var i = 0; i < n; i++)
for (var a = parent[i]; a >= 0; a = parent[a]) depth[i]++;
var roots = new List<int>();
for (var i = 0; i < n; i++)
if (parent[i] < 0) roots.Add(i);
tree._roots = roots.ToArray();
var deepestFirst = new int[n];
for (var i = 0; i < n; i++) deepestFirst[i] = i;
Array.Sort(deepestFirst, (x, y) => depth[y] != depth[x] ? depth[y].CompareTo(depth[x]) : x.CompareTo(y));
var body = new float[n]; // the body's own mass
var family = new float[n]; // the family's mass M_F
for (var i = 0; i < n; i++) family[i] = body[i] = nodes[i].Mass;
foreach (var i in deepestFirst)
if (parent[i] >= 0) family[parent[i]] += family[i];
// Massive roots are the zone's subzones; with two or more, their centre of mass is the zone node's position.
var massive = new List<int>();
foreach (var r in roots)
if (family[r] > 0) massive.Add(r);
tree._massiveRoots = massive.ToArray();
var centre = float2.zero;
var total = 0f;
foreach (var r in massive)
{
centre += tree._anchor[r] * family[r];
total += family[r];
}
tree.Centre = total > 0f ? centre / total : float2.zero;
tree._kids = new int[n][];
var ring = tree._ring = new int[n][];
for (var i = 0; i < n; i++)
{
tree._kids[i] = kids[i].ToArray();
ring[i] = RingOf(kids[i], family, tree._distance);
}
BuildReaches(tree, parent, body, family, ring, deepestFirst, out var lower, out var group);
for (var i = 0; i < n; i++)
if (body[i] > 0) ring[i] = null; // a ring stands in for a body only on a bodyless node
tree.EmitTerms(parent, body, family, ring, lower, group, settings, packRadius);
return tree;
}
// The massive children of a node at its innermost massive distance, within RingWindow of it; null when there are none.
private static int[] RingOf(List<int> kids, float[] family, float[] distance)
{
var inner = float.PositiveInfinity;
foreach (var k in kids)
if (family[k] > 0) inner = Math.Min(inner, distance[k]);
if (float.IsPositiveInfinity(inner)) return null;
var members = new List<int>();
foreach (var k in kids)
if (family[k] > 0 && distance[k] - inner <= RingWindow) members.Add(k);
return members.ToArray();
}
// The massive children of a bodyless node that stand in for its body: their mass, the ring's distance and its support.
private readonly struct RingGroup
{
public readonly float Mass, Inner, Support;
public RingGroup(float mass, float inner, float support)
{
Mass = mass;
Inner = inner;
Support = support;
}
}
// Children before parents: a family's support needs its children's reaches. lower = L_F, the outermost child's extent.
private static void BuildReaches(WellTree t, int[] parent, float[] body, float[] family, int[][] ring, int[] order,
out float[] lower, out RingGroup[] group)
{
var n = parent.Length;
var d = t._distance;
lower = new float[n];
group = new RingGroup[n];
t._support = new float[n];
t._reach = new float[n];
foreach (var i in order)
{
var low = 0f;
foreach (var k in t._kids[i]) low = Math.Max(low, t._extent[k]);
var support = GravityField.Beta * low;
foreach (var k in t._kids[i]) support = Math.Max(support, t._extent[k] + t._reach[k]);
if (family[i] > 0 && body[i] <= 0)
{
// Bodyless: its ring group stands in for a body (it has massive children, or it would carry no mass).
var inner = float.PositiveInfinity;
var mass = 0f;
foreach (var j in ring[i])
{
inner = Math.Min(inner, d[j]);
mass += family[j];
}
var ringSupport = GravityField.Beta * inner;
foreach (var j in ring[i]) ringSupport = Math.Max(ringSupport, d[j] + t._reach[j]);
low = Math.Max(low, ringSupport);
support = Math.Max(support, GravityField.Beta * low);
group[i] = new RingGroup(mass, inner, ringSupport);
}
lower[i] = low;
t._support[i] = support;
if (family[i] <= 0) continue;
var floor = GravityField.Beta * support;
var p = parent[i];
float reach;
if (p < 0) reach = t.SubzoneReach(i, floor);
else if (body[p] <= 0 && ring[p] != null && ring[p].Length >= 2 && Array.IndexOf(ring[p], i) >= 0)
reach = Math.Max(GravityField.RingReach * d[i], floor);
else reach = Math.Max(GravityField.Reach * d[i], floor);
t._reach[i] = reach > 0f ? reach : 0f;
}
}
// A root's reach. The one place the subzone shape is decided (question multi-root-zone-shape): a sole massive root is the
// whole zone (+inf, so its global term is the far field); with several, isolated-subzones bounds each root to Reach x the
// distance to the nearest other massive root, or its floor (Beta x support) if more. The superposed-inside option
// (Beta x max(pack radius, max over massive roots of |anchor - Centre| + support)) is not a change to this return: it needs
// every massive root's support, and roots are processed one at a time here, so it would take a second pass over the roots
// once their supports exist.
private float SubzoneReach(int root, float floor)
{
if (_massiveRoots.Length < 2) return float.PositiveInfinity;
var nearest = float.PositiveInfinity;
foreach (var other in _massiveRoots)
if (other != root) nearest = Math.Min(nearest, length(_anchor[other] - _anchor[root]));
return Math.Max(GravityField.Reach * nearest, floor);
}
private void EmitTerms(int[] parent, float[] body, float[] family, int[][] ring, float[] lower, RingGroup[] group,
PlanetSettings settings, float packRadius)
{
var n = parent.Length;
var terms = new List<Term>();
_termStart = new int[n];
_termCount = new int[n];
for (var i = 0; i < n; i++)
{
_termStart[i] = terms.Count;
if (family[i] <= 0 || !(_reach[i] > 0f)) continue;
var global = float.IsPositiveInfinity(_reach[i]);
var fadeFrom = global ? float.PositiveInfinity : Math.Max(_support[i], GravityField.Gamma * _reach[i]);
void Add(float strength, float softening, float inner, float outer) => terms.Add(new Term(_key[i], strength, softening, inner,
outer, float.IsPositiveInfinity(outer) ? strength * GravityField.Potential(softening, inner, outer, packRadius) : 0f));
if (body[i] > 0 && family[i] <= body[i])
{
var (c, eps) = Curves(settings, body[i]);
Add(c, eps, fadeFrom, _reach[i]);
_termCount[i] = terms.Count - _termStart[i];
continue;
}
if (body[i] > 0)
{
var (c, eps) = Curves(settings, body[i]);
Add(c, eps, lower[i], _support[i]);
}
else
{
var (c, eps) = Curves(settings, group[i].Mass);
Add(c, eps, lower[i], _support[i]);
Add(-c, eps, group[i].Inner, group[i].Support);
}
var (cf, ef) = Curves(settings, family[i]);
Add(cf, ef, fadeFrom, _reach[i]);
Add(-cf, ef, lower[i], _support[i]);
_termCount[i] = terms.Count - _termStart[i];
}
_zoneStart = terms.Count;
if (_massiveRoots.Length >= 2)
{
var total = 0f;
var span = 0f;
foreach (var r in _massiveRoots)
{
total += family[r];
span = Math.Max(span, length(_anchor[r] - Centre) + _reach[r]);
}
var (c, eps) = Curves(settings, total);
terms.Add(new Term(default, c, eps, float.PositiveInfinity, float.PositiveInfinity,
c * GravityField.Potential(eps, float.PositiveInfinity, float.PositiveInfinity, packRadius)));
if (span > 0f) terms.Add(new Term(default, -c, eps, span, GravityField.Beta * span, 0f));
}
_zoneCount = terms.Count - _zoneStart;
_terms = terms.ToArray();
}
// One node's own terms at distance r from it, toward its centre.
private float PullOf(int node, float r) => PullOf(_termStart[node], _termCount[node], r);
private float PullOf(int start, int count, float r)
{
var pull = 0f;
for (var t = start; t < start + count; t++)
{
var term = _terms[t];
pull += term.Strength * GravityField.Pull(term.Softening, term.Inner, term.Outer, r);
}
return pull;
}
/// <summary>
/// The field's force at p (toward the wells, u/s^2): every node's terms, walked from the roots; a subtree is skipped
/// the moment p is outside its reach, which is exact. `positionAt(orbit, time)` is the caller's, so this never owns time.
/// </summary>
public float2 Force(float2 p, Func<CultRecordKey, double, float2> positionAt, double time)
{
var force = float2.zero;
var toZone = Centre - p;
var rz = length(toZone);
if (_zoneCount > 0 && rz > 0f) force += toZone / rz * PullOf(_zoneStart, _zoneCount, rz);
foreach (var root in _roots) ForceFrom(root, p, positionAt, time, ref force);
return force;
}
private void ForceFrom(int node, float2 p, Func<CultRecordKey, double, float2> positionAt, double time, ref float2 force)
{
var offset = positionAt(_key[node], time) - p;
var r = length(offset);
if (r > 0f)
force += offset / r * PullOf(node, r);
foreach (var kid in _kids[node])
if (length(positionAt(_key[kid], time) - p) < _reach[kid]) ForceFrom(kid, p, positionAt, time, ref force);
}
/// <summary>The field's potential at p: every node's terms less their references, walked as Force walks.</summary>
public float Potential(float2 p, Func<CultRecordKey, double, float2> positionAt, double time)
{
var sum = 0f;
var rz = length(Centre - p);
for (var t = _zoneStart; t < _zoneStart + _zoneCount; t++)
sum += _terms[t].Strength * GravityField.Potential(_terms[t].Softening, _terms[t].Inner, _terms[t].Outer, rz) - _terms[t].Reference;
foreach (var root in _roots) sum += PotentialFrom(root, p, positionAt, time);
return sum;
}
private float PotentialFrom(int node, float2 p, Func<CultRecordKey, double, float2> positionAt, double time)
{
var r = length(positionAt(_key[node], time) - p);
var sum = 0f;
for (var t = _termStart[node]; t < _termStart[node] + _termCount[node]; t++)
{
var term = _terms[t];
sum += term.Strength * GravityField.Potential(term.Softening, term.Inner, term.Outer, r) - term.Reference;
}
foreach (var kid in _kids[node])
if (length(positionAt(_key[kid], time) - p) < _reach[kid]) sum += PotentialFrom(kid, p, positionAt, time);
return sum;
}
/// <summary>
/// The pull (toward the node, u/s^2) a rail at `distance` round `parentOrbit` is given: the node's own terms at that
/// distance. A ring member of a bodyless node adds each other member's own terms at their chord (positions at t = 0),
/// projected toward the node; a strengthless orbit on that ring (a station, a turret) takes the mean of what the ring's
/// members feel, so it keeps their period. Infinite period for a pull of 0 is the caller's.
/// </summary>
public float PullAround(CultRecordKey parentOrbit, float distance, CultRecordKey self)
{
var node = _index[parentOrbit];
var pull = PullOf(node, distance);
var members = _ring[node];
if (members == null || members.Length < 2) return pull;
if (_index.TryGetValue(self, out var me) && Array.IndexOf(members, me) >= 0) return pull + PartnersPull(members, me);
var inner = float.PositiveInfinity;
foreach (var m in members) inner = Math.Min(inner, _distance[m]);
if (Math.Abs(distance - inner) > RingWindow) return pull;
var mean = 0f;
foreach (var m in members) mean += PartnersPull(members, m) / members.Length;
return pull + mean;
}
// What the other members of a ring pull member m with, along the direction from m to the node.
private float PartnersPull(int[] members, int m)
{
var at = _distance[m] * OrbitData.Evaluate(_phase[m]);
var sum = 0f;
foreach (var q in members)
{
if (q == m) continue;
var chord = _distance[q] * OrbitData.Evaluate(_phase[q]) - at;
var r = length(chord);
if (r <= 0f || _distance[m] <= 0f) continue;
sum += PullOf(q, r) * dot(chord / r, -at / _distance[m]);
}
return sum;
}
}

2
Assets/Scripts/ServerShared/WellTree.cs.meta Normal file
View file

@ -0,0 +1,2 @@
fileFormatVersion: 2
guid: da20befcb0f606b4ab1165988b836096

131
Assets/Scripts/ServerShared/Zone.cs
View file

@ -141,6 +141,7 @@ public class Zone
Pack = pack;
_itemManager = itemManager;
Settings = settings;
_orbitPositionAt = OrbitPositionAt;
CombatSeed = galaxyZone?.Name.StableHash() ?? 1337u;
var cache = itemManager.ItemData;
@ -187,18 +188,29 @@ public class Zone
}
}
// Orbits, bodies and belts exist: every rail's period is derived once, from the field it orbits in.
foreach (var orbit in Orbits)
{
var parent = orbit.Value.Data.Parent;
if (!parent.IsSet()) continue;
if (!_childOrbits.TryGetValue(parent.Key, out var siblings))
_childOrbits[parent.Key] = siblings = new List<CultRecordKey>();
siblings.Add(orbit.Key);
}
// Orbits, bodies and belts exist: the field's terms are built once, and every rail's period is derived from them.
var wellMass = new Dictionary<CultRecordKey, float>();
foreach (var body in PlanetInstances)
if (body.Value.WellStrength > 0)
_wellOn.TryAdd(Planets[body.Key].Orbit.Key, body.Value);
wellMass.TryAdd(Planets[body.Key].Orbit.Key, body.Value.BodyData.Mass);
// A belt's rail carries its asteroids, which spread beyond the rail's own distance: its extent is the outermost.
var beltExtent = new Dictionary<CultRecordKey, float>();
foreach (var belt in AsteroidBelts)
{
var extent = 0f;
foreach (var asteroid in belt.Value.Data.Asteroids) extent = Math.Max(extent, asteroid.Distance);
beltExtent[Planets[belt.Key].Orbit.Key] = extent;
}
var wellNodes = new List<WellNode>();
foreach (var orbit in Orbits)
{
var data = orbit.Value.Data;
wellNodes.Add(new WellNode(orbit.Key, data.Parent.Key, data.Distance, data.Phase,
wellMass.TryGetValue(orbit.Key, out var mass) ? mass : 0f,
beltExtent.TryGetValue(orbit.Key, out var reach) ? Math.Max(reach, data.Distance) : data.Distance,
data.FixedPosition + OrbitData.Evaluate(frac(data.Phase)) * data.Distance)); // a root's stationary place: Offset at an infinite period
}
Wells = WellTree.Build(wellNodes, settings, pack.Radius);
foreach (var orbit in Orbits)
orbit.Value.Period = orbit.Value.Data.Parent.IsSet()
? PeriodAround(orbit.Value.Data.Parent.Key, orbit.Value.Data.Distance, orbit.Key)
@ -561,67 +573,19 @@ public class Zone
return time == _time ? GetOrbitPosition(orbitID) : RailAt(orbitID, time);
}
// orbit -> the well-bearing body on it, and orbit -> the orbits whose parent it is. Derived from the pack in the
// constructor, read only by PeriodAround and Strength.
private readonly Dictionary<CultRecordKey, Planet> _wellOn = new Dictionary<CultRecordKey, Planet>();
private readonly Dictionary<CultRecordKey, List<CultRecordKey>> _childOrbits = new Dictionary<CultRecordKey, List<CultRecordKey>>();
private readonly Dictionary<CultRecordKey, float> _strengths = new Dictionary<CultRecordKey, float>();
/// <summary>The zone's gravity: the terms of its orbit tree (WellTree owns which exist and the walk that sums them).</summary>
public WellTree Wells { get; }
// The one owner of a rail's period (Orbit.Period, AsteroidBelt.Periods): 2 pi sqrt(d / a), a being the pull the
// field gives at distance d from the node being orbited. Round a body, that is the body's own well through
// GravityField. Round a bodyless node it is the ring law: the strength of the node's orbits inside d, plus the
// ring term for the orbits on d itself (within RingWindow of it), taken as a point at the node. A bodied self
// feels half the strength of its ring partners (as a binary does); a strengthless self (a station, a turret)
// takes the period its bodied ring partners have, the mean of what each of them feels. No pull gives an
// infinite period.
private const float RingWindow = 0.1f;
private readonly Func<CultRecordKey, double, float2> _orbitPositionAt;
// The one owner of a rail's period (Orbit.Period, AsteroidBelt.Periods): 2 pi sqrt(d / a), a being the pull the field
// gives at distance d from the node being orbited (WellTree.PullAround). No pull gives an infinite period.
private float PeriodAround(CultRecordKey parentOrbit, float distance, CultRecordKey self)
{
float a;
if (_wellOn.TryGetValue(parentOrbit, out var body))
a = length(GravityField.Acceleration(body.WellStrength, body.RocheLimit, float2(distance, 0)));
else
{
var inside = 0f;
var ring = new List<(CultRecordKey Orbit, float Strength)>();
if (_childOrbits.TryGetValue(parentOrbit, out var children))
foreach (var child in children)
{
var childDistance = Orbits[child].Data.Distance;
var childStrength = Strength(child);
if (childDistance < distance - RingWindow) inside += childStrength;
else if (childDistance <= distance + RingWindow && childStrength > 0) ring.Add((child, childStrength));
}
// Whose pull self feels: its own if it has strength, else each bodied member of the ring.
var ringSum = 0f;
foreach (var member in ring) ringSum += member.Strength;
var selfBodied = Strength(self) > 0;
var feelings = 0f;
var feelers = 0;
foreach (var member in ring)
if (!selfBodied || member.Orbit.Equals(self))
{
feelings += (ringSum - member.Strength) / 2;
feelers++;
}
var strength = inside + (feelers > 0 ? feelings / feelers : 0f);
a = length(GravityField.Acceleration(strength, 0f, float2(distance, 0)));
}
var a = Wells.PullAround(parentOrbit, distance, self);
return a > 0 && distance > 0 ? 2 * PI * sqrt(distance / a) : float.PositiveInfinity;
}
// The well strength an orbit carries as a node: its own body's plus every orbit round it, recursively.
private float Strength(CultRecordKey orbit)
{
if (_strengths.TryGetValue(orbit, out var known)) return known;
var strength = _wellOn.TryGetValue(orbit, out var body) ? body.WellStrength : 0f;
if (_childOrbits.TryGetValue(orbit, out var children))
foreach (var child in children)
strength += Strength(child);
return _strengths[orbit] = strength;
}
public float2 GetOrbitVelocity(CultRecordKey orbit)
{
if (Orbits.ContainsKey(orbit))
@ -734,10 +698,6 @@ public class Zone
return security;
}
// The potential of a body's well at the zone's edge, subtracted so every well is zero there: the one owner of the
// reference. GetHeight sums with it and the renderer's lowering (GravityLowering) hands it to the shader unchanged.
public float WellReference(Planet body) => GravityField.Potential(body.WellStrength, body.RocheLimit, Pack.Radius);
// The extent of the populated area, from the zone's origin: owned here; presenters read it, nothing authors it. A static
// bound on where anything can be, so it holds at every time and reads nothing the pack's Radius decides: for each body,
// its orbit chain's reach plus its Roche limit; for each belt, the chain's reach to the orbit its asteroids circle (the
@ -769,19 +729,9 @@ public class Zone
public float GetHeight(float2 position) => GetHeight(position, _time);
// The field at a time: each well sits where its body's rail puts it then. The sim reads it at _time, a presenter at
// The field at a time: each term sits where its node's rail puts it then. The sim reads it at _time, a presenter at
// the drawn time, so a ship rides the same slope the drawn planet's well makes.
public float GetHeight(float2 position, double time)
{
var result = 0f;
foreach (var body in PlanetInstances.Values)
{
var offset = position - OrbitPositionAt(body.BodyData.Orbit.Key, time);
result += GravityField.Potential(body.WellStrength, body.RocheLimit, length(offset)) - WellReference(body);
}
return result / Settings.HeightExaggeration;
}
public float GetHeight(float2 position, double time) => Wells.Potential(position, _orbitPositionAt, time) / Settings.HeightExaggeration;
public float GetLight(float2 position)
{
@ -800,16 +750,10 @@ public class Zone
internal float SchematicCellSize => _itemManager.GameplaySettings.SchematicCellSize;
internal float FieldCrossSection(AsteroidBeltData belt) => _itemManager.ItemData.Get(belt.Kind)?.CrossSection ?? 0f;
// The one source of gravitational acceleration (u/s^2): the sum of every body's well (GravityField). Ships,
// floating items, clouds and mines all read it, and so does Flow. It is the gradient of the potential GetHeight
// draws, so the two cannot disagree.
public float2 GetForce(float2 position)
{
var force = float2.zero;
foreach (var body in PlanetInstances.Values)
force += GravityField.Acceleration(body.WellStrength, body.RocheLimit, position - body.Orbit.Position);
return force;
}
// The one source of gravitational acceleration (u/s^2): the walk of the zone's terms (WellTree). Ships, floating items,
// clouds and mines all read it, and so does Flow. It is the gradient of the potential GetHeight draws, so the two
// cannot disagree.
public float2 GetForce(float2 position) => Wells.Force(position, _orbitPositionAt, _time);
// The global flow velocity at `position` in zone time, in u/s (FlowField is the law). The slope term reads the
// downhill direction, which is the force.
@ -861,12 +805,7 @@ public class Planet
public void CalculateProperties()
{
BodyRadius = Settings.BodyRadius.Evaluate(BodyData.Mass) * BodyData.BodyRadiusMultiplier;
WellStrength = Settings.WellStrength.Evaluate(BodyData.Mass);
RocheLimit = Settings.RocheLimit.Evaluate(BodyData.Mass);
if (!(WellStrength >= 0f) || float.IsInfinity(WellStrength))
throw new InvalidOperationException("A body's well strength must be finite and not negative (PlanetSettings.WellStrength).");
if (WellStrength > 0 && !(RocheLimit > 0f && !float.IsInfinity(RocheLimit)))
throw new InvalidOperationException("A body's Roche limit, which is its core radius, must be positive and finite (PlanetSettings.RocheLimit).");
(WellStrength, RocheLimit) = WellTree.Curves(Settings, BodyData.Mass);
}
}

51
Assets/Scripts/Zone Display/GravityLowering.cs
View file

@ -3,8 +3,8 @@ using UnityEngine;
/// <summary>
/// The renderer's lowering of the zone's potential (aetheria-gravity, cut grid-lowering). Zone.GetHeight owns the drawn
/// height: this hands its inputs (each well's position, strength, core radius and reference, and the exaggeration) to
/// Assets/Shaders/Gravity.cginc, which sums them per texel, and makes the one quad that carries the sum into every
/// height: this hands its inputs (each term's position, strength, softening, band edges and reference, and the
/// exaggeration) to Assets/Shaders/Gravity.cginc, which sums them per texel, and makes the one quad that carries the sum into every
/// camera that culls layer 9. Nothing here decides a height.
/// </summary>
public static class GravityLowering
@ -14,30 +14,32 @@ public static class GravityLowering
/// <summary>The quad's side in u: larger than any zone (10000) plus any view (4096); texels outside a camera cost nothing.</summary>
public const float FieldExtent = 100000f;
// Mirrors the well record of Gravity.cginc: float2 position, strength, coreRadius, reference.
// Mirrors the term record of Gravity.cginc: float2 position, strength, softening, inner, outer, reference.
[StructLayout(LayoutKind.Sequential)]
struct WellRecord
struct TermRecord
{
public Vector2 Position;
public float Strength;
public float CoreRadius;
public float Softening;
public float Inner;
public float Outer;
public float Reference;
}
static readonly int Stride = Marshal.SizeOf<WellRecord>();
static readonly int Stride = Marshal.SizeOf<TermRecord>();
static readonly int WellsId = Shader.PropertyToID("_GravityWells");
static readonly int CountId = Shader.PropertyToID("_GravityWellCount");
static readonly int ExaggerationId = Shader.PropertyToID("_HeightExaggeration");
static WellRecord[] _records = new WellRecord[0];
static TermRecord[] _records = new TermRecord[0];
/// <summary>
/// Fills one record per body of the zone (a null zone has none) and sets the shader globals. The buffer is the
/// caller's, (re)allocated to the body count (at least one, so the global is always bound). Wells are lowered at the
/// drawn time (DrawAhead.Time), the same T the planets and a ship's ride height use.
/// Fills one record per term of the zone's WellTree (a null zone has none) and sets the shader globals. The buffer is
/// the caller's, (re)allocated to the term count (at least one, so the global is always bound). Terms are lowered at
/// their node's rail at the drawn time (DrawAhead.Time), the same T the planets and a ship's ride height use.
/// </summary>
public static void Lower(Zone zone, float lead, ref ComputeBuffer wells)
{
var count = zone?.PlanetInstances.Count ?? 0;
var count = zone?.Wells.Terms.Count ?? 0;
var capacity = Mathf.Max(count, 1);
if (wells == null || wells.count != capacity)
{
@ -45,20 +47,21 @@ public static class GravityLowering
wells = new ComputeBuffer(capacity, Stride);
}
if (_records.Length != capacity) _records = new WellRecord[capacity];
var i = 0;
if (zone != null)
foreach (var body in zone.PlanetInstances.Values)
if (_records.Length != capacity) _records = new TermRecord[capacity];
for (var t = 0; t < count; t++)
{
var term = zone.Wells.Terms[t];
var at = term.Node.IsSet() ? DrawAhead.Orbit(zone, term.Node, lead) : zone.Wells.Centre; // an unset node is the zone node, stationary at Centre
_records[t] = new TermRecord
{
var well = DrawAhead.Well(zone, body, lead);
_records[i++] = new WellRecord
{
Position = new Vector2(well.Position.x, well.Position.y),
Strength = well.Strength,
CoreRadius = well.CoreRadius,
Reference = well.Reference
};
}
Position = new Vector2(at.x, at.y),
Strength = term.Strength,
Softening = term.Softening,
Inner = term.Inner,
Outer = term.Outer,
Reference = term.Reference
};
}
wells.SetData(_records);
Shader.SetGlobalBuffer(WellsId, wells);

79
Assets/Shaders/Gravity.cginc
View file

@ -4,38 +4,91 @@
// GENERATED by ServerShared/GravityShader.cs from GravityField; do not edit (GravityShaderTests fails on any difference).
// The lowering of the zone's potential (aetheria-gravity, cut grid-lowering). The law is ServerShared's GravityField
// and the height is Zone.GetHeight; this file decides nothing. GravityLowering.Lower fills the buffer and the globals
// from the zone (positions from Orbit.Position, reference from Zone.WellReference) each frame.
// from the zone's WellTree (positions from the nodes' rails at the drawn time) each frame.
// GravityFieldParityVerify renders it and compares every texel with -Zone.GetHeight.
struct GravityWell
struct GravityTerm
{
float2 position;
float strength;
float coreRadius;
float softening;
float inner;
float outer;
float reference;
};
StructuredBuffer<GravityWell> _GravityWells;
StructuredBuffer<GravityTerm> _GravityWells;
int _GravityWellCount;
float _HeightExaggeration;
// GravityField.Potential: c is the well's strength, R its core radius, r the distance.
float gravityPotential(float c, float R, float r)
// GravityField.OneMinusS and GravityField.S: the band profile B(t) = 1 - S(t) and the order-4 smoothstep, lowest first.
static const float GravityOneMinusS[8] = { 1, -0, -0, -0, -35, 84, -70, 20 };
static const float GravityS[8] = { 0, 0, 0, 0, 35, -84, 70, -20 };
// GravityField.BandClosed: the band's potential at t for tau <= 0.4, per unit strength.
float gravityBandClosed(float tau, float t)
{
return c == 0 ? 0 :
r < R
? c * r * r / (2 * R * R)
: c * (log(r / R) + 0.5);
float q[7];
q[6] = GravityOneMinusS[7];
[unroll] for (int i = 6; i > 0; i--) q[i - 1] = GravityOneMinusS[i] - tau * q[i];
float remainder = GravityOneMinusS[0] - tau * q[0];
float whole = 0;
float part = 0;
[unroll] for (int j = 6; j >= 0; j--)
{
whole = whole + q[j] / (j + 1);
part = (part + q[j] / (j + 1)) * t;
}
return 0.5 * (remainder * log((t + tau) / (1 + tau)) - (whole - part));
}
// Zone.GetHeight at world xz: each well's potential minus its own reference, summed, divided by the exaggeration.
// GravityField.BandSeries: the band's potential at t for tau > 0.4, a series in (1 - t) / (1 + tau).
float gravityBandSeries(float tau, float t)
{
float y = 1 - t;
float yj[8];
yj[0] = 1;
[unroll] for (int j = 1; j < 8; j++) yj[j] = yj[j - 1] * y;
float w = 1 / (1 + tau);
float wk = w;
float yk = y;
float acc = 0;
for (int k = 0; k < 40; k++)
{
float inner = 0;
[unroll] for (int m = 4; m < 8; m++) inner += GravityS[m] * yj[m] / (m + k + 1);
acc += wk * yk * inner;
wk *= w;
yk *= y;
}
return -0.5 * acc;
}
// GravityField.Potential: one band term, per unit strength; the global term when b is infinite.
float gravityTermPotential(float softening, float a, float b, float r)
{
float s = r * r;
float e2 = softening * softening;
if (isinf(b)) return 0.5 * log(s + e2);
float a2 = a * a;
float b2 = b * b;
if (s >= b2) return 0;
float tau = (a2 + e2) / (b2 - a2);
if (s <= a2) return 0.5 * log((s + e2) / (a2 + e2)) + (tau <= 0.4 ? gravityBandClosed(tau, 0) : gravityBandSeries(tau, 0));
float t = (s - a2) / (b2 - a2);
return tau <= 0.4 ? gravityBandClosed(tau, t) : gravityBandSeries(tau, t);
}
// Zone.GetHeight at world xz: each term's potential minus its own reference, summed, divided by the exaggeration.
float gravityHeight(float2 p)
{
float result = 0;
for (int i = 0; i < _GravityWellCount; i++)
{
GravityWell well = _GravityWells[i];
result += gravityPotential(well.strength, well.coreRadius, length(p - well.position)) - well.reference;
GravityTerm term = _GravityWells[i];
float r = length(p - term.position);
if (r < term.outer)
result += term.strength * gravityTermPotential(term.softening, term.inner, term.outer, r) - term.reference;
}
return result / _HeightExaggeration;
}

46
tests/Aetheria.Shared.Tests/DrawnPoseTests.cs
View file

@ -16,7 +16,8 @@ using float3 = CultMath.float3;
// Every drawn pose is evaluated at one drawn time, T = Zone.ExactTime + Clock.Lead (DrawAhead.Time). The operator's report:
// coasting beside a planet and looking at it, the planet jumped back and forth every frame, because the planet was drawn at
// the last step's orbit while the ship was drawn a lead ahead of its last step. Each rail below is spelled here from the
// field law outside the core (a = c / d, so T = 2 pi d / sqrt(c)), never by calling the production rail; the zone is built
// field law inside a node's outermost child (the node pulls as its body alone: a = c d / (d^2 + eps^2), T = 2 pi sqrt(d / a)),
// never by calling the production rail; the zone is built
// by the production constructor and stepped by Zone.Update, and the clock is the production SimClock.
public sealed class DrawnPoseTests : IDisposable
{
@ -124,9 +125,12 @@ public sealed class DrawnPoseTests : IDisposable
// ---- The rails, spelled from the field law ----
// T = 2 pi sqrt(d / a), a the Plummer pull of a body alone (softening 2) at d.
private static float PeriodAt(float distance, float strength) => 2 * PI * sqrt(distance / (strength * distance / (distance * distance + 4f)));
private static float2 Around(float2 center, double time, float distance, float phase, float strength)
{
var period = 2 * PI * distance / sqrt(strength);
var period = PeriodAt(distance, strength);
var wrapped = (float) frac(time / period);
return center + OrbitData.Evaluate(frac(wrapped + phase)) * distance;
}
@ -135,13 +139,14 @@ public sealed class DrawnPoseTests : IDisposable
private static float2 MoonAt(double time) => Around(PlanetAt(time), time, MoonDistance, MoonPhase, PlanetMass);
private static float2 StationAt(double time) => Around(PlanetAt(time), time, StationDistance, StationPhase, PlanetMass);
// The field at a point at a time, with each well where its rail puts it then.
// The field at a point at a time: the zone's terms (WellTree owns which exist), each at the position this test spells
// for its node's rail at that time.
private static float FieldAt(World w, float2 point, double time)
{
var wells = new (float2 At, float Strength)[] { (SunCentre, SunMass), (PlanetAt(time), PlanetMass), (MoonAt(time), MoonMass) };
float2 At(CultRecordKey node) => node.Equals(w.PlanetOrbit) ? PlanetAt(time) : node.Equals(w.MoonOrbit) ? MoonAt(time) : SunCentre;
var sum = 0f;
foreach (var (at, strength) in wells)
sum += GravityField.Potential(strength, 2f, length(point - at)) - GravityField.Potential(strength, 2f, w.Zone.Pack.Radius);
foreach (var term in w.Zone.Wells.Terms)
sum += term.Strength * GravityField.Potential(term.Softening, term.Inner, term.Outer, length(point - At(term.Node))) - term.Reference;
return sum / w.HeightExaggeration;
}
@ -332,7 +337,7 @@ public sealed class DrawnPoseTests : IDisposable
for (var i = 0; i < w.Asteroids.Length; i++)
{
var asteroid = w.Asteroids[i];
var period = 2 * PI * asteroid.Distance / sqrt(PlanetMass);
var period = PeriodAt(asteroid.Distance, PlanetMass);
var expected = OrbitData.Evaluate((float) frac(time / period + asteroid.Phase)) * asteroid.Distance + parent;
AssertNear(expected, drawn[i].xy, 1e-3f, $"chunk {i} drawn");
Assert.Equal(w.Zone.ChunkPose(w.Belt, i, time), drawn[i]);
@ -342,34 +347,27 @@ public sealed class DrawnPoseTests : IDisposable
}
}
// The well the shader takes is its body's rail at the drawn time: the lowering copies DrawAhead.Well and decides
// nothing, so the drawn well and the drawn planet share one pose at a non-zero lead, and the strength, core and
// reference are the field law's (spelled here from the masses, the core 2 and the zone's edge), not the step's.
// The term the shader takes sits on its node's rail at the drawn time: the lowering copies DrawAhead.Orbit of the term's
// node and decides nothing, so a drawn term and the drawn planet share one pose at a non-zero lead.
[Fact]
public void AWellIsLoweredAtTheDrawnTime()
public void ATermIsLoweredAtItsNodesDrawnRail()
{
var w = Build();
Step(w, 1f / 60f, 803);
Assert.Contains(w.Zone.Wells.Terms, t => t.Node.Equals(w.PlanetOrbit));
Assert.Contains(w.Zone.Wells.Terms, t => t.Node.Equals(w.MoonOrbit));
foreach (var lead in new[] { 0f, .004f, .0166f })
{
var time = w.Zone.ExactTime + lead;
foreach (var body in w.Zone.PlanetInstances.Values)
foreach (var term in w.Zone.Wells.Terms)
{
var key = body.BodyData.Orbit.Key;
var (at, strength) = key.Equals(w.PlanetOrbit) ? (PlanetAt(time), PlanetMass)
: key.Equals(w.MoonOrbit) ? (MoonAt(time), MoonMass)
: (SunCentre, SunMass);
var well = DrawAhead.Well(w.Zone, body, lead);
AssertNear(at, well.Position, 1e-3f, $"well of the {strength} body at lead {lead}");
Assert.Equal(strength, well.Strength, 3);
Assert.Equal(2f, well.CoreRadius, 3);
Assert.Equal(GravityField.Potential(strength, 2f, w.Zone.Pack.Radius), well.Reference, 3);
var at = term.Node.Equals(w.PlanetOrbit) ? PlanetAt(time) : term.Node.Equals(w.MoonOrbit) ? MoonAt(time) : SunCentre;
AssertNear(at, DrawAhead.Orbit(w.Zone, term.Node, lead), 1e-3f, $"term of a body at lead {lead}");
}
}
var planet = w.Zone.PlanetInstances.Values.Single(b => b.BodyData.Orbit.Key.Equals(w.PlanetOrbit));
Assert.True(length(DrawAhead.Well(w.Zone, planet, .0166f).Position - DrawAhead.Well(w.Zone, planet, 0f).Position) > .1f, "the well moves over a lead");
AssertNear(w.Zone.Orbits[w.PlanetOrbit].Position, DrawAhead.Well(w.Zone, planet, 0f).Position, 0f, "at lead 0 the well is the step's memo");
Assert.True(length(DrawAhead.Orbit(w.Zone, w.PlanetOrbit, .0166f) - DrawAhead.Orbit(w.Zone, w.PlanetOrbit, 0f)) > .1f, "the term moves over a lead");
AssertNear(w.Zone.Orbits[w.PlanetOrbit].Position, DrawAhead.Orbit(w.Zone, w.PlanetOrbit, 0f), 0f, "at lead 0 the term is the step's memo");
}
// A HUD marker round the ship the camera follows (an incoming lock, a contact's bearing, a ping's flash) is anchored on

6
tests/Aetheria.Shared.Tests/FlowFieldTests.cs
View file

@ -124,8 +124,8 @@ public sealed class FlowFieldTests : IDisposable
}
}
// Zone.Flow's slope term reads the owner's superposed force, not any one body's: two equal wells, probes off the
// axis between them where the nearest body's direction and the sum's differ.
// Zone.Flow's slope term reads the owner's force at the point: two equal wells 600 apart (two subzones, each felt within
// its reach of 300), probes inside each off the axis, where the direction is not along the axis.
[Fact]
public void TheSlopeTermFollowsTheSuperposedForce()
{
@ -135,7 +135,7 @@ public sealed class FlowFieldTests : IDisposable
WellFixture.Dig(_items.ItemData, pack, planets, 300f, 100f, new float2(300f, 0f));
var zone = new Zone(_items, planets, pack, new GalaxyZone { Name = "Two Wells", Owner = null }, null);
var checkedPoints = 0;
foreach (var at in new[] { float2(100f, 400f), float2(-50f, -300f), float2(250f, 220f), float2(-420f, 90f) })
foreach (var at in new[] { float2(-200f, 150f), float2(-350f, -120f), float2(250f, 220f), float2(420f, -90f) })
{
var flow = zone.Flow(float3(at.x, 3f, at.y));
var expected = 1.5f * normalize(zone.GetForce(at));

4
tests/Aetheria.Shared.Tests/GravityShaderTests.cs
View file

@ -4,7 +4,7 @@ using Xunit;
// potential-owns-gravity (aetheria-gravity, cut grid-lowering): Assets/Shaders/Gravity.cginc is the output of
// GravityShader.Source(), which reads GravityField's constants and writes the one lowering of the potential, the
// zero-strength guard, each well's reference and gravityHeight's sum. The committed file must equal it, character for
// band term's closed form and series, each term's reference and gravityHeight's sum. The committed file must equal it, character for
// character (line ends aside): any added statement, respelling, flipped sign, swapped argument or macro fails here.
// To change the shader, change GravityShader and run once with AETHERIA_REGENERATE_GRAVITY_CGINC=1.
public sealed class GravityShaderTests
@ -17,7 +17,7 @@ public sealed class GravityShaderTests
File.WriteAllText(path, GravityShader.Source());
var generated = GravityShader.Source();
// The law was generated, not an empty template agreeing with an empty file.
Assert.Contains("log(r / R) + 0.5", generated);
Assert.Contains("static const float GravityS[8] = { 0, 0, 0, 0, 35, -84, 70, -20 };", generated);
Assert.Equal(generated, File.ReadAllText(path).Replace("\r\n", "\n"));
}
}

8
tests/Aetheria.Shared.Tests/LockOnScenarioTests.cs
View file

@ -35,6 +35,7 @@ public sealed partial class RunStartTests
var lockedAt = -1;
var bearingDot = 0f;
var stagedDrift = 0f;
var strengths = new List<float>();
for (var tick = 0; tick < 600 && (lockedAt < 0 || tick < lockedAt + 40); tick++)
{
@ -44,13 +45,16 @@ public sealed partial class RunStartTests
if (lockedAt < 0)
{
lockedAt = tick;
bearingDot = dot(normalize(paint.Bearing), toPlayer * -1f);
// The warning points at the hostile where it is at the lock; the hostile has left the staged line by then.
var live = normalize(hostile.Position.xz - player.Position.xz);
bearingDot = dot(normalize(paint.Bearing), live);
stagedDrift = degrees(acos(clamp(dot(live, -toPlayer), -1f, 1f)));
}
strengths.Add(paint.Strength);
}
}
Assert.True(lockedAt >= 0, "the hostile never locked the player within 600 ticks of .1 s");
Assert.True(bearingDot >= cos(radians(5f)), "the warning points at the hostile");
Assert.True(bearingDot >= cos(radians(.5f)), $"the warning points at the hostile where it is (the hostile is {stagedDrift} deg off its staged line at the lock)");
Assert.True(strengths.First() < strengths.Last(), $"the lock tightens: Strength went {strengths.First()} to {strengths.Last()}");
}
}

4
tests/Aetheria.Shared.Tests/MiningCut3Tests.cs
View file

@ -134,7 +134,7 @@ public sealed partial class MiningCut3Tests : IDisposable
var items = new ItemManager(cache, new ProvenanceLedger(), GameSettings(), _ => { });
var sunOrbit = cache.Upsert(new OrbitData { FixedPosition = sunAt }).Key;
var sun = cache.Upsert(new SunData { Orbit = new CultRecordRef<OrbitData>(sunOrbit) }).Key;
var sun = cache.Upsert(new SunData { Orbit = new CultRecordRef<OrbitData>(sunOrbit), Mass = 1f }).Key;
var pack = new ZonePack
{
Orbits = { new CultRecordRef<OrbitData>(sunOrbit) },
@ -147,7 +147,7 @@ public sealed partial class MiningCut3Tests : IDisposable
{
var orbit = cache.Upsert(new OrbitData { FixedPosition = at }).Key;
pack.Orbits.Add(new CultRecordRef<OrbitData>(orbit));
pack.Planets.Add(new CultRecordRef<BodyData>(cache.Upsert(new SunData { Orbit = new CultRecordRef<OrbitData>(orbit) }).Key));
pack.Planets.Add(new CultRecordRef<BodyData>(cache.Upsert(new SunData { Orbit = new CultRecordRef<OrbitData>(orbit), Mass = 1f }).Key));
}
var scene = new Scene { Cache = cache, Items = items, Hull = hull, Eye = eye, DimEye = dimEye, Pack = pack };
scene.Kinds.AddRange(kinds ?? new FieldKindData[0]);

146
tests/Aetheria.Shared.Tests/OrbitsFromFieldTests.cs
View file

@ -9,76 +9,31 @@ using static CultMath.math;
using float2 = CultMath.float2;
using float3 = CultMath.float3;
// Aetheria gravity, orbits-from-field: every rail's period is derived from the field it orbits in. Round a body it is
// that body's own well; round a bodyless node it is the ring law (the strength inside the rail, plus half of what shares
// it). The zones here are hand-built packs read by the production constructor, with the shipped well curves
// (c = 45.08 m^.2, R = 7.32 m^.25). Each expected period is spelled here from the field law outside the core, T = 2 pi d / sqrt(c),
// with c and R read from the zone's own bodies and every distance outside the core.
// Aetheria gravity, orbits-from-field: every rail's period is derived from the field it orbits in (WellTree.PullAround).
// Round a body that is the body's own pull inside its outermost child; round a bodyless node a ring member also feels each
// other member's own pull at their chord, projected toward the node. The zones here are hand-built packs read by the
// production constructor, with the shipped well curves (c = 45.08 m^.2, eps = 7.32 m^.25). Each expected period is
// spelled here as T = 2 pi sqrt(d / a), a composed from GravityField.Pull and c and eps read from the zone's own bodies.
public sealed class OrbitsFromFieldTests : IDisposable
{
private readonly string _root = Path.Combine(Path.GetTempPath(), "aetheria-orbitsfield-" + Guid.NewGuid().ToString("N"));
private readonly List<CultCache> _caches = new List<CultCache>();
private readonly ZoneWorlds _worlds = new ZoneWorlds();
public OrbitsFromFieldTests() => Directory.CreateDirectory(_root);
public void Dispose() => _worlds.Dispose();
public void Dispose()
{
foreach (var cache in _caches) cache.Dispose();
Directory.Delete(_root, true);
}
private static GameplaySettings GameSettings() => new GameplaySettings
{
DefaultEntitySettings = new EntitySettings(),
Tiers = new[] { new RarityTier { Name = "Common", Quality = .5f, Rarity = 0, Color = new float3(1, 1, 1) } },
QualityPriceModifier = new ExponentialLerp()
};
private static PlanetSettings Curves() => new PlanetSettings
{
WellStrength = new ExponentialCurve { Exponent = .2f, Multiplier = 45.08f },
RocheLimit = new ExponentialCurve { Exponent = .25f, Multiplier = 7.32f },
BodyRadius = new ExponentialCurve { Exponent = .25f, Multiplier = 3f },
LightRadius = new ExponentialCurve()
};
private sealed class World
{
public CultCache Cache;
public ItemManager Items;
public ZonePack Pack = new ZonePack { Radius = 20000f, Mass = 10000f };
public CultRecordKey Orbit(CultRecordKey parent, float distance, float phase = 0f)
{
var data = new OrbitData { Distance = distance, Phase = phase };
if (parent.IsSet()) data.Parent = new CultRecordRef<OrbitData>(parent);
var key = Cache.Upsert(data).Key;
Pack.Orbits.Add(new CultRecordRef<OrbitData>(key));
return key;
}
public CultRecordKey Body(CultRecordKey orbit, float mass)
{
var key = Cache.Upsert(new PlanetData { Name = "Body", Orbit = new CultRecordRef<OrbitData>(orbit), Mass = mass }).Key;
Pack.Planets.Add(new CultRecordRef<BodyData>(key));
return key;
}
public Zone Build() => new Zone(Items, Curves(), Pack, new GalaxyZone { Name = "OrbitsFromField", Owner = null }, null);
}
private World NewWorld()
{
var dir = Path.Combine(_root, $"world{_caches.Count}");
Directory.CreateDirectory(dir);
var cache = AetheriaStores.Open(Path.Combine(dir, "Aetheria.cc"), Path.Combine(dir, "run.cc"), catalogWritable: true);
_caches.Add(cache);
return new World { Cache = cache, Items = new ItemManager(cache, new ProvenanceLedger(), GameSettings(), _ => { }) };
}
private ZoneWorlds.World NewWorld() => _worlds.New();
private static PlanetSettings Curves() => ZoneWorlds.Curves();
private static float C(Zone zone, CultRecordKey body) => zone.PlanetInstances[body].WellStrength;
private static float Core(Zone zone, CultRecordKey body) => zone.PlanetInstances[body].RocheLimit;
private static float Outside(float distance, float c) => 2 * PI * distance / sqrt(c);
// The pull of a body's well as one body, at distance d: c times the band term's pull with no band (GravityField.Pull).
private static float Pull(float c, float roche, float d) => c * GravityField.Pull(roche, float.PositiveInfinity, float.PositiveInfinity, d);
private static float Pull(Zone zone, CultRecordKey body, float d) => Pull(C(zone, body), Core(zone, body), d);
// A mass the zone has no body of (a family's total): its well from the authored curves, as the law reads them.
private static float FamilyPull(float mass, float d) => Pull(Curves().WellStrength.Evaluate(mass), Curves().RocheLimit.Evaluate(mass), d);
private static float PeriodUnder(float pull, float distance) => 2 * PI * sqrt(distance / pull);
private static void AssertPeriod(float expected, float actual, string what) =>
Assert.True(Math.Abs(actual / expected - 1f) <= 1e-4f, $"{what}: expected {expected}, got {actual}");
@ -95,13 +50,12 @@ public sealed class OrbitsFromFieldTests : IDisposable
w.Body(moonOrbit, 1f);
var zone = w.Build();
Assert.True(Core(zone, planet) < 150f && Core(zone, sun) < 3000f, "fixture: both rails run outside the core");
AssertPeriod(Outside(150f, C(zone, planet)), zone.Orbits[moonOrbit].Period, "the moon, in its planet's field");
AssertPeriod(Outside(3000f, C(zone, sun)), zone.Orbits[planetOrbit].Period, "the planet, in its sun's field");
AssertPeriod(PeriodUnder(Pull(zone, planet, 150f), 150f), zone.Orbits[moonOrbit].Period, "the moon, in its planet's field");
AssertPeriod(PeriodUnder(Pull(zone, sun, 3000f), 3000f), zone.Orbits[planetOrbit].Period, "the planet, in its sun's field");
}
[Fact]
public void AChildInsideItsParentsCoreOrbitsAsASolidBody()
public void AChildInsideItsParentsRocheLimitOrbitsInThePlummerPull()
{
var w = NewWorld();
var sunOrbit = w.Orbit(default, 0f);
@ -113,9 +67,8 @@ public sealed class OrbitsFromFieldTests : IDisposable
var far = w.Orbit(planetOrbit, coreRadius * .7f);
var zone = w.Build();
var solid = 2 * PI * coreRadius / sqrt(C(zone, planet)); // a = c d / R^2, so T = 2 pi R / sqrt(c) at every d inside
AssertPeriod(solid, zone.Orbits[near].Period, "the inner child");
AssertPeriod(solid, zone.Orbits[far].Period, "the outer child");
AssertPeriod(PeriodUnder(Pull(zone, planet, coreRadius * .3f), coreRadius * .3f), zone.Orbits[near].Period, "the inner child");
AssertPeriod(PeriodUnder(Pull(zone, planet, coreRadius * .7f), coreRadius * .7f), zone.Orbits[far].Period, "the outer child");
}
[Fact]
@ -130,11 +83,13 @@ public sealed class OrbitsFromFieldTests : IDisposable
w.Body(outerOrbit, 100f);
var zone = w.Build();
AssertPeriod(Outside(3000f, C(zone, sun)), zone.Orbits[outerOrbit].Period, "the outer planet, with a heavy planet inside it");
AssertPeriod(PeriodUnder(Pull(zone, sun, 3000f), 3000f), zone.Orbits[outerOrbit].Period, "the outer planet, with a heavy planet inside it");
}
// A binary on a bodyless node: a member feels its partner's own pull at their chord (twice the ring's radius), toward
// the node, and nothing of the node itself.
[Fact]
public void ABinaryOrbitsUnderHalfItsPartnersPull()
public void ABinaryMemberFeelsItsPartnerAtTheChord()
{
var w = NewWorld();
var node = w.Orbit(default, 0f);
@ -145,10 +100,11 @@ public sealed class OrbitsFromFieldTests : IDisposable
var zone = w.Build();
Assert.True(C(zone, bodyHeavy) > 1.5f * C(zone, bodyA), "fixture: unequal members");
AssertPeriod(Outside(400f, C(zone, bodyHeavy) / 2), zone.Orbits[a].Period, "the light member");
AssertPeriod(Outside(400f, C(zone, bodyA) / 2), zone.Orbits[heavy].Period, "the heavy member");
AssertPeriod(PeriodUnder(Pull(zone, bodyHeavy, 800f), 400f), zone.Orbits[a].Period, "the light member");
AssertPeriod(PeriodUnder(Pull(zone, bodyA, 800f), 400f), zone.Orbits[heavy].Period, "the heavy member");
}
// Conservation: the outer body feels one well of the whole mass inside it, c(M_pair + moon), not the sum of strengths.
[Fact]
public void AnOuterBodyFeelsTheWholePairAndItsMoons()
{
@ -157,15 +113,14 @@ public sealed class OrbitsFromFieldTests : IDisposable
var m1 = w.Orbit(node, 400f, 0f);
var m2 = w.Orbit(node, 400f, .5f);
var moon = w.Orbit(m1, 40f);
var body1 = w.Body(m1, 100f);
var body2 = w.Body(m2, 300f);
var moonBody = w.Body(moon, 5f);
w.Body(m1, 100f);
w.Body(m2, 300f);
w.Body(moon, 5f);
var outer = w.Orbit(node, 3000f);
w.Body(outer, 50f);
var zone = w.Build();
var whole = C(zone, body1) + C(zone, moonBody) + C(zone, body2);
AssertPeriod(Outside(3000f, whole), zone.Orbits[outer].Period, "the outer body, round the pair and the first one's moon");
AssertPeriod(PeriodUnder(FamilyPull(405f, 3000f), 3000f), zone.Orbits[outer].Period, "the outer body, round the pair and the first one's moon");
}
[Fact]
@ -191,7 +146,7 @@ public sealed class OrbitsFromFieldTests : IDisposable
w.Body(planet2, 100f);
var companion2 = w.Orbit(node, 500f, 1f / 6f);
zone = w.Build();
AssertPeriod(Outside(500f, C(zone, inner)), zone.Orbits[planet2].Period, "the planet");
AssertPeriod(PeriodUnder(Pull(zone, inner, 500f), 500f), zone.Orbits[planet2].Period, "the planet");
Assert.Equal(zone.Orbits[planet2].Period, zone.Orbits[companion2].Period);
}
@ -213,10 +168,10 @@ public sealed class OrbitsFromFieldTests : IDisposable
var zone = w.Build();
Assert.True(C(zone, heavyBody) > 1.5f * C(zone, lightBody), "fixture: unequal members");
AssertPeriod(Outside(400f + offset, (C(zone, lightBody) + C(zone, heavyBody)) / 4), zone.Orbits[companion].Period,
AssertPeriod(PeriodUnder((Pull(zone, lightBody, 800f) + Pull(zone, heavyBody, 800f)) / 2, 400f + offset), zone.Orbits[companion].Period,
"the companion, between its two members");
// An equal pair: the companion has exactly the members' period.
// An equal pair: the companion has the members' period.
w = NewWorld();
node = w.Orbit(default, 0f);
var a = w.Orbit(node, 400f, 0f);
@ -226,12 +181,13 @@ public sealed class OrbitsFromFieldTests : IDisposable
companion = w.Orbit(node, 400f, .25f);
zone = w.Build();
Assert.True(float.IsFinite(zone.Orbits[a].Period));
Assert.Equal(zone.Orbits[a].Period, zone.Orbits[companion].Period);
Assert.Equal(zone.Orbits[b].Period, zone.Orbits[companion].Period);
AssertPeriod(zone.Orbits[a].Period, zone.Orbits[companion].Period, "the companion, against the first member");
AssertPeriod(zone.Orbits[b].Period, zone.Orbits[companion].Period, "the companion, against the second member");
}
// The ring is the 0.1 u window about a rail's distance: a member just inside it is a ring partner (pulls half), one just
// outside is inside or beyond the rail (pulls whole, or not at all).
// The ring is the 0.1 u window about the innermost massive distance: a member just inside it is a ring partner (felt at
// the chord), one just outside it is not a partner and the rail has the node's own terms alone, which give a
// ring's own member next to no pull.
[Theory]
[InlineData(-.09f, true)]
[InlineData(.09f, true)]
@ -246,8 +202,11 @@ public sealed class OrbitsFromFieldTests : IDisposable
var partnerBody = w.Body(partner, 2000f);
var zone = w.Build();
var c = C(zone, partnerBody);
AssertPeriod(Outside(400f, onRing ? c / 2 : c), zone.Orbits[rail].Period, "the rail, with a partner at " + partnerOffset);
if (onRing)
AssertPeriod(PeriodUnder(Pull(zone, partnerBody, 800f + partnerOffset), 400f), zone.Orbits[rail].Period,
"the rail, with a partner at " + partnerOffset);
else
Assert.True(zone.Orbits[rail].Period > 1e6f, "a body outside the window is no partner: the rail has no pull of any size");
}
[Fact]
@ -264,7 +223,8 @@ public sealed class OrbitsFromFieldTests : IDisposable
Assert.True(float.IsPositiveInfinity(zone.Orbits[rail].Period), "a body beyond the rail does not pull it");
}
// A bodyless node's children read only that node's ring: the well of the body the node itself orbits does not enter.
// A bodyless node's children read only that node's own terms and their ring partners: the well of the body the node
// itself orbits does not enter, and a lone body on the node has no pull.
[Fact]
public void ABodylessNodesChildrenNeverReadItsParentsWell()
{
@ -286,8 +246,8 @@ public sealed class OrbitsFromFieldTests : IDisposable
var bodyA = w.Body(a, 100f);
var bodyB = w.Body(b, 300f);
zone = w.Build();
AssertPeriod(Outside(200f, C(zone, bodyB) / 2), zone.Orbits[a].Period, "the first member, under its partner alone");
AssertPeriod(Outside(200f, C(zone, bodyA) / 2), zone.Orbits[b].Period, "the second member, under its partner alone");
AssertPeriod(PeriodUnder(Pull(zone, bodyB, 400f), 200f), zone.Orbits[a].Period, "the first member, under its partner alone");
AssertPeriod(PeriodUnder(Pull(zone, bodyA, 400f), 200f), zone.Orbits[b].Period, "the second member, under its partner alone");
}
[Fact]
@ -317,12 +277,12 @@ public sealed class OrbitsFromFieldTests : IDisposable
var binaries = 0;
for (uint seed = 1; seed <= 12; seed++)
{
var dir = Path.Combine(_root, $"galaxy{seed}");
var dir = Path.Combine(_worlds.Root, $"galaxy{seed}");
Directory.CreateDirectory(dir);
var catalog = Path.Combine(dir, "Aetheria.cc");
File.Copy(TestCatalog.Repo, catalog);
var cache = new CultCache(TestCatalog.Registry());
_caches.Add(cache);
_worlds.Caches.Add(cache);
cache.AddBackingStore(new GameCult.Caching.MessagePack.SingleFileMessagePackBackingStore(catalog), AetheriaStores.CatalogTypes);
cache.AddBackingStore(new GameCult.Caching.MessagePack.SingleFileMessagePackBackingStore(Path.Combine(dir, "run.cc")), AetheriaStores.RunTypes);
var authored = AuthoredSettings.Load(RunStartTests.FindRepoRoot());

118
tests/Aetheria.Shared.Tests/PotentialOwnerTests.cs
View file

@ -73,8 +73,8 @@ public sealed class PotentialOwnerTests : IDisposable
private static void Near(float expected, float actual, float relative) =>
Assert.InRange(actual, expected - Math.Abs(expected) * relative - 1e-5f, expected + Math.Abs(expected) * relative + 1e-5f);
// The law written out: potential c r^2 / 2R^2 inside the core, c (ln(r / R) + 1/2) beyond.
private static float Phi(float c, float core, float r) => r < core ? c * r * r / (2 * core * core) : c * (log(r / core) + .5f);
// The law of a lone body written out (its global term): the Plummer potential c ln(r^2 + eps^2) / 2, eps the Roche limit.
private static float Phi(float c, float eps, float r) => c * .5f * log(r * r + eps * eps);
// A zone whose bodies have no authored gravity curves (most test zones) has no pull and no dip, and no NaN either.
[Fact]
@ -103,12 +103,13 @@ public sealed class PotentialOwnerTests : IDisposable
var zone = ZoneOf(planets, (mass, float2.zero));
var core = planets.RocheLimit.Evaluate(mass);
foreach (var r in new[] { 3f * core, 10f * core, PackRadius / 2f })
Near(planets.WellStrength.Evaluate(mass), length(zone.GetForce(new float2(r, 0))) * r, 1e-3f);
Near(planets.WellStrength.Evaluate(mass) * r * r / (r * r + core * core), length(zone.GetForce(new float2(r, 0))) * r, 1e-3f);
}
}
// The steepest pull of the Plummer well is at its softening, which is the Roche curve: it rises up to it and falls beyond.
[Fact]
public void TheCoreRadiusIsTheRocheCurve()
public void TheSteepestPullIsAtTheRocheCurve()
{
var planets = Curves();
foreach (var mass in new[] { 29.2f, 14152.4f })
@ -116,10 +117,8 @@ public sealed class PotentialOwnerTests : IDisposable
var zone = ZoneOf(planets, (mass, float2.zero));
var core = planets.RocheLimit.Evaluate(mass);
float Pull(float r) => length(zone.GetForce(new float2(0, r)));
var slope = Pull(.25f * core) / (.25f * core);
Near(slope, Pull(.5f * core) / (.5f * core), 1e-3f);
Near(slope, Pull(.9f * core) / (.9f * core), 1e-3f);
Assert.True(Pull(core) > Pull(.98f * core), "the pull still rises up to the core radius");
Near(planets.WellStrength.Evaluate(mass) / (2f * core), Pull(core), 1e-3f);
Assert.True(Pull(core) > Pull(.98f * core), "the pull still rises up to the Roche limit");
Assert.True(Pull(core) > Pull(1.02f * core), "and falls beyond it");
}
}
@ -146,10 +145,10 @@ public sealed class PotentialOwnerTests : IDisposable
ZoneOf(Constants(50f), (100f, float2.zero));
}
// core-at-roche-and-exaggeration-up: growing the core to the Roche limit (3 m^.25 to 7.32 m^.25, 2.44x) and dividing the
// global K by that factor leaves the steepest drawn slope (|grad phi| / K, at the core) where it was, and draws every
// slope outside the old core 2.44x steeper. The "before" side is the previous authoring, stated here from the ruling:
// core 3 m^.25 and K 10. The "after" side is the shipped Settings.asset read through the production zone.
// core-at-roche-and-exaggeration-up: the Plummer well's steepest pull is 0.5 c / eps at the Roche limit (eps = 7.32 m^.25),
// half the harmonic core's c / R, so K halves with the law and the steepest drawn slope stays where the previous authoring
// had it. The "before" side is stated from the ruling: core 3 m^.25 and K 10. The "after" side is the shipped
// Settings.asset read through the production zone. Outside the old core the funnel draws steeper than it did.
[Fact]
public void ThePeakDrawnSlopeStaysAndTheFunnelDrawsSteeper()
{
@ -165,10 +164,8 @@ public sealed class PotentialOwnerTests : IDisposable
var peakBefore = c / oldCore / oldExaggeration;
Assert.InRange(Drawn(planets.RocheLimit.Evaluate(mass)), peakBefore * .98f, peakBefore * 1.02f);
var rocheOverOldCore = planets.RocheLimit.Evaluate(mass) / oldCore;
var farBefore = c / (3f * oldCore) / oldExaggeration;
Near(rocheOverOldCore, Drawn(3f * oldCore) / farBefore, 1e-2f);
Assert.InRange(Drawn(3f * oldCore) / farBefore, 2.44f * .99f, 2.44f * 1.01f);
Assert.InRange(Drawn(3f * oldCore) / farBefore, 2.9f, 2.98f);
}
}
@ -247,8 +244,9 @@ public sealed class PotentialOwnerTests : IDisposable
}
}
// A lone body is a global term: it pulls as the Plummer well at every range, with no rim and no gate.
[Fact]
public void TheTailIsAPointMassWithNoRimAndNoGate()
public void TheGlobalTermIsPlummerAtEveryRange()
{
var planets = Curves();
var zone = ZoneOf(planets, (273.9f, float2.zero));
@ -259,34 +257,11 @@ public sealed class PotentialOwnerTests : IDisposable
{
var at = new float2(r * cos(angle), r * sin(angle));
var force = zone.GetForce(at);
Near(c, length(force) * r, 1e-3f);
Near(c * r * r / (r * r + core * core), length(force) * r, 1e-3f);
Assert.True(dot(normalize(force), normalize(-at)) > .9999f, "the pull points at the body");
}
}
[Fact]
public void WellsSuperpose()
{
var planets = Curves();
var a = (mass: 29.2f, at: new float2(-300f, 0f));
var b = (mass: 14152.4f, at: new float2(300f, 0f));
var both = ZoneOf(planets, a, b);
var onlyA = ZoneOf(planets, a);
var onlyB = ZoneOf(planets, b);
var probes = new[]
{
float2(-300f, 0f), float2(300f, 0f), float2(-310f, 4f), float2(305f, -9f), float2(0f, 0f), float2(0f, 250f),
float2(-150f, 100f), float2(150f, -100f), float2(900f, 700f), float2(-1200f, -300f), float2(40f, 3000f), float2(-300f, -40f)
};
foreach (var p in probes)
{
var sum = onlyA.GetForce(p) + onlyB.GetForce(p);
Near(sum.x, both.GetForce(p).x, 1e-4f);
Near(sum.y, both.GetForce(p).y, 1e-4f);
Near(onlyA.GetHeight(p) + onlyB.GetHeight(p), both.GetHeight(p), 1e-4f);
}
}
[Fact]
public void HeightIsThePotentialOverOneK()
{
@ -297,7 +272,7 @@ public sealed class PotentialOwnerTests : IDisposable
var c = planets.WellStrength.Evaluate(273.9f);
var core = planets.RocheLimit.Evaluate(273.9f);
Assert.Equal(0f, lone.GetHeight(new float2(PackRadius, 0)), 4);
Near(-c * (log(PackRadius / core) + .5f) / k, lone.GetHeight(float2.zero), 1e-4f);
Near(-c * .5f * log((PackRadius * PackRadius + core * core) / (core * core)) / k, lone.GetHeight(float2.zero), 1e-4f);
var two = ZoneOf(planets, (29.2f, new float2(-300f, 0f)), (14152.4f, new float2(300f, 0f)));
var c0 = planets.WellStrength.Evaluate(29.2f);
@ -305,65 +280,34 @@ public sealed class PotentialOwnerTests : IDisposable
var c1 = planets.WellStrength.Evaluate(14152.4f);
var core1 = planets.RocheLimit.Evaluate(14152.4f);
const float step = .5f;
// Inside the two cores as well, out to their outer bands (0.5 R to R), where the potential is the quadratic bowl.
// Inside the two softenings as well, where the Plummer potential is smooth through the centre.
var inCore = new[]
{
float2(-300f + .5f * core0, 0f), float2(-300f, .65f * core0), float2(-300f - .86f * core0, 0f), float2(-300f + .6f * core0, .6f * core0),
float2(300f + .55f * core1, 0f), float2(300f, -.7f * core1), float2(300f - .86f * core1, 0f), float2(300f + .9f * core1, 0f)
};
foreach (var p in inCore.Concat(new[] { float2(-200f, 120f), float2(0f, 0f), float2(100f, -300f), float2(500f, 40f), float2(-290f, 12f), float2(310f, 20f), float2(700f, 700f), float2(-80f, -500f) }))
// Two suns 600 apart are two subzones: each is felt only inside its own reach (half the distance, 300), so the height
// there is that well's alone (plus the zone node's constant, which differences cancel) and the force is its alone.
var w0 = float2(-300f, 0f);
var w1 = float2(300f, 0f);
foreach (var p in inCore.Concat(new[] { float2(-200f, 120f), float2(100f, -300f), float2(500f, 40f), float2(-290f, 12f), float2(310f, 20f), float2(700f, 700f), float2(-80f, -500f) }))
{
var expected = Phi(c0, core0, length(p - float2(-300f, 0f))) - Phi(c0, core0, PackRadius)
+ Phi(c1, core1, length(p - float2(300f, 0f))) - Phi(c1, core1, PackRadius);
Near(expected, two.GetHeight(p) * k, 1e-3f);
var near0 = length(p - w0) < 140f; // inside Gamma x reach (150) the well is its Plummer pull whole
var near1 = length(p - w1) < 140f;
if (near0 || near1)
{
var (at, wc, wcore) = near0 ? (w0, c0, core0) : (w1, c1, core1);
var expected = Phi(wc, wcore, length(p - at)) - Phi(wc, wcore, 0f);
Near(expected, (two.GetHeight(p) - two.GetHeight(at)) * k, 1e-3f);
}
// The downhill slope of the potential is the force: the central difference of height * K is -GetForce.
var gradient = float2(
(two.GetHeight(p + float2(step, 0)) - two.GetHeight(p - float2(step, 0))) * k / (2 * step),
(two.GetHeight(p + float2(0, step)) - two.GetHeight(p - float2(0, step))) * k / (2 * step));
var force = two.GetForce(p);
Assert.InRange(length(gradient + force), 0f, length(force) * 1e-2f + 1e-4f);
Assert.InRange(length(gradient + force), 0f, length(force) * 1e-2f + 1e-4f + 3e-5f * k); // flat ground between the subzones reads the float32 resolution of the height, times K
}
}
}
// Inside the core the potential is the bowl c r^2 / 2R^2 all the way out to R: the drop from a lone body's centre to a point
// at r is exactly that, from the middle of the core to a hair inside its edge (the log form agrees with the bowl to second
// order at R, so only a tight tolerance sees a join that sits just inside it).
[Fact]
public void TheCoreIsTheBowlUpToItsEdge()
{
foreach (var k in new[] { 1f, 10f })
{
var planets = Curves(k);
foreach (var mass in new[] { 29.2f, 14152.4f })
{
var zone = ZoneOf(planets, (mass, float2.zero));
var c = planets.WellStrength.Evaluate(mass);
var core = planets.RocheLimit.Evaluate(mass);
foreach (var fraction in new[] { .5f, .75f, .9f, .95f, .98f, .99f, .996f })
foreach (var angle in new[] { 0f, 2.3f })
{
var r = fraction * core;
var drop = (zone.GetHeight(new float2(r * cos(angle), r * sin(angle))) - zone.GetHeight(float2.zero)) * k;
Near(c * fraction * fraction / 2f, drop, 2e-5f);
}
}
}
}
[Fact]
public void TheCoreMeetsTheTailAtItsRadius()
{
var planets = Curves();
foreach (var mass in new[] { 29.2f, 14152.4f })
{
var zone = ZoneOf(planets, (mass, float2.zero));
var core = planets.RocheLimit.Evaluate(mass);
var inside = new float2(core * (1f - 1e-3f), 0);
var outside = new float2(core * (1f + 1e-3f), 0);
Assert.InRange(abs(zone.GetHeight(inside) - zone.GetHeight(outside)), 0f, abs(zone.GetHeight(outside)) * 1e-3f);
Assert.InRange(abs(length(zone.GetForce(inside)) - length(zone.GetForce(outside))), 0f, length(zone.GetForce(outside)) * 1e-3f);
}
}
}

4
tests/Aetheria.Shared.Tests/WellFixture.cs
View file

@ -2,7 +2,7 @@ using CultMath;
using GameCult.Caching;
// The one fixture that gives a test zone a slope: a body on a parentless orbit (the origin unless `at` says otherwise; call again for more bodies), with the
// curves a body needs. Gravity is the sum of body wells (GravityField), so a test that wants ground to slide on
// curves a body needs. Gravity is the zone's WellTree over its bodies' masses, so a body carries a mass (any: the constant curves ignore it), and a test that wants ground to slide on
// digs it here and nowhere else. `strength` is the well's c (u^2/s^2) and `coreRadius` its softened core in u, which is
// the Roche limit curve (the one curve of a body's core).
internal static class WellFixture
@ -13,7 +13,7 @@ internal static class WellFixture
settings.WellStrength = new ExponentialCurve { Constant = strength };
settings.RocheLimit = new ExponentialCurve { Constant = coreRadius };
var orbit = cache.Upsert(new OrbitData { FixedPosition = at }).Key;
var body = cache.Upsert(new PlanetData { Name = $"Well {at.x},{at.y}", Orbit = new CultRecordRef<OrbitData>(orbit) }).Key;
var body = cache.Upsert(new PlanetData { Name = $"Well {at.x},{at.y}", Orbit = new CultRecordRef<OrbitData>(orbit), Mass = 1f }).Key;
pack.Orbits.Add(new CultRecordRef<OrbitData>(orbit));
pack.Planets.Add(new CultRecordRef<BodyData>(body));
}

533
tests/Aetheria.Shared.Tests/WellLawTests.cs Normal file
View file

@ -0,0 +1,533 @@
using System;
using System.Linq;
using CultMath;
using GameCult.Caching;
using Xunit;
using static CultMath.math;
using float2 = CultMath.float2;
// aetheria-gravity, cut well-law: the conserving compact law. GravityField is the law of one band term; WellTree turns a zone's
// orbit tree into terms and walks them. The oracles here are written out from the law (the order-4 smoothstep, Simpson
// integration of the pull in double, merging a family into one body), never by calling the code under test twice.
public sealed class WellLawTests : IDisposable
{
private readonly ZoneWorlds _worlds = new ZoneWorlds();
public void Dispose() => _worlds.Dispose();
// The Plummer pull of a term without a band: r / (r^2 + eps^2).
private static double Plummer(double eps, double r) => r / (r * r + eps * eps);
// The order-4 smoothstep, spelled: S(t) = 35 t^4 - 84 t^5 + 70 t^6 - 20 t^7.
private static double S(double t) => t * t * t * t * (35 - 84 * t + 70 * t * t - 20 * t * t * t);
private static double Simpson(Func<double, double> f, double from, double to, int n = 4000)
{
var h = (to - from) / n;
var sum = f(from) + f(to);
for (var i = 1; i < n; i++) sum += f(from + i * h) * (i % 2 == 1 ? 4 : 2);
return sum * h / 3;
}
public static readonly object[][] Terms =
{
new object[] { 79.8f, 1000f, 1500f },
new object[] { 17f, 40f, 60f },
new object[] { 79.8f, 0f, 1525f },
new object[] { 1658f, 1000f, 1500f }, // tau 3: the series branch
new object[] { 59f, 11f, 22f }, // tau 10
new object[] { 79.8f, float.PositiveInfinity, float.PositiveInfinity }
};
// The pull is the slope of the potential: the potential's drop between two radii is the integral of the pull between
// them, by Simpson in double, for terms on both branches, inside the band, across its edges and beyond it.
[Theory]
[MemberData(nameof(Terms))]
public void TheTermIsItsPotentialsSlope(float eps, float inner, float outer)
{
var span = float.IsPositiveInfinity(outer) ? 3000f : outer;
foreach (var (from, to) in new[] { (.2f, .55f), (.4f, .99f), (.9f, 1.3f), (.05f, 1.0f) })
{
float r0 = from * span, r1 = Math.Min(to * span, float.IsPositiveInfinity(outer) ? 6000f : outer * 1.1f);
var drop = GravityField.Potential(eps, inner, outer, r1) - GravityField.Potential(eps, inner, outer, r0);
var integral = Simpson(r => GravityField.Pull(eps, inner, outer, (float) r), r0, r1);
Assert.InRange(drop, integral - 3e-5 - 1e-4 * Math.Abs(integral), integral + 3e-5 + 1e-4 * Math.Abs(integral));
}
if (!float.IsPositiveInfinity(outer))
{
Assert.Equal(0f, GravityField.Potential(eps, inner, outer, outer));
Assert.Equal(0f, GravityField.Potential(eps, inner, outer, outer * 1.5f));
Assert.Equal(0f, GravityField.Pull(eps, inner, outer, outer * 1.5f));
}
}
// The pull is the Plummer pull faded by 1 - S(t) across the band, S the order-4 smoothstep: a smoothstep of order 2 has a
// different profile, and the profile's first three derivatives vanish at both edges (the deficit goes as t^4).
[Fact]
public void ABandEdgeIsSmooth()
{
const float eps = 79.8f, a = 1000f, b = 1500f;
foreach (var t in new[] { 0f, .02f, .05f, .1f, .3f, .5f, .7f, .9f, .98f, 1f })
{
var r = (float) Math.Sqrt(a * a + t * (b * b - a * a));
var pull = GravityField.Pull(eps, a, b, r);
var expected = Plummer(eps, r) * (1 - S(t));
Assert.InRange(pull, expected - 1e-5 * Plummer(eps, r), expected + 1e-5 * Plummer(eps, r));
}
Assert.Equal((double) GravityField.Pull(eps, a, b, a), Plummer(eps, a), 6);
Assert.Equal(0f, GravityField.Pull(eps, a, b, b));
var deficit = new Func<float, double>(t => 1 - GravityField.Pull(eps, a, b, (float) Math.Sqrt(a * a + t * (b * b - a * a))) / Plummer(eps, Math.Sqrt(a * a + t * (b * b - a * a))));
Assert.InRange(deficit(.1f) / deficit(.05f), 13.0, 15.0); // doubling t multiplies the deficit by about 14 (t^4), not by 4
}
// The band's potential is the integral of its pull, in single precision, on both sides of the series threshold: tau =
// (a^2 + eps^2) / (b^2 - a^2) is swept from 0.05 to 10 (every (a, b) where that tau has a real softening), and the
// potential stays within 1.8e-6 per unit strength of the double-precision Simpson integral of the pull. A closed form kept
// past tau 0.6 or a series cut short misses this by 1e-5 to 1e-3. Both branches are within the bound between tau 0.3 and 0.6,
// so the exact seam (SeriesFrom 0.4) inside that interval is pinned by the committed shader being the generated one, not here.
[Fact]
public void TheBandPotentialIsExactInSinglePrecision()
{
var checkedCases = 0;
foreach (var tau in new[] { .05, .3, .35, .39, .41, .45, .5, .55, .6, .7, 1, 1.5, 2.5, 10 })
foreach (var (a, b) in new[] { (1000.0, 1500.0), (40.0, 60.0), (11.0, 22.0), (0.0, 1525.0), (150.0, 4000.0) })
{
var e2 = tau * (b * b - a * a) - a * a;
if (e2 < 0) continue;
var eps = (float) Math.Sqrt(e2);
e2 = (double) eps * eps; // the oracle integrates the softening the code is given
for (var k = 0; k < 15; k++)
{
var r = b * (.02 + .979 * k / 14);
var oracle = -Simpson(x =>
{
var t = x <= a ? 0.0 : (x * x - a * a) / (b * b - a * a);
return x / (x * x + e2) * (1 - S(t));
}, r, b, 8000);
var actual = GravityField.Potential(eps, (float) a, (float) b, (float) r);
Assert.InRange(actual, oracle - 1.8e-6, oracle + 1.8e-6);
checkedCases++;
}
}
Assert.True(checkedCases > 400, "the sweep ran its cases");
}
private static float Plummer(PlanetSettings curves, float mass, float d) =>
curves.WellStrength.Evaluate(mass) * d / (d * d + curves.RocheLimit.Evaluate(mass) * curves.RocheLimit.Evaluate(mass));
// A sun, a planet at 3000 and a moon at 150 round the planet, the moon's mass as given.
private Zone SunPlanetMoon(float moonMass, out CultRecordKey planetOrbit, out CultRecordKey moonOrbit, float planetMass = 100f)
{
var w = _worlds.New();
var sunOrbit = w.Orbit(default, 0f);
w.Body(sunOrbit, 5000f);
planetOrbit = w.Orbit(sunOrbit, 3000f);
w.Body(planetOrbit, planetMass);
moonOrbit = w.Orbit(planetOrbit, 150f);
if (moonMass > 0f) w.Body(moonOrbit, moonMass);
return w.Build();
}
// Inside its outermost child a node pulls exactly as its body alone: the family terms cancel there, so a child's rail
// never feels its own mass through its parent.
[Fact]
public void InsideAFamilyItsParentPullsAsItsBodyAlone()
{
var curves = ZoneWorlds.Curves();
foreach (var moonMass in new[] { 1f, 2f })
{
var zone = SunPlanetMoon(moonMass, out var planetOrbit, out var moonOrbit);
foreach (var d in new[] { 20f, 60f, 150f })
{
var expected = Plummer(curves, 100f, d);
Assert.InRange(zone.Wells.PullAround(planetOrbit, d, moonOrbit), expected * (1 - 1e-5f), expected * (1 + 1e-5f));
}
}
}
// Merging a family into one body of its mass in the same slot changes nothing beyond the family's support, and does
// inside it.
[Fact]
public void AFamilyMergedIntoOneBodyIsTheSameFieldOutsideIt()
{
var family = SunPlanetMoon(1f, out var planetOrbit, out _);
var merged = SunPlanetMoon(0f, out _, out _, planetMass: 101f);
var support = family.Wells.Support(planetOrbit);
Assert.True(support > 150f, "fixture: the support clears the moon");
var planetAt = family.Orbits[planetOrbit].Position;
foreach (var (distance, angle) in new[] { (1.05f, .3f), (1.4f, 2.2f), (3f, 4.0f), (9f, 1.0f), (20f, 5.5f) })
{
var p = planetAt + float2(support * distance * cos(angle), support * distance * sin(angle));
var f = family.GetForce(p);
var m = merged.GetForce(p);
Assert.InRange(length(f - m), 0f, 1e-5f * length(f) + 1e-9f);
Assert.InRange(abs(family.GetHeight(p) - merged.GetHeight(p)), 0f, 1e-5f * (1f + abs(family.GetHeight(p))));
}
var inside = planetAt + float2(100f, 0f);
Assert.True(length(family.GetForce(inside) - merged.GetForce(inside)) > 1e-3f * length(merged.GetForce(inside)), "inside the support the moon is a body of its own");
}
// A sole root has no reach: far from everything the zone is one body of its total mass, however the generator split it.
[Fact]
public void TheZonesFarFieldIsItsTotalMass()
{
var curves = ZoneWorlds.Curves();
var one = SunPlanetMoon(1f, out _, out _);
var w = _worlds.New();
var sunOrbit = w.Orbit(default, 0f);
w.Body(sunOrbit, 5000f);
for (var i = 0; i < 5; i++) w.Body(w.Orbit(sunOrbit, 800f + 450f * i, i / 5f), 20.2f);
var five = w.Build();
foreach (var zone in new[] { one, five })
{
var far = 10f * zone.Pack.Radius;
foreach (var angle in new[] { 0f, 1.7f, 4.1f })
Near(curves.WellStrength.Evaluate(5101f), length(zone.GetForce(float2(far * cos(angle), far * sin(angle)))) * far, 1e-3f);
}
}
// Three subzones: a sun of 5000 with a planet, a sun of 900 alone, and a bodyless pair; each is a root (a parentless orbit
// stands where its distance and phase put it). Masses are the test's own, the zone's total 6100.
private Zone ThreeSubzones(out CultRecordKey[] roots)
{
var w = _worlds.New();
var big = w.Orbit(default, 3000f, .85f);
w.Body(big, 5000f);
w.Body(w.Orbit(big, 800f), 100f);
var small = w.Orbit(default, 9000f, .3f);
w.Body(small, 900f);
var pair = w.Orbit(default, 12000f, .7f);
w.Body(w.Orbit(pair, 300f, 0f), 60f);
w.Body(w.Orbit(pair, 300f, .5f), 40f);
roots = new[] { big, small, pair };
return w.Build();
}
private static readonly float[] ThreeSubzoneMass = { 5100f, 900f, 100f };
// The zone node sits at the family-mass-weighted mean of its roots' places and pulls the zone's total mass from there.
[Fact]
public void AZoneOfSeveralSunsPullsAsItsTotalMassFromTheirBarycentre()
{
var zone = ThreeSubzones(out var roots);
var centre = float2.zero;
for (var i = 0; i < roots.Length; i++) centre += zone.Orbits[roots[i]].Position * ThreeSubzoneMass[i] / 6100f;
Assert.InRange(length(zone.Wells.Centre - centre), 0f, 1e-3f);
// Each subzone reaches half the way to its nearest neighbour, or 1.5 of its support if more.
for (var i = 0; i < roots.Length; i++)
{
var nearest = Enumerable.Range(0, roots.Length).Where(j => j != i).Min(j => length(zone.Orbits[roots[j]].Position - zone.Orbits[roots[i]].Position));
Assert.Equal(Math.Max(.5f * nearest, 1.5f * zone.Wells.Support(roots[i])), zone.Wells.Reach(roots[i]), 2);
}
var far = 10f * zone.Pack.Radius;
foreach (var angle in new[] { 0f, 2.1f, 4.4f })
Near(ZoneWorlds.Curves().WellStrength.Evaluate(6100f),
length(zone.GetForce(zone.Wells.Centre + float2(far * cos(angle), far * sin(angle)))) * far, 1e-3f);
}
// One sun with a planet, alone in its zone at the place a root of `distance` and `phase` stands.
private Zone OneSubzone(float sun, float distance, float phase, out CultRecordKey root)
{
var w = _worlds.New();
root = w.Orbit(default, distance, phase);
w.Body(root, sun);
w.Body(w.Orbit(root, 700f, .2f), 120f);
return w.Build();
}
// Isolated subzones: inside its own support a sun is felt exactly as if it were alone in the zone, the other suns and the
// zone node adding nothing; a root's reach is half the distance to the nearest other sun, or 1.5 of its support if more.
[Fact]
public void ASubzoneFeelsOnlyItsOwnSun()
{
var w = _worlds.New();
var a = w.Orbit(default, 0f);
w.Body(a, 4000f);
w.Body(w.Orbit(a, 700f, .2f), 120f);
var b = w.Orbit(default, 9000f, 0f);
w.Body(b, 2500f);
w.Body(w.Orbit(b, 700f, .2f), 120f);
var zone = w.Build();
var rng = new System.Random(3);
foreach (var (root, sun, distance) in new[] { (a, 4000f, 0f), (b, 2500f, 9000f) })
{
var alone = OneSubzone(sun, distance, 0f, out var aloneRoot);
Assert.True(float.IsPositiveInfinity(alone.Wells.Reach(aloneRoot)), "a sole massive root is the whole zone");
var support = zone.Wells.Support(root);
Assert.Equal(Math.Max(.5f * 9000f, 1.5f * support), zone.Wells.Reach(root), 2);
var at = zone.Orbits[root].Position;
for (var i = 0; i < 50; i++)
{
var angle = (float) (rng.NextDouble() * 2 * Math.PI);
var radius = support * (.2f + .8f * (float) rng.NextDouble());
var p = at + float2(radius * cos(angle), radius * sin(angle));
var expected = alone.GetForce(p);
Assert.InRange(length(zone.GetForce(p) - expected), 0f, 1e-5f * length(expected) + 1e-9f);
}
}
}
// Between the subzones the ground is flat; beyond the populated area the zone is one body of its total mass at its centre.
[Fact]
public void BeyondEverySubzoneTheZoneIsOneBody()
{
var zone = ThreeSubzones(out var roots);
var curves = ZoneWorlds.Curves();
var centre = zone.Wells.Centre;
var span = roots.Max(r => length(zone.Orbits[r].Position - centre) + zone.Wells.Reach(r)); // A, as the law defines it
var flat = 0;
for (var x = -14000f; x <= 14000f; x += 700f)
for (var y = -14000f; y <= 14000f; y += 700f)
{
var p = float2(x, y);
if (length(p - centre) >= span || roots.Any(r => length(p - zone.Orbits[r].Position) <= zone.Wells.Reach(r))) continue;
flat++;
Assert.InRange(length(zone.GetForce(p)), 0f, 1e-7f);
}
Assert.True(flat > 100, "fixture: there is open ground between the subzones");
// The rim: from A to 1.5 A the zone's pull rises from nothing to its whole Plummer pull as the order-4 smoothstep of
// (d^2 - A^2) / ((1.5 A)^2 - A^2), A being the reach-derived radius above.
foreach (var (radius, angle) in new[] { (1.05f, .7f), (1.15f, 3.1f), (1.25f, 5.0f), (1.35f, 1.9f), (1.45f, 4.2f) })
{
var p = centre + float2(span * radius * cos(angle), span * radius * sin(angle));
var d = length(centre - p);
var t = ((double) d * d - (double) span * span) / (2.25 * span * span - (double) span * span);
var expected = (centre - p) / d * (float) (Plummer(curves, 6100f, d) * S(t));
Assert.InRange(length(zone.GetForce(p) - expected), 0f, 2e-5f * Plummer(curves, 6100f, d));
}
foreach (var (radius, angle) in new[] { (1.6f, .4f), (3f, 2.5f), (8f, 5.1f) })
{
var p = centre + float2(span * radius * cos(angle), span * radius * sin(angle));
var expected = (centre - p) / length(centre - p) * Plummer(curves, 6100f, length(centre - p));
Assert.InRange(length(zone.GetForce(p) - expected), 0f, 1e-5f * length(expected));
}
}
// A belt is a child of its parent's family: the family's support covers the belt's outermost asteroid, so every asteroid
// circles the planet's body alone at its own distance (not at the belt orbit's, 43, where the fade used to start).
[Fact]
public void ABeltRidesInsideItsParentsBody()
{
var w = _worlds.New();
var sun = w.Orbit(default, 0f);
w.Body(sun, 5000f);
var planet = w.Orbit(sun, 3000f);
w.Body(planet, 172f);
var belt = w.Belt(w.Orbit(planet, 43f, .4f), 10f, 87f, 43f); // the outermost asteroid is in the middle: neither the first nor the last listed
var zone = w.Build();
Assert.Equal(1.5f * 87f, zone.Wells.Support(planet), 3);
var curves = ZoneWorlds.Curves();
var periods = zone.AsteroidBelts[belt].Periods;
Assert.Equal(3, periods.Length);
for (var i = 0; i < 3; i++)
{
var d = zone.AsteroidBelts[belt].Data.Asteroids[i].Distance;
var expected = 2 * PI * sqrt(d / Plummer(curves, 172f, d));
Assert.InRange(periods[i], expected * (1 - 1e-5f), expected * (1 + 1e-5f));
}
}
// The fade of a tight family begins at its support: the planet (moon at 150, support 225, reach at its floor 337.5) pulls
// as the whole family out to 225 and then gives way as the order-4 smoothstep of (r^2 - 225^2) / (337.5^2 - 225^2).
[Fact]
public void ATightFamilysPullFadesFromItsSupport()
{
var w = _worlds.New();
var sun = w.Orbit(default, 0f);
w.Body(sun, 5000f);
var planet = w.Orbit(sun, 400f);
w.Body(planet, 100f);
var moon = w.Orbit(planet, 150f);
w.Body(moon, 1f);
var zone = w.Build();
var curves = ZoneWorlds.Curves();
Assert.Equal(225f, zone.Wells.Support(planet), 3);
Assert.Equal(337.5f, zone.Wells.Reach(planet), 3);
foreach (var r in new[] { 232f, 250f, 270f, 290f, 310f, 325f, 333f, 336f })
{
var t = (r * r - 225.0 * 225) / (337.5 * 337.5 - 225.0 * 225);
var expected = Plummer(curves, 101f, r) * (1 - S(t));
var scale = Plummer(curves, 101f, r); // float rounding of the four cancelling terms is relative to the Plummer pull
Assert.InRange(zone.Wells.PullAround(planet, r, moon), expected - 1e-5 * scale, expected + 1e-5 * scale);
}
}
// A wide family fades from Gamma x its reach, not its support: a planet at 2000 with a moon at 100 has support 150 and reach
// 1000 (half its distance), so it pulls as the whole family out to 500 and gives way from there as the order-4 smoothstep
// of (r^2 - 500^2) / (1000^2 - 500^2).
[Fact]
public void AWideFamilysPullFadesFromGammaTimesItsReach()
{
var w = _worlds.New();
var sun = w.Orbit(default, 0f);
w.Body(sun, 5000f);
var planet = w.Orbit(sun, 2000f);
w.Body(planet, 100f);
var moon = w.Orbit(planet, 100f);
w.Body(moon, 1f);
var zone = w.Build();
var curves = ZoneWorlds.Curves();
Assert.Equal(150f, zone.Wells.Support(planet), 3);
Assert.Equal(1000f, zone.Wells.Reach(planet), 3);
foreach (var r in new[] { 480f, 520f, 600f, 700f, 800f, 900f, 980f })
{
var t = Math.Max(0.0, (r * r - 500.0 * 500) / (1000.0 * 1000 - 500.0 * 500));
var expected = Plummer(curves, 101f, r) * (1 - S(t));
var scale = Plummer(curves, 101f, r);
Assert.InRange(zone.Wells.PullAround(planet, r, moon), expected - 1e-5 * scale, expected + 1e-5 * scale);
}
}
// A ring group's pull rises from its ring: a pair at 300 round a bodyless node pulls from the node as the order-4 smoothstep
// of (r^2 - 300^2) / (1500^2 - 300^2) of the pair's Plummer, 1500 being the ring's distance and the members' RingReach.
[Fact]
public void ARingGroupsPullRisesFromItsRing()
{
var w = _worlds.New();
var sun = w.Orbit(default, 0f);
w.Body(sun, 5000f);
var node = w.Orbit(sun, 3000f);
w.Body(w.Orbit(node, 300f, 0f), 80f);
w.Body(w.Orbit(node, 300f, .5f), 60f);
var zone = w.Build();
var curves = ZoneWorlds.Curves();
foreach (var r in new[] { 310f, 400f, 600f, 900f, 1200f, 1400f, 1490f })
{
var t = (r * r - 300.0 * 300) / (1500.0 * 1500 - 300.0 * 300);
var expected = Plummer(curves, 140f, r) * S(t);
var scale = Plummer(curves, 140f, r); // float rounding of the four cancelling terms is relative to the Plummer pull
Assert.InRange(zone.Wells.PullAround(node, r, default), expected - 1e-5 * scale, expected + 1e-5 * scale);
}
}
// The drawn height is the force's potential everywhere, however far: its slope along the radius is the pull, for a zone of
// one sun and for one of three subzones, out to ten pack radii (a potential flattened at some radius reads zero here).
[Fact]
public void FarOutTheHeightIsTheForcesPotential()
{
var single = SunPlanetMoon(1f, out _, out _);
var several = ThreeSubzones(out _);
foreach (var (zone, centre) in new[] { (single, float2.zero), (several, several.Wells.Centre) })
foreach (var radius in new[] { 2f, 4f, 7f, 10f })
foreach (var angle in new[] { .3f, 2.4f, 4.6f })
{
var outward = float2(cos(angle), sin(angle));
var p = centre + outward * (radius * zone.Pack.Radius);
const float h = 10f;
var slope = (zone.GetHeight(p + outward * h) - zone.GetHeight(p - outward * h)) * zone.Settings.HeightExaggeration / (2 * h);
var pull = -dot(zone.GetForce(p), outward);
Assert.InRange(slope, pull * .98f, pull * 1.02f);
}
}
private static void Near(float expected, float actual, float relative) =>
Assert.InRange(actual, expected * (1 - relative), expected * (1 + relative));
// A node's reach is half its distance, or one and a half times its support if that is more; the root's is infinite.
[Fact]
public void AReachIsHalfTheDistanceOrOneAndAHalfTheSupport()
{
var w = _worlds.New();
var sunOrbit = w.Orbit(default, 0f);
w.Body(sunOrbit, 5000f);
var wide = w.Orbit(sunOrbit, 3000f);
w.Body(wide, 100f);
var tight = w.Orbit(sunOrbit, 100f);
w.Body(tight, 100f);
w.Body(w.Orbit(wide, 150f), 1f);
w.Body(w.Orbit(tight, 150f), 1f);
var zone = w.Build();
Assert.True(float.IsPositiveInfinity(zone.Wells.Reach(sunOrbit)));
Assert.Equal(1500f, zone.Wells.Reach(wide), 3);
Assert.Equal(225f, zone.Wells.Support(wide), 3);
Assert.Equal(337.5f, zone.Wells.Reach(tight), 3);
Assert.Equal(225f, zone.Wells.Support(tight), 3);
}
// A lone body on a bodyless node is pulled by nothing of its own: the node's terms cancel where the ring is.
[Fact]
public void ALoneBodyOnABodylessNodeFeelsNoPullOfItsOwn()
{
var w = _worlds.New();
var sunOrbit = w.Orbit(default, 0f);
w.Body(sunOrbit, 5000f);
var node = w.Orbit(sunOrbit, 3000f);
var alone = w.Orbit(node, 200f);
w.Body(alone, 100f);
var zone = w.Build();
Assert.Equal(0f, zone.Wells.PullAround(node, 200f, alone));
Assert.True(float.IsPositiveInfinity(zone.Orbits[alone].Period));
}
// The walk skips a subtree outside its reach, which is exact: at points all over a tree with families, a ring and a belt of
// moons, the force and potential of the walk equal the flat sum over every term to float rounding of the terms.
[Fact]
public void TheWalkIsTheSum()
{
var w = _worlds.New();
var sunOrbit = w.Orbit(default, 0f);
w.Body(sunOrbit, 5000f);
for (var i = 0; i < 5; i++)
{
var planet = w.Orbit(sunOrbit, 700f + 800f * i, frac(i * .3819f));
w.Body(planet, 40f + 30f * i);
for (var j = 0; j < 3; j++) w.Body(w.Orbit(planet, 90f + 70f * j, frac(j * .6180339887f)), 1f + j);
}
var node = w.Orbit(sunOrbit, 5200f, .13f);
w.Body(w.Orbit(node, 300f, 0f), 80f);
w.Body(w.Orbit(node, 300f, .5f), 60f);
w.Body(w.Orbit(default, 9000f, .3f), 900f); // a second massive root: the zone node joins the walk
var zone = w.Build();
var rng = new System.Random(7);
for (var n = 0; n < 200; n++)
{
var p = float2((float) (rng.NextDouble() * 24000 - 6000), (float) (rng.NextDouble() * 24000 - 12000));
var force = float2.zero;
var potential = 0f;
var scale = 0f;
foreach (var term in zone.Wells.Terms)
{
var offset = (term.Node.IsSet() ? zone.Orbits[term.Node].Position : zone.Wells.Centre) - p; // an unset node is the zone node
var r = length(offset);
var pull = term.Strength * GravityField.Pull(term.Softening, term.Inner, term.Outer, r);
if (r > 0f) force += offset / r * pull;
var level = term.Strength * GravityField.Potential(term.Softening, term.Inner, term.Outer, r) - term.Reference;
potential += level;
scale += abs(pull) + abs(level);
}
Assert.InRange(length(zone.GetForce(p) - force), 0f, 2e-6f * scale + 1e-9f);
Assert.InRange(abs(zone.GetHeight(p) * zone.Settings.HeightExaggeration - potential), 0f, 2e-6f * scale + 1e-9f);
}
}
// A ring of three on a bodyless node: each member feels both partners' own pulls at their chord (the side of the inscribed
// triangle), each projected toward the node (cos 30 degrees), so the projection and the chord both decide the period.
[Fact]
public void ARosetteMemberFeelsItsPartnersProjectedTowardTheNode()
{
var w = _worlds.New();
var node = w.Orbit(default, 0f);
var members = new[] { w.Orbit(node, 400f, 0f), w.Orbit(node, 400f, 1f / 3f), w.Orbit(node, 400f, 2f / 3f) };
var bodies = members.Select(m => w.Body(m, 120f)).ToArray();
var zone = w.Build();
var chord = 400f * sqrt(3f);
var pull = 2f * zone.PlanetInstances[bodies[0]].WellStrength * GravityField.Pull(zone.PlanetInstances[bodies[0]].RocheLimit,
float.PositiveInfinity, float.PositiveInfinity, chord) * cos(PI / 6f);
foreach (var member in members)
Assert.InRange(zone.Orbits[member].Period, 2 * PI * sqrt(400f / pull) * (1 - 1e-4f), 2 * PI * sqrt(400f / pull) * (1 + 1e-4f));
}
}

85
tests/Aetheria.Shared.Tests/ZoneWorlds.cs Normal file
View file

@ -0,0 +1,85 @@
using System;
using System.Collections.Generic;
using System.IO;
using System.Linq;
using CultMath;
using GameCult.Caching;
using Xunit;
using static CultMath.math;
using float2 = CultMath.float2;
using float3 = CultMath.float3;
// The hand-built zones of the tree tests: packs of orbits and bodies read by the production Zone constructor with the shipped
// well curves (c = 45.08 m^.2, eps = 7.32 m^.25), each in its own store under one temporary root.
internal sealed class ZoneWorlds : IDisposable
{
public readonly string Root = Path.Combine(Path.GetTempPath(), "aetheria-orbitsfield-" + Guid.NewGuid().ToString("N"));
public readonly List<CultCache> Caches = new List<CultCache>();
public ZoneWorlds() => Directory.CreateDirectory(Root);
public void Dispose()
{
foreach (var cache in Caches) cache.Dispose();
Directory.Delete(Root, true);
}
private static GameplaySettings GameSettings() => new GameplaySettings
{
DefaultEntitySettings = new EntitySettings(),
Tiers = new[] { new RarityTier { Name = "Common", Quality = .5f, Rarity = 0, Color = new float3(1, 1, 1) } },
QualityPriceModifier = new ExponentialLerp()
};
public static PlanetSettings Curves() => new PlanetSettings
{
WellStrength = new ExponentialCurve { Exponent = .2f, Multiplier = 45.08f },
RocheLimit = new ExponentialCurve { Exponent = .25f, Multiplier = 7.32f },
BodyRadius = new ExponentialCurve { Exponent = .25f, Multiplier = 3f },
LightRadius = new ExponentialCurve(),
AsteroidSize = new ExponentialLerp { Minimum = 2f, Maximum = 9f, Exponent = 1.7f }
};
public sealed class World
{
public CultCache Cache;
public ItemManager Items;
public ZonePack Pack = new ZonePack { Radius = 20000f, Mass = 10000f };
public CultRecordKey Orbit(CultRecordKey parent, float distance, float phase = 0f)
{
var data = new OrbitData { Distance = distance, Phase = phase };
if (parent.IsSet()) data.Parent = new CultRecordRef<OrbitData>(parent);
var key = Cache.Upsert(data).Key;
Pack.Orbits.Add(new CultRecordRef<OrbitData>(key));
return key;
}
public CultRecordKey Body(CultRecordKey orbit, float mass)
{
var key = Cache.Upsert(new PlanetData { Name = "Body", Orbit = new CultRecordRef<OrbitData>(orbit), Mass = mass }).Key;
Pack.Planets.Add(new CultRecordRef<BodyData>(key));
return key;
}
// A belt of asteroids at the given distances riding `orbit`; returns the belt's body key.
public CultRecordKey Belt(CultRecordKey orbit, params float[] distances)
{
var asteroids = distances.Select((d, i) => new Asteroid { Distance = d, Phase = frac(i * .6180339887f), Size = .5f, RotationSpeed = .1f }).ToArray();
var key = Cache.Upsert(new AsteroidBeltData { Orbit = new CultRecordRef<OrbitData>(orbit), Asteroids = asteroids }).Key;
Pack.Planets.Add(new CultRecordRef<BodyData>(key));
return key;
}
public Zone Build() => new Zone(Items, Curves(), Pack, new GalaxyZone { Name = "OrbitsFromField", Owner = null }, null);
}
public World New()
{
var dir = Path.Combine(Root, $"world{Caches.Count}");
Directory.CreateDirectory(dir);
var cache = AetheriaStores.Open(Path.Combine(dir, "Aetheria.cc"), Path.Combine(dir, "run.cc"), catalogWritable: true);
Caches.Add(cache);
return new World { Cache = cache, Items = new ItemManager(cache, new ProvenanceLedger(), GameSettings(), _ => { }) };
}
}