3.0.0 SIMD

.github/workflows/ci.yml

@@ -0,0 +1,34 @@
+name: NuGet Pack Check
+
+on:
+  push:
+    branches: [main]
+  pull_request:
+    branches: [main]
+  workflow_dispatch:
+
+jobs:
+  pack:
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v7
+
+      - uses: actions/setup-dotnet@v6
+        with:
+          dotnet-version: 10.0.x
+
+      - name: Pack
+        run: dotnet pack VSOP2013.NET/VSOP2013.NET.csproj -c Release -o artifacts
+
+      - name: Package smoke test
+        run: |
+          VERSION=$(dotnet msbuild VSOP2013.NET/VSOP2013.NET.csproj -getProperty:Version)
+          dotnet run --project PackageSmokeTest -c Release -p:VsopVersion=$VERSION
+
+      - uses: actions/upload-artifact@v7
+        with:
+          name: nuget-package
+          retention-days: 7
+          path: |
+            artifacts/*.nupkg
+            artifacts/*.snupkg

DataConverter/DataConverter.csproj

@@ -52,7 +52,7 @@
   </ItemGroup>
 
   <ItemGroup>
-    <ProjectReference Include="..\VSOP2013.NET.DataStructures\VSOP2013.NET.DataStructures.csproj" />
+    <ProjectReference Include="..\VSOP2013.NET\VSOP2013.NET.csproj" />
   </ItemGroup>
 
 

DataConverter/DataReader.cs

@@ -191,17 +191,12 @@ private static Term ReadTerm(ReadOnlySpan<char> line, ref Term T)
                 ptr += 3;
             }
 
-            ci = double.Parse(line[ptr..(ptr + 20)].Trim());
-            ptr += 20;
-            ptr++;
-            ci *= Math.Pow(10, double.Parse(line[ptr..(ptr + 3)].Trim()));
-            ptr += 3;
+            // Parse mantissa and exponent together: mantissa * 10^exp would round twice.
+            ci = ParseScaled(line[ptr..(ptr + 20)], line[(ptr + 21)..(ptr + 24)]);
+            ptr += 24;
             T.ss = ci;
 
-            ci = double.Parse(line[ptr..(ptr + 20)].Trim());
-            ptr += 20;
-            ptr++;
-            ci *= Math.Pow(10, double.Parse(line[ptr..(ptr + 3)].Trim()));
+            ci = ParseScaled(line[ptr..(ptr + 20)], line[(ptr + 21)..(ptr + 24)]);
 
             T.cc = ci;
 
@@ -212,5 +207,11 @@ private static Term ReadTerm(ReadOnlySpan<char> line, ref Term T)
             }
             return T;
         }
+
+        private static double ParseScaled(ReadOnlySpan<char> mantissa, ReadOnlySpan<char> exponent)
+        {
+            string text = $"{mantissa.Trim()}e{int.Parse(exponent.Trim())}";
+            return double.Parse(text, System.Globalization.CultureInfo.InvariantCulture);
+        }
     }
 }

Demo/PerfTest.cs

@@ -1,36 +1,50 @@
+// Licensed to the .NET Foundation under one or more agreements.
+// The .NET Foundation licenses this file to you under the MIT license.
+
 using System;
 using BenchmarkDotNet.Attributes;
+using BenchmarkDotNet.Configs;
 using BenchmarkDotNet.Jobs;
 using VSOP2013;
 
 namespace Demo
 {
+    /// <summary>
+    /// GetPlanetPosition(EMB) for {Scalar, Vector256} x {Calculator (double), CalculatorF (float)}.
+    /// Ratios are computed against the double calculator within the same SIMD mode.
+    /// </summary>
     [SimpleJob(RuntimeMoniker.Net10_0)]
     [MemoryDiagnoser]
+    [GroupBenchmarksBy(BenchmarkLogicalGroupRule.ByParams)]
     public class PerfTest
     {
-        private Calculator _vsop;
+        private Calculator _double = null!;
+        private CalculatorF _float = null!;
 
-        private DateTime _dt;
         private VSOPTime _vTime;
 
-        public PerfTest()
-        { _vsop = new Calculator(); }
+        [Params(SimdMode.Scalar, SimdMode.Vector256)]
+        public SimdMode Mode { get; set; }
 
         [GlobalSetup]
         public void init()
         {
-            _dt = DateTime.Now.ToUniversalTime();
-            _dt.ToUniversalTime();
-            _dt = _dt.AddSeconds(-69.184);
-            _vTime = new VSOPTime(_dt, TimeFrame.UTC);
+            // Double: full series. Float: default truncation (1e-8), the intended configuration.
+            _double = new Calculator(Mode);
+            _float = new CalculatorF(Mode);
+            _vTime = new VSOPTime(DateTime.UtcNow.AddSeconds(-69.184), TimeFrame.UTC);
         }
 
         [Benchmark(Baseline = true)]
-        public VSOPResult Compute()
+        public VSOPResult Double()
+        {
+            return _double.GetPlanetPosition(VSOPBody.EMB, _vTime);
+        }
+
+        [Benchmark]
+        public VSOPResult Float()
         {
-            var ell = _vsop.GetPlanetPosition(VSOPBody.EMB, _vTime);
-            return ell;
+            return _float.GetPlanetPosition(VSOPBody.EMB, _vTime);
         }
     }
 }

Demo/Program.cs

@@ -1,7 +1,11 @@
+// Licensed to the .NET Foundation under one or more agreements.
+// The .NET Foundation licenses this file to you under the MIT license.
+
 using System;
 using System.Globalization;
 using System.Text.Json;
 using System.Threading.Tasks;
+using BenchmarkDotNet.Configs;
 using BenchmarkDotNet.Running;
 using VSOP2013;
 
@@ -38,7 +42,7 @@ private static async Task Main(string[] args)
 
             Console.WriteLine("Press Enter to Start Performance Test...");
             Console.ReadLine();
-            return;
+            //return;
 #if DEBUG
             var summary = BenchmarkRunner.Run<PerfTest>(new DebugBuildConfig());
 #else
@@ -131,12 +135,21 @@ public static void FormattedPrint(VSOPResult Result, VSOPTime vtime)
         private static void WriteColorLine(params object[] flow)
         {
             foreach (var element in flow)
+            {
                 if (element == null)
+                {
                     Console.ResetColor();
+                }
                 else if (element is ConsoleColor color)
+                {
                     Console.ForegroundColor = color;
+                }
                 else
+                {
                     Console.Write(element.ToString());
+                }
+            }
+
             Console.WriteLine();
             Console.ResetColor();
         }

Documentation/README.md

@@ -28,43 +28,33 @@ This is the best VSOP2013 library ever.
 
 ## Performance
 
-```
+`GetPlanetPosition(EMB)` on an i9-12950HX, .NET 10:
 
-BenchmarkDotNet v0.13.11, Windows 11 (10.0.22631.4751/23H2/2023Update/SunValley3)
-12th Gen Intel Core i9-12950HX, 1 CPU, 16 logical and 16 physical cores
-.NET SDK 9.0.102
-  [Host]   : .NET 8.0.12 (8.0.1224.60305), X64 RyuJIT AVX2 [AttachedDebugger]
-  .NET 8.0 : .NET 8.0.12 (8.0.1224.60305), X64 RyuJIT AVX2
-  .NET 9.0 : .NET 9.0.1 (9.0.124.61010), X64 RyuJIT AVX2
+| Calculator | Scalar | Vector256 |
+| ---------- | -----: | --------: |
+| Calculator (double, full series) | 1,241 μs | 243 μs |
+| CalculatorF (float, truncated) | 12.6 μs | 2.4 μs |
 
-
-```
-
-| Method  | Job      | Runtime  |     Mean |     Error |    StdDev | Ratio | Allocated | Alloc Ratio |
-| ------- | -------- | -------- | -------: | --------: | --------: | ----: | --------: | ----------: |
-| Compute | .NET 8.0 | .NET 8.0 | 1.473 ms | 0.0150 ms | 0.0133 ms |  1.00 |   2.91 KB |        1.00 |
-|         |          |          |          |           |           |       |           |             |
-| Compute | .NET 9.0 | .NET 9.0 | 1.431 ms | 0.0196 ms | 0.0174 ms |  1.00 |    2.9 KB |        1.00 |
-
-Note: .NET 8 occurs performance regression due to RyuJIT bugs. [Detail Here](https://github.com/dotnet/runtime/issues/95954#issuecomment-1956661569)
+See [Performance](vsop2013.net/performance.md) for details.
 
 ## Features
 
 1. Use VSOPResult class to manage calculate results.
 2. Use VSOPTime class to manage time.\
    Easy to convert time by calling `VSOPTime.UTC`, `VSOPTime.TAI`, `VSOPTime.TDB`
-3. Very high performance per solution
+3. Very high performance per solution. SIMD (`Vector128/256/512`) evaluation with automatic width detection, plus multi-thread across variables.
 4. Useful Utility class. Convert Elliptic coordinates to cartesian and spherical
 5. Async Api included.
 6. precalculation on **φ** in data loading, which gives 20%+ speed up of calculation.
-7. Use  [MessagePack](https://github.com/neuecc/MessagePack-CSharp) for binary serialize.\
+7. Optional term truncation (`Calculator(mode, truncation)`) to trade precision for speed.
+8. `CalculatorF`: single precision calculator for visualization, 100 to 500 times faster.
+9. Validated against the official `VSOP2013.ctl` reference values (max error 1e-10).
+10. Use [MemoryPack](https://github.com/Cysharp/MemoryPack) for binary serialize.\
    Initialization time becomes less than 10% of previous version.
-8. FastLZMA2 compression on data. \~300Mb -> \~50MB.
+11. Brotli compression on data. \~300Mb -> \~43MB.
 
 
 
 ### Reference
 
-[https://github.com/neuecc/MessagePack-CSharp](https://github.com/neuecc/MessagePack-CSharp)
-
-[https://github.com/kingsznhone/FastLZMA2Net](https://github.com/kingsznhone/FastLZMA2Net)
+[https://github.com/Cysharp/MemoryPack](https://github.com/Cysharp/MemoryPack)

Documentation/SUMMARY.md

@@ -4,6 +4,7 @@
   * [Quick Start](vsop2013.net/quick-start.md)
   * [Change Logs](vsop2013.net/change-logs.md)
   * [Requirements](vsop2013.net/requirements.md)
+  * [Performance](vsop2013.net/performance.md)
 * [VSOP2013](vsop2013/README.md)
   * [Original Readme](vsop2013/original-readme/README.md)
     * [Introduction](vsop2013/original-readme/introduction.md)
@@ -43,5 +44,6 @@
     * [VSOPResult\_XYZ Class](api-reference/vsop2013-namespace/vsopresult-class/vsopresult_xyz-class.md)
     * [VSOPResult\_LBR Class](api-reference/vsop2013-namespace/vsopresult-class/vsopresult_lbr-class.md)
   * [VSOPTime Class](api-reference/vsop2013-namespace/vsoptime-class.md)
+  * [CalculatorF Class](api-reference/vsop2013-namespace/calculatorf-class.md)
   * [Data Structs](api-reference/vsop2013-namespace/data-structs.md)
   * [Enums](api-reference/vsop2013-namespace/enums.md)

Documentation/api-reference/vsop2013-namespace/README.md

@@ -3,6 +3,7 @@
 ## Classes
 
 * [Calculator](calculator-class/)
+* [CalculatorF](calculatorf-class.md)
 * [VSOPResult](vsopresult-class/)
 * [VSOPResult\_ELL](vsopresult-class/vsopresult\_ell-class.md)
 * [VSOPResult\_LBR](vsopresult-class/vsopresult\_lbr-class.md)
@@ -23,6 +24,7 @@
 * [CoordinatesReference](enums.md#coordinatesreference)
 * [CoordinatesType](enums.md#coordinatestype)
 * [ReferenceFrame](enums.md#referenceframe)
+* [SimdMode](enums.md#simdmode)
 * [VSOPBody](enums.md#fields)
 * [VSOPVariable](enums.md#fields-2)
 * [TimeFrame](enums.md#fields-4)

Documentation/api-reference/vsop2013-namespace/calculator-class/README.md

@@ -2,24 +2,66 @@
 
 ## Definition
 
-Provide methods to calculate planet position. Full 6 variable or one of them.
+Full precision VSOP2013 calculator (double precision). It matches the official `VSOP2013.ctl` reference values to 1e-10.
+
+Full 6 variable or one of them.
 
 ```csharp
-public class Calculator
+public partial class Calculator
 ```
 
-## Constructors <a href="#constructors" id="constructors"></a>
+An instance is read-only after construction and safe to share across threads. Create it once and reuse it; construction loads and decompresses the data file.
+
+For fast, low precision visualization use [CalculatorF](../calculatorf-class.md).
 
-Initializes a new instance of the Calculator class.
+## Constructors <a href="#constructors" id="constructors"></a>
 
 ```csharp
-public VSOP2013 ();
+public Calculator();
+public Calculator(SimdMode mode, double truncation = 0d);
 ```
 
+**Parameters**
+
+`mode` [SimdMode](../enums.md#simdmode)
+
+Vector width used to evaluate the series. `SimdMode.Auto` picks the widest hardware accelerated width.
+
+`truncation` double
+
+Truncation level ρ as defined by VSOP2013. Terms with `sqrt(C² + S²) < ρ` are dropped when the data is loaded (au for A, rad for L, unitless for K/H/Q/P). `0` (default) keeps the full series.
+
+**Exceptions**
+
+`PlatformNotSupportedException`: the requested vector width is not hardware accelerated.
+
+`ArgumentOutOfRangeException`: `truncation` is negative or not finite.
+
+## Properties <a href="#properties" id="properties"></a>
+
+`Mode` [SimdMode](../enums.md#simdmode)
+
+Vector width actually used by this instance after hardware detection. It is never `Auto`.
+
+`Truncation` double
+
+Truncation level applied to this instance.
+
 ## Methods <a href="#methods" id="methods"></a>
 
 * [GetPlanetPosition(VSOPBody body, VSOPTime time)](getplanetposition.md)
 * [GetPlanetPositionAsync(VSOPBody body, VSOPTime time)](getplanetpositionasync.md)
-* [GetVariable(VSOPBody body, int iv, VSOPTime time)](getvariable.md)
-* [GetVariableAsync(VSOPBody body, int iv, VSOPTime time)](getvariableasync.md)
+* [GetVariable(VSOPBody body, VSOPVariable variable, VSOPTime time)](getvariable.md)
+* [GetVariableAsync(VSOPBody body, VSOPVariable variable, VSOPTime time)](getvariableasync.md)
+
+## Example
+
+```csharp
+// Default: full series, widest SIMD width.
+Calculator vsop = new Calculator();
+
+// Scalar evaluation, drop terms smaller than 1e-14.
+Calculator scalar = new Calculator(SimdMode.Scalar, truncation: 1e-14);
+Console.WriteLine(scalar.Mode); // Scalar
+```
 

Documentation/api-reference/vsop2013-namespace/calculator-class/getplanetposition.md

@@ -1,6 +1,6 @@
 # GetPlanetPosition
 
-Calculate all 6 variable of a planet.
+Calculate all 6 variable of a planet. The six variables are evaluated in parallel.
 
 ```csharp
 public VSOPResult_ELL GetPlanetPosition(VSOPBody body, VSOPTime time)

Documentation/api-reference/vsop2013-namespace/calculator-class/getplanetpositionasync.md

@@ -3,7 +3,7 @@
 Calculate all 6 variable of a planet.
 
 ```csharp
-public Task<VSOPResult_ELL> GetPlanetPositionAsync(VSOPBody body, VSOPTime time)
+public ValueTask<VSOPResult_ELL> GetPlanetPositionAsync(VSOPBody body, VSOPTime time)
 ```
 
 **Parameters**
@@ -22,7 +22,7 @@ Exclusive time class in VSOP2013.
 
 **Return**
 
-`Task<`[`VSOPResult_ELL`](../vsopresult-class/vsopresult_ell-class.md)`>`
+`ValueTask<`[`VSOPResult_ELL`](../vsopresult-class/vsopresult_ell-class.md)`>`
 
 Full Result with 6 variable of Elliptic Coordinates.