1 /**
2     Inochi2D Physics Subsystem
3 
4     Copyright: 
5         Copyright © 2020-2026, Inochi2D Project
6     
7     License:
8         $(LINK2 https://github.com/Inochi2D/inochi2d/blob/main/LICENSE, BSD 2-clause License)
9     
10     Authors:
11         Hoshino Lina
12 */
13 module inochi2d.core.phys.system;
14 import inochi2d;
15 import numath;
16 import numem;
17 import nulib;
18 
19 abstract
20 class PhysicsSystem : NuObject {
21 private:
22 @nogc:
23     map!(float*, size_t) variableMap;
24     vector!(float*) refs;
25 
26     // Derivatives
27     float[] kstate;
28     float[] k0;
29     float[] k1;
30     float[] k2;
31     float[] k3;
32     float[] k4;
33 
34     float t;
35 
36     /// Helper that resizes all of the temporary arrays.
37     void resize(size_t length) {
38         kstate = kstate.nu_resize(length);
39         k0 = k0.nu_resize(length);
40         k1 = k1.nu_resize(length);
41         k2 = k2.nu_resize(length);
42         k3 = k3.nu_resize(length);
43         k4 = k4.nu_resize(length);
44     }
45 
46 protected:
47     /**
48         Add a float variable to the simulation
49     */
50     size_t addVariable(float* var) {
51         size_t index = refs.length;
52 
53         variableMap[var] = index;
54         refs ~= var;
55 
56         this.resize(refs.length);
57         return index;
58     }
59 
60     /**
61         Add a vec2 variable to the simulation
62     */
63     size_t addVariable(vec2* var) {
64         size_t index = addVariable(&var.data[0]);
65         addVariable(&var.data[1]);
66         return index;
67     }
68 
69     /**
70         Set the derivative of a variable (solver input) by index
71     */
72     void setD(size_t index, float value) {
73         k0[index] = value;
74     }
75 
76     /**
77         Set the derivative of a float variable (solver input)
78     */
79     void setD(ref float var, float value) {
80         setD(variableMap[&var], value);
81     }
82 
83     /**
84         Set the derivative of a vec2 variable (solver input)
85     */
86     void setD(ref vec2 var, vec2 value) {
87         setD(var.data[0], value.x);
88         setD(var.data[1], value.y);
89     }
90 
91     float[] getState() {
92         foreach (idx, ptr; refs)
93             kstate[idx] = *ptr;
94         return kstate;
95     }
96 
97     void setState(float[] vals) {
98         foreach (idx, ptr; refs) {
99             *ptr = vals[idx];
100         }
101     }
102 
103     /**
104         Evaluate the simulation at a given time
105     */
106     abstract void eval(float t);
107 
108 public:
109 
110      ~this() {
111         variableMap.clear();
112         refs.clear();
113 
114         // Free derivatives.
115         nu_freea(kstate);
116         nu_freea(k0);
117         nu_freea(k1);
118         nu_freea(k2);
119         nu_freea(k3);
120         nu_freea(k4);
121     }
122 
123     /**
124         Run a simulation tick (Runge-Kutta method)
125     */
126     void tick(float h) {
127         float[] cur = getState();
128         k0[0 .. $] = 0;
129 
130         eval(t);
131         k1[0 .. $] = k0[0 .. $];
132 
133         foreach (i; 0 .. cur.length)
134             *refs[i] = cur[i] + h * k1[i] / 2f;
135         eval(t + h / 2f);
136         k2[0 .. $] = k0[0 .. $];
137 
138         foreach (i; 0 .. cur.length)
139             *refs[i] = cur[i] + h * k2[i] / 2f;
140         eval(t + h / 2f);
141         k3[0 .. $] = k0[0 .. $];
142 
143         foreach (i; 0 .. cur.length)
144             *refs[i] = cur[i] + h * k3[i];
145         eval(t + h);
146         k4[0 .. $] = k0[0 .. $];
147 
148         foreach (i; 0 .. cur.length) {
149             *refs[i] = cur[i] + h * (k1[i] + 2 * k2[i] + 2 * k3[i] + k4[i]) / 6f;
150             if (!isFinite(*refs[i])) {
151                 // Simulation failed, revert
152                 foreach (j; 0 .. cur.length)
153                     *refs[j] = cur[j];
154                 break;
155             }
156         }
157 
158         t += h;
159     }
160 
161     /**
162         Updates the anchor for the physics system
163     */
164     abstract void updateAnchor();
165 }