diff options
Diffstat (limited to 'OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs')
-rw-r--r-- | OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs | 307 |
1 files changed, 0 insertions, 307 deletions
diff --git a/OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs b/OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs deleted file mode 100644 index 8f12637..0000000 --- a/OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs +++ /dev/null | |||
@@ -1,307 +0,0 @@ | |||
1 | /* | ||
2 | * Copyright (c) Contributors, http://opensimulator.org/ | ||
3 | * See CONTRIBUTORS.TXT for a full list of copyright holders. | ||
4 | * | ||
5 | * Redistribution and use in source and binary forms, with or without | ||
6 | * modification, are permitted provided that the following conditions are met: | ||
7 | * * Redistributions of source code must retain the above copyright | ||
8 | * notice, this list of conditions and the following disclaimer. | ||
9 | * * Redistributions in binary form must reproduce the above copyright | ||
10 | * notice, this list of conditions and the following disclaimer in the | ||
11 | * documentation and/or other materials provided with the distribution. | ||
12 | * * Neither the name of the OpenSim Project nor the | ||
13 | * names of its contributors may be used to endorse or promote products | ||
14 | * derived from this software without specific prior written permission. | ||
15 | * | ||
16 | * THIS SOFTWARE IS PROVIDED BY THE DEVELOPERS ``AS IS'' AND ANY | ||
17 | * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED | ||
18 | * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE | ||
19 | * DISCLAIMED. IN NO EVENT SHALL THE CONTRIBUTORS BE LIABLE FOR ANY | ||
20 | * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES | ||
21 | * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; | ||
22 | * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND | ||
23 | * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT | ||
24 | * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS | ||
25 | * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. | ||
26 | */ | ||
27 | |||
28 | namespace libTerrain | ||
29 | { | ||
30 | partial class Channel | ||
31 | { | ||
32 | // Navier Stokes Algorithms ported from | ||
33 | // "Real-Time Fluid Dynamics for Games" by Jos Stam. | ||
34 | // presented at GDC 2003. | ||
35 | |||
36 | // Poorly ported from C++. (I gave up making it properly native somewhere after nsSetBnd) | ||
37 | |||
38 | private static int nsIX(int i, int j, int N) | ||
39 | { | ||
40 | return ((i) + (N + 2)*(j)); | ||
41 | } | ||
42 | |||
43 | // TODO: unused | ||
44 | // private static void nsSwap(ref double x0, ref double x) | ||
45 | // { | ||
46 | // double tmp = x0; | ||
47 | // x0 = x; | ||
48 | // x = tmp; | ||
49 | // } | ||
50 | |||
51 | private static void nsSwap(ref double[] x0, ref double[] x) | ||
52 | { | ||
53 | double[] tmp = x0; | ||
54 | x0 = x; | ||
55 | x = tmp; | ||
56 | } | ||
57 | |||
58 | private void nsAddSource(int N, ref double[] x, ref double[] s, double dt) | ||
59 | { | ||
60 | int i; | ||
61 | int size = (N + 2)*(N + 2); | ||
62 | for (i = 0; i < size; i++) | ||
63 | { | ||
64 | x[i] += dt*s[i]; | ||
65 | } | ||
66 | } | ||
67 | |||
68 | private void nsSetBnd(int N, int b, ref double[] x) | ||
69 | { | ||
70 | int i; | ||
71 | for (i = 0; i <= N; i++) | ||
72 | { | ||
73 | x[nsIX(0, i, N)] = b == 1 ? -x[nsIX(1, i, N)] : x[nsIX(1, i, N)]; | ||
74 | x[nsIX(0, N + 1, N)] = b == 1 ? -x[nsIX(N, i, N)] : x[nsIX(N, i, N)]; | ||
75 | x[nsIX(i, 0, N)] = b == 2 ? -x[nsIX(i, 1, N)] : x[nsIX(i, 1, N)]; | ||
76 | x[nsIX(i, N + 1, N)] = b == 2 ? -x[nsIX(i, N, N)] : x[nsIX(i, N, N)]; | ||
77 | } | ||
78 | x[nsIX(0, 0, N)] = 0.5f*(x[nsIX(1, 0, N)] + x[nsIX(0, 1, N)]); | ||
79 | x[nsIX(0, N + 1, N)] = 0.5f*(x[nsIX(1, N + 1, N)] + x[nsIX(0, N, N)]); | ||
80 | x[nsIX(N + 1, 0, N)] = 0.5f*(x[nsIX(N, 0, N)] + x[nsIX(N + 1, 1, N)]); | ||
81 | x[nsIX(N + 1, N + 1, N)] = 0.5f*(x[nsIX(N, N + 1, N)] + x[nsIX(N + 1, N, N)]); | ||
82 | } | ||
83 | |||
84 | private void nsLinSolve(int N, int b, ref double[] x, ref double[] x0, double a, double c) | ||
85 | { | ||
86 | int i, j; | ||
87 | for (i = 1; i <= N; i++) | ||
88 | { | ||
89 | for (j = 1; j <= N; j++) | ||
90 | { | ||
91 | x[nsIX(i, j, N)] = (x0[nsIX(i, j, N)] + a* | ||
92 | (x[nsIX(i - 1, j, N)] + | ||
93 | x[nsIX(i + 1, j, N)] + | ||
94 | x[nsIX(i, j - 1, N)] + x[nsIX(i, j + 1, N)]) | ||
95 | )/c; | ||
96 | } | ||
97 | } | ||
98 | |||
99 | nsSetBnd(N, b, ref x); | ||
100 | } | ||
101 | |||
102 | private void nsDiffuse(int N, int b, ref double[] x, ref double[] x0, double diff, double dt) | ||
103 | { | ||
104 | double a = dt*diff*N*N; | ||
105 | nsLinSolve(N, b, ref x, ref x0, a, 1 + 4*a); | ||
106 | } | ||
107 | |||
108 | private void nsAdvect(int N, int b, ref double[] d, ref double[] d0, ref double[] u, ref double[] v, double dt) | ||
109 | { | ||
110 | int i, j, i0, j0, i1, j1; | ||
111 | double x, y, s0, t0, s1, t1, dt0; | ||
112 | |||
113 | dt0 = dt*N; | ||
114 | |||
115 | for (i = 1; i <= N; i++) | ||
116 | { | ||
117 | for (j = 1; j <= N; j++) | ||
118 | { | ||
119 | x = i - dt0*u[nsIX(i, j, N)]; | ||
120 | y = j - dt0*v[nsIX(i, j, N)]; | ||
121 | |||
122 | if (x < 0.5) | ||
123 | x = 0.5; | ||
124 | if (x > N + 0.5) | ||
125 | x = N + 0.5; | ||
126 | i0 = (int) x; | ||
127 | i1 = i0 + 1; | ||
128 | |||
129 | if (y < 0.5) | ||
130 | y = 0.5; | ||
131 | if (y > N + 0.5) | ||
132 | y = N + 0.5; | ||
133 | j0 = (int) y; | ||
134 | j1 = j0 + 1; | ||
135 | |||
136 | s1 = x - i0; | ||
137 | s0 = 1 - s1; | ||
138 | t1 = y - j0; | ||
139 | t0 = 1 - t1; | ||
140 | |||
141 | d[nsIX(i, j, N)] = s0*(t0*d0[nsIX(i0, j0, N)] + t1*d0[nsIX(i0, j1, N)]) + | ||
142 | s1*(t0*d0[nsIX(i1, j0, N)] + t1*d0[nsIX(i1, j1, N)]); | ||
143 | } | ||
144 | } | ||
145 | |||
146 | nsSetBnd(N, b, ref d); | ||
147 | } | ||
148 | |||
149 | public void nsProject(int N, ref double[] u, ref double[] v, ref double[] p, ref double[] div) | ||
150 | { | ||
151 | int i, j; | ||
152 | |||
153 | for (i = 1; i <= N; i++) | ||
154 | { | ||
155 | for (j = 1; j <= N; j++) | ||
156 | { | ||
157 | div[nsIX(i, j, N)] = -0.5* | ||
158 | (u[nsIX(i + 1, j, N)] - u[nsIX(i - 1, j, N)] + v[nsIX(i, j + 1, N)] - | ||
159 | v[nsIX(i, j - 1, N)])/N; | ||
160 | p[nsIX(i, j, N)] = 0; | ||
161 | } | ||
162 | } | ||
163 | |||
164 | nsSetBnd(N, 0, ref div); | ||
165 | nsSetBnd(N, 0, ref p); | ||
166 | |||
167 | nsLinSolve(N, 0, ref p, ref div, 1, 4); | ||
168 | |||
169 | for (i = 1; i <= N; i++) | ||
170 | { | ||
171 | for (j = 1; j <= N; j++) | ||
172 | { | ||
173 | u[nsIX(i, j, N)] -= 0.5*N*(p[nsIX(i + 1, j, N)] - p[nsIX(i - 1, j, N)]); | ||
174 | v[nsIX(i, j, N)] -= 0.5*N*(p[nsIX(i, j + 1, N)] - p[nsIX(i, j - 1, N)]); | ||
175 | } | ||
176 | } | ||
177 | |||
178 | nsSetBnd(N, 1, ref u); | ||
179 | nsSetBnd(N, 2, ref v); | ||
180 | } | ||
181 | |||
182 | private void nsDensStep(int N, ref double[] x, ref double[] x0, ref double[] u, ref double[] v, double diff, | ||
183 | double dt) | ||
184 | { | ||
185 | nsAddSource(N, ref x, ref x0, dt); | ||
186 | nsSwap(ref x0, ref x); | ||
187 | nsDiffuse(N, 0, ref x, ref x0, diff, dt); | ||
188 | nsSwap(ref x0, ref x); | ||
189 | nsAdvect(N, 0, ref x, ref x0, ref u, ref v, dt); | ||
190 | } | ||
191 | |||
192 | private void nsVelStep(int N, ref double[] u, ref double[] v, ref double[] u0, ref double[] v0, double visc, | ||
193 | double dt) | ||
194 | { | ||
195 | nsAddSource(N, ref u, ref u0, dt); | ||
196 | nsAddSource(N, ref v, ref v0, dt); | ||
197 | nsSwap(ref u0, ref u); | ||
198 | nsDiffuse(N, 1, ref u, ref u0, visc, dt); | ||
199 | nsSwap(ref v0, ref v); | ||
200 | nsDiffuse(N, 2, ref v, ref v0, visc, dt); | ||
201 | nsProject(N, ref u, ref v, ref u0, ref v0); | ||
202 | nsSwap(ref u0, ref u); | ||
203 | nsSwap(ref v0, ref v); | ||
204 | nsAdvect(N, 1, ref u, ref u0, ref u0, ref v0, dt); | ||
205 | nsAdvect(N, 2, ref v, ref v0, ref u0, ref v0, dt); | ||
206 | nsProject(N, ref u, ref v, ref u0, ref v0); | ||
207 | } | ||
208 | |||
209 | private void nsBufferToDoubles(ref double[] dens, int N, ref double[,] doubles) | ||
210 | { | ||
211 | int i; | ||
212 | int j; | ||
213 | |||
214 | for (i = 1; i <= N; i++) | ||
215 | { | ||
216 | for (j = 1; j <= N; j++) | ||
217 | { | ||
218 | doubles[i - 1, j - 1] = dens[nsIX(i, j, N)]; | ||
219 | } | ||
220 | } | ||
221 | } | ||
222 | |||
223 | private void nsDoublesToBuffer(double[,] doubles, int N, ref double[] dens) | ||
224 | { | ||
225 | int i; | ||
226 | int j; | ||
227 | |||
228 | for (i = 1; i <= N; i++) | ||
229 | { | ||
230 | for (j = 1; j <= N; j++) | ||
231 | { | ||
232 | dens[nsIX(i, j, N)] = doubles[i - 1, j - 1]; | ||
233 | } | ||
234 | } | ||
235 | } | ||
236 | |||
237 | private void nsSimulate(int N, int rounds, double dt, double diff, double visc) | ||
238 | { | ||
239 | int size = (N*2)*(N*2); | ||
240 | |||
241 | double[] u = new double[size]; // Force, X axis | ||
242 | double[] v = new double[size]; // Force, Y axis | ||
243 | double[] u_prev = new double[size]; | ||
244 | double[] v_prev = new double[size]; | ||
245 | double[] dens = new double[size]; | ||
246 | double[] dens_prev = new double[size]; | ||
247 | |||
248 | nsDoublesToBuffer(map, N, ref dens); | ||
249 | nsDoublesToBuffer(map, N, ref dens_prev); | ||
250 | |||
251 | for (int i = 0; i < rounds; i++) | ||
252 | { | ||
253 | u_prev = u; | ||
254 | v_prev = v; | ||
255 | dens_prev = dens; | ||
256 | |||
257 | nsVelStep(N, ref u, ref v, ref u_prev, ref v_prev, visc, dt); | ||
258 | nsDensStep(N, ref dens, ref dens_prev, ref u, ref v, diff, dt); | ||
259 | } | ||
260 | |||
261 | nsBufferToDoubles(ref dens, N, ref map); | ||
262 | } | ||
263 | |||
264 | /// <summary> | ||
265 | /// Performs computational fluid dynamics on a channel | ||
266 | /// </summary> | ||
267 | /// <param name="rounds">The number of steps to perform (Recommended: 20)</param> | ||
268 | /// <param name="dt">Delta Time - The time between steps (Recommended: 0.1)</param> | ||
269 | /// <param name="diff">Fluid diffusion rate (Recommended: 0.0)</param> | ||
270 | /// <param name="visc">Fluid viscosity (Recommended: 0.0)</param> | ||
271 | public void navierStokes(int rounds, double dt, double diff, double visc) | ||
272 | { | ||
273 | nsSimulate(h, rounds, dt, diff, visc); | ||
274 | } | ||
275 | |||
276 | public void navierStokes(int rounds, double dt, double diff, double visc, ref double[,] uret, ref double[,] vret) | ||
277 | { | ||
278 | int N = h; | ||
279 | |||
280 | int size = (N*2)*(N*2); | ||
281 | |||
282 | double[] u = new double[size]; // Force, X axis | ||
283 | double[] v = new double[size]; // Force, Y axis | ||
284 | double[] u_prev = new double[size]; | ||
285 | double[] v_prev = new double[size]; | ||
286 | double[] dens = new double[size]; | ||
287 | double[] dens_prev = new double[size]; | ||
288 | |||
289 | nsDoublesToBuffer(map, N, ref dens); | ||
290 | nsDoublesToBuffer(map, N, ref dens_prev); | ||
291 | |||
292 | for (int i = 0; i < rounds; i++) | ||
293 | { | ||
294 | u_prev = u; | ||
295 | v_prev = v; | ||
296 | dens_prev = dens; | ||
297 | |||
298 | nsVelStep(N, ref u, ref v, ref u_prev, ref v_prev, visc, dt); | ||
299 | nsDensStep(N, ref dens, ref dens_prev, ref u, ref v, diff, dt); | ||
300 | } | ||
301 | |||
302 | nsBufferToDoubles(ref u, N, ref uret); | ||
303 | nsBufferToDoubles(ref v, N, ref vret); | ||
304 | nsBufferToDoubles(ref dens, N, ref map); | ||
305 | } | ||
306 | } | ||
307 | } | ||