diff options
Diffstat (limited to '')
-rw-r--r-- | OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs | 214 |
1 files changed, 214 insertions, 0 deletions
diff --git a/OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs b/OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs new file mode 100644 index 0000000..3eae46b --- /dev/null +++ b/OpenSim/Region/Terrain.BasicTerrain/libTerrainBSD/Channel/Manipulators/NavierStokes.cs | |||
@@ -0,0 +1,214 @@ | |||
1 | /* | ||
2 | * Copyright (c) Contributors, http://www.openmetaverse.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 | |||
29 | using System; | ||
30 | using System.Collections.Generic; | ||
31 | using System.Text; | ||
32 | |||
33 | namespace libTerrain | ||
34 | { | ||
35 | partial class Channel | ||
36 | { | ||
37 | // Navier Stokes Algorithms ported from | ||
38 | // "Real-Time Fluid Dynamics for Games" by Jos Stam. | ||
39 | // presented at GDC 2003. | ||
40 | |||
41 | // Poorly ported from C++. (I gave up making it properly native somewhere after nsSetBnd) | ||
42 | |||
43 | private static int nsIX(int i, int j, int N) | ||
44 | { | ||
45 | return ((i) + (N + 2) * (j)); | ||
46 | } | ||
47 | |||
48 | private static void nsSwap(ref double x0, ref double x) | ||
49 | { | ||
50 | double tmp = x0; | ||
51 | x0 = x; | ||
52 | x = tmp; | ||
53 | } | ||
54 | |||
55 | private static void nsSwap(ref double[] x0, ref double[] x) | ||
56 | { | ||
57 | double[] tmp = x0; | ||
58 | x0 = x; | ||
59 | x = tmp; | ||
60 | } | ||
61 | |||
62 | private void nsAddSource(int N, ref double[] x, ref double[] s, double dt) | ||
63 | { | ||
64 | int i; | ||
65 | int size = (N + 2) * (N + 2); | ||
66 | for (i = 0; i < size; i++) | ||
67 | { | ||
68 | x[i] += dt * s[i]; | ||
69 | } | ||
70 | } | ||
71 | |||
72 | private void nsSetBnd(int N, int b, ref double[] x) | ||
73 | { | ||
74 | int i; | ||
75 | for (i = 0; i <= N; i++) | ||
76 | { | ||
77 | x[nsIX(0, i, N)] = b == 1 ? -x[nsIX(1, i, N)] : x[nsIX(1, i, N)]; | ||
78 | x[nsIX(0, N + 1, N)] = b == 1 ? -x[nsIX(N, i, N)] : x[nsIX(N, i, N)]; | ||
79 | x[nsIX(i, 0, N)] = b == 2 ? -x[nsIX(i, 1, N)] : x[nsIX(i, 1, N)]; | ||
80 | x[nsIX(i, N + 1, N)] = b == 2 ? -x[nsIX(i, N, N)] : x[nsIX(i, N, N)]; | ||
81 | } | ||
82 | x[nsIX(0, 0, N)] = 0.5f * (x[nsIX(1, 0, N)] + x[nsIX(0, 1, N)]); | ||
83 | x[nsIX(0, N + 1, N)] = 0.5f * (x[nsIX(1, N + 1, N)] + x[nsIX(0, N, N)]); | ||
84 | x[nsIX(N + 1, 0, N)] = 0.5f * (x[nsIX(N, 0, N)] + x[nsIX(N + 1, 1, N)]); | ||
85 | x[nsIX(N + 1, N + 1, N)] = 0.5f * (x[nsIX(N, N + 1, N)] + x[nsIX(N + 1, N, N)]); | ||
86 | } | ||
87 | |||
88 | private void nsLinSolve(int N, int b, ref double[] x, ref double[] x0, double a, double c) | ||
89 | { | ||
90 | int i, j; | ||
91 | for (i = 1; i <= N; i++) | ||
92 | { | ||
93 | for (j = 1; j <= N; j++) | ||
94 | { | ||
95 | x[nsIX(i, j, N)] = (x0[nsIX(i, j, N)] + a * | ||
96 | (x[nsIX(i - 1, j, N)] + | ||
97 | x[nsIX(i + 1, j, N)] + | ||
98 | x[nsIX(i, j - 1, N)] + x[nsIX(i, j + 1, N)]) | ||
99 | ) / c; | ||
100 | } | ||
101 | } | ||
102 | |||
103 | nsSetBnd(N, b, ref x); | ||
104 | } | ||
105 | |||
106 | private void nsDiffuse(int N, int b, ref double[] x, ref double[] x0, double diff, double dt) | ||
107 | { | ||
108 | double a = dt * diff * N * N; | ||
109 | nsLinSolve(N, b, ref x, ref x0, a, 1 + 4 * a); | ||
110 | } | ||
111 | |||
112 | private void nsAdvect(int N, int b, ref double[] d, ref double[] d0, ref double[] u, ref double[] v, double dt) | ||
113 | { | ||
114 | int i, j, i0, j0, i1, j1; | ||
115 | double x, y, s0, t0, s1, t1, dt0; | ||
116 | |||
117 | dt0 = dt * N; | ||
118 | |||
119 | for (i = 1; i <= N; i++) | ||
120 | { | ||
121 | for (j = 1; j <= N; j++) | ||
122 | { | ||
123 | x = i - dt0 * u[nsIX(i, j, N)]; | ||
124 | y = j - dt0 * v[nsIX(i, j, N)]; | ||
125 | |||
126 | if (x < 0.5) | ||
127 | x = 0.5; | ||
128 | if (x > N + 0.5) | ||
129 | x = N + 0.5; | ||
130 | i0 = (int)x; | ||
131 | i1 = i0 + 1; | ||
132 | |||
133 | if (y < 0.5) | ||
134 | y = 0.5; | ||
135 | if (y > N + 0.5) | ||
136 | y = N + 0.5; | ||
137 | j0 = (int)y; | ||
138 | j1 = j0 + 1; | ||
139 | |||
140 | s1 = x - i0; | ||
141 | s0 = 1 - s1; | ||
142 | t1 = y - j0; | ||
143 | t0 = 1 - t1; | ||
144 | |||
145 | d[nsIX(i, j, N)] = s0 * (t0 * d0[nsIX(i0, j0, N)] + t1 * d0[nsIX(i0, j1, N)]) + | ||
146 | s1 * (t0 * d0[nsIX(i1, j0, N)] + t1 * d0[nsIX(i1, j1, N)]); | ||
147 | } | ||
148 | } | ||
149 | |||
150 | nsSetBnd(N, b, ref d); | ||
151 | } | ||
152 | |||
153 | public void nsProject(int N, ref double[] u, ref double[] v, ref double[] p, ref double[] div) | ||
154 | { | ||
155 | int i, j; | ||
156 | |||
157 | for (i = 1; i <= N; i++) | ||
158 | { | ||
159 | for (j = 1; j <= N; j++) | ||
160 | { | ||
161 | div[nsIX(i, j, N)] = -0.5 * (u[nsIX(i + 1, j, N)] - u[nsIX(i - 1, j, N)] + v[nsIX(i, j + 1, N)] - v[nsIX(i, j - 1, N)]) / N; | ||
162 | p[nsIX(i, j, N)] = 0; | ||
163 | } | ||
164 | } | ||
165 | |||
166 | nsSetBnd(N, 0, ref div); | ||
167 | nsSetBnd(N, 0, ref p); | ||
168 | |||
169 | nsLinSolve(N, 0, ref p, ref div, 1, 4); | ||
170 | |||
171 | for (i = 1; i <= N; i++) | ||
172 | { | ||
173 | for (j = 1; j <= N; j++) | ||
174 | { | ||
175 | u[nsIX(i, j, N)] -= 0.5 * N * (p[nsIX(i + 1, j, N)] - p[nsIX(i - 1, j, N)]); | ||
176 | v[nsIX(i, j, N)] -= 0.5 * N * (p[nsIX(i, j + 1, N)] - p[nsIX(i, j - 1, N)]); | ||
177 | } | ||
178 | } | ||
179 | |||
180 | nsSetBnd(N, 1, ref u); | ||
181 | nsSetBnd(N, 2, ref v); | ||
182 | } | ||
183 | |||
184 | private void nsDensStep(int N, ref double[] x, ref double[] x0, ref double[] u, ref double[] v, double diff, double dt) | ||
185 | { | ||
186 | nsAddSource(N, ref x, ref x0, dt); | ||
187 | nsSwap(ref x0, ref x); | ||
188 | nsDiffuse(N, 0, ref x, ref x0, diff, dt); | ||
189 | nsSwap(ref x0, ref x); | ||
190 | nsAdvect(N, 0, ref x, ref x0, ref u, ref v, dt); | ||
191 | } | ||
192 | |||
193 | private void nsVelStep(int N, ref double[] u, ref double[] v, ref double[] u0, ref double[] v0, double visc, 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 | public void navierSimulate() | ||
210 | { | ||
211 | |||
212 | } | ||
213 | } | ||
214 | } \ No newline at end of file | ||