diff --git a/Assets/Shaders/Gravity.cginc b/Assets/Shaders/Gravity.cginc index 47a419518..b6b4d3300 100644 --- a/Assets/Shaders/Gravity.cginc +++ b/Assets/Shaders/Gravity.cginc @@ -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 _GravityWells; +StructuredBuffer _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 <= 1, 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 > 1, 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 <= 1 ? gravityBandClosed(tau, 0) : gravityBandSeries(tau, 0)); + float t = (s - a2) / (b2 - a2); + return tau <= 1 ? 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; } diff --git a/tests/Aetheria.Shared.Tests/DrawnPoseTests.cs b/tests/Aetheria.Shared.Tests/DrawnPoseTests.cs index 355c57432..acc7f8f72 100644 --- a/tests/Aetheria.Shared.Tests/DrawnPoseTests.cs +++ b/tests/Aetheria.Shared.Tests/DrawnPoseTests.cs @@ -125,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 * sqrt(distance / (strength * distance / (distance * distance + 4f))); + var period = PeriodAt(distance, strength); var wrapped = (float) frac(time / period); return center + OrbitData.Evaluate(frac(wrapped + phase)) * distance; } @@ -334,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]); diff --git a/tests/Aetheria.Shared.Tests/MiningCut3Tests.cs b/tests/Aetheria.Shared.Tests/MiningCut3Tests.cs index fdabaf6bd..fd522e357 100644 --- a/tests/Aetheria.Shared.Tests/MiningCut3Tests.cs +++ b/tests/Aetheria.Shared.Tests/MiningCut3Tests.cs @@ -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(sunOrbit) }).Key; + var sun = cache.Upsert(new SunData { Orbit = new CultRecordRef(sunOrbit), Mass = 1f }).Key; var pack = new ZonePack { Orbits = { new CultRecordRef(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(orbit)); - pack.Planets.Add(new CultRecordRef(cache.Upsert(new SunData { Orbit = new CultRecordRef(orbit) }).Key)); + pack.Planets.Add(new CultRecordRef(cache.Upsert(new SunData { Orbit = new CultRecordRef(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]); diff --git a/tests/Aetheria.Shared.Tests/OrbitsFromFieldTests.cs b/tests/Aetheria.Shared.Tests/OrbitsFromFieldTests.cs index 10de2c62e..d20fe1436 100644 --- a/tests/Aetheria.Shared.Tests/OrbitsFromFieldTests.cs +++ b/tests/Aetheria.Shared.Tests/OrbitsFromFieldTests.cs @@ -9,11 +9,11 @@ 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")); @@ -78,7 +78,15 @@ public sealed class OrbitsFromFieldTests : IDisposable 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 +103,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 +120,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 +136,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 +153,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 +166,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 +199,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 +221,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 +234,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 +255,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 +276,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 +299,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] diff --git a/tests/Aetheria.Shared.Tests/WellFixture.cs b/tests/Aetheria.Shared.Tests/WellFixture.cs index 08bed1ad4..60e2ea29e 100644 --- a/tests/Aetheria.Shared.Tests/WellFixture.cs +++ b/tests/Aetheria.Shared.Tests/WellFixture.cs @@ -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(orbit) }).Key; + var body = cache.Upsert(new PlanetData { Name = $"Well {at.x},{at.y}", Orbit = new CultRecordRef(orbit), Mass = 1f }).Key; pack.Orbits.Add(new CultRecordRef(orbit)); pack.Planets.Add(new CultRecordRef(body)); }