use psys->seed for smoke random number generator, increase size of MATHUTILS_TOT_CB and reduce float->double conversions.
This commit is contained in:
@@ -29,7 +29,7 @@ void Eigentred2(sEigenvalue& eval) {
|
|||||||
for (int k = 0; k < i; k++) {
|
for (int k = 0; k < i; k++) {
|
||||||
scale = scale + fabs(eval.d[k]);
|
scale = scale + fabs(eval.d[k]);
|
||||||
}
|
}
|
||||||
if (scale == 0.0) {
|
if (scale == 0.0f) {
|
||||||
eval.e[i] = eval.d[i-1];
|
eval.e[i] = eval.d[i-1];
|
||||||
for (int j = 0; j < i; j++) {
|
for (int j = 0; j < i; j++) {
|
||||||
eval.d[j] = eval.V[i-1][j];
|
eval.d[j] = eval.V[i-1][j];
|
||||||
@@ -96,7 +96,7 @@ void Eigentred2(sEigenvalue& eval) {
|
|||||||
eval.V[n-1][i] = eval.V[i][i];
|
eval.V[n-1][i] = eval.V[i][i];
|
||||||
eval.V[i][i] = 1.0;
|
eval.V[i][i] = 1.0;
|
||||||
float h = eval.d[i+1];
|
float h = eval.d[i+1];
|
||||||
if (h != 0.0) {
|
if (h != 0.0f) {
|
||||||
for (int k = 0; k <= i; k++) {
|
for (int k = 0; k <= i; k++) {
|
||||||
eval.d[k] = eval.V[k][i+1] / h;
|
eval.d[k] = eval.V[k][i+1] / h;
|
||||||
}
|
}
|
||||||
@@ -181,7 +181,7 @@ void Eigentql2 (sEigenvalue& eval) {
|
|||||||
// Compute implicit shift
|
// Compute implicit shift
|
||||||
|
|
||||||
float g = eval.d[l];
|
float g = eval.d[l];
|
||||||
float p = (eval.d[l+1] - g) / (2.0 * eval.e[l]);
|
float p = (eval.d[l+1] - g) / (2.0f * eval.e[l]);
|
||||||
float r = hypot(p,1.0);
|
float r = hypot(p,1.0);
|
||||||
if (p < 0) {
|
if (p < 0) {
|
||||||
r = -r;
|
r = -r;
|
||||||
@@ -280,7 +280,7 @@ void Eigenorthes (sEigenvalue& eval) {
|
|||||||
for (int i = m; i <= high; i++) {
|
for (int i = m; i <= high; i++) {
|
||||||
scale = scale + fabs(eval.H[i][m-1]);
|
scale = scale + fabs(eval.H[i][m-1]);
|
||||||
}
|
}
|
||||||
if (scale != 0.0) {
|
if (scale != 0.0f) {
|
||||||
|
|
||||||
// Compute Householder transformation.
|
// Compute Householder transformation.
|
||||||
|
|
||||||
@@ -334,7 +334,7 @@ void Eigenorthes (sEigenvalue& eval) {
|
|||||||
}
|
}
|
||||||
|
|
||||||
for (int m = high-1; m >= low+1; m--) {
|
for (int m = high-1; m >= low+1; m--) {
|
||||||
if (eval.H[m][m-1] != 0.0) {
|
if (eval.H[m][m-1] != 0.0f) {
|
||||||
for (int i = m+1; i <= high; i++) {
|
for (int i = m+1; i <= high; i++) {
|
||||||
eval.ort[i] = eval.H[i][m-1];
|
eval.ort[i] = eval.H[i][m-1];
|
||||||
}
|
}
|
||||||
@@ -406,7 +406,7 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
int l = n;
|
int l = n;
|
||||||
while (l > low) {
|
while (l > low) {
|
||||||
s = fabs(eval.H[l-1][l-1]) + fabs(eval.H[l][l]);
|
s = fabs(eval.H[l-1][l-1]) + fabs(eval.H[l][l]);
|
||||||
if (s == 0.0) {
|
if (s == 0.0f) {
|
||||||
s = norm;
|
s = norm;
|
||||||
}
|
}
|
||||||
if (fabs(eval.H[l][l-1]) < eps * s) {
|
if (fabs(eval.H[l][l-1]) < eps * s) {
|
||||||
@@ -429,7 +429,7 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
|
|
||||||
} else if (l == n-1) {
|
} else if (l == n-1) {
|
||||||
w = eval.H[n][n-1] * eval.H[n-1][n];
|
w = eval.H[n][n-1] * eval.H[n-1][n];
|
||||||
p = (eval.H[n-1][n-1] - eval.H[n][n]) / 2.0;
|
p = (eval.H[n-1][n-1] - eval.H[n][n]) / 2.0f;
|
||||||
q = p * p + w;
|
q = p * p + w;
|
||||||
z = sqrt(fabs(q));
|
z = sqrt(fabs(q));
|
||||||
eval.H[n][n] = eval.H[n][n] + exshift;
|
eval.H[n][n] = eval.H[n][n] + exshift;
|
||||||
@@ -446,7 +446,7 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
}
|
}
|
||||||
eval.d[n-1] = x + z;
|
eval.d[n-1] = x + z;
|
||||||
eval.d[n] = eval.d[n-1];
|
eval.d[n] = eval.d[n-1];
|
||||||
if (z != 0.0) {
|
if (z != 0.0f) {
|
||||||
eval.d[n] = x - w / z;
|
eval.d[n] = x - w / z;
|
||||||
}
|
}
|
||||||
eval.e[n-1] = 0.0;
|
eval.e[n-1] = 0.0;
|
||||||
@@ -516,21 +516,21 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
eval.H[i][i] -= x;
|
eval.H[i][i] -= x;
|
||||||
}
|
}
|
||||||
s = fabs(eval.H[n][n-1]) + fabs(eval.H[n-1][n-2]);
|
s = fabs(eval.H[n][n-1]) + fabs(eval.H[n-1][n-2]);
|
||||||
x = y = 0.75 * s;
|
x = y = 0.75f * s;
|
||||||
w = -0.4375 * s * s;
|
w = -0.4375f * s * s;
|
||||||
}
|
}
|
||||||
|
|
||||||
// MATLAB's new ad hoc shift
|
// MATLAB's new ad hoc shift
|
||||||
|
|
||||||
if (iter == 30) {
|
if (iter == 30) {
|
||||||
s = (y - x) / 2.0;
|
s = (y - x) / 2.0f;
|
||||||
s = s * s + w;
|
s = s * s + w;
|
||||||
if (s > 0) {
|
if (s > 0) {
|
||||||
s = sqrt(s);
|
s = sqrt(s);
|
||||||
if (y < x) {
|
if (y < x) {
|
||||||
s = -s;
|
s = -s;
|
||||||
}
|
}
|
||||||
s = x - w / ((y - x) / 2.0 + s);
|
s = x - w / ((y - x) / 2.0f + s);
|
||||||
for (int i = low; i <= n; i++) {
|
for (int i = low; i <= n; i++) {
|
||||||
eval.H[i][i] -= s;
|
eval.H[i][i] -= s;
|
||||||
}
|
}
|
||||||
@@ -580,15 +580,15 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
if (k != m) {
|
if (k != m) {
|
||||||
p = eval.H[k][k-1];
|
p = eval.H[k][k-1];
|
||||||
q = eval.H[k+1][k-1];
|
q = eval.H[k+1][k-1];
|
||||||
r = (notlast ? eval.H[k+2][k-1] : 0.0);
|
r = (notlast ? eval.H[k+2][k-1] : 0.0f);
|
||||||
x = fabs(p) + fabs(q) + fabs(r);
|
x = fabs(p) + fabs(q) + fabs(r);
|
||||||
if (x != 0.0) {
|
if (x != 0.0f) {
|
||||||
p = p / x;
|
p = p / x;
|
||||||
q = q / x;
|
q = q / x;
|
||||||
r = r / x;
|
r = r / x;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
if (x == 0.0) {
|
if (x == 0.0f) {
|
||||||
break;
|
break;
|
||||||
}
|
}
|
||||||
s = sqrt(p * p + q * q + r * r);
|
s = sqrt(p * p + q * q + r * r);
|
||||||
@@ -651,7 +651,7 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
|
|
||||||
// Backsubstitute to find vectors of upper triangular form
|
// Backsubstitute to find vectors of upper triangular form
|
||||||
|
|
||||||
if (norm == 0.0) {
|
if (norm == 0.0f) {
|
||||||
return;
|
return;
|
||||||
}
|
}
|
||||||
|
|
||||||
@@ -670,13 +670,13 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
for (int j = l; j <= n; j++) {
|
for (int j = l; j <= n; j++) {
|
||||||
r = r + eval.H[i][j] * eval.H[j][n];
|
r = r + eval.H[i][j] * eval.H[j][n];
|
||||||
}
|
}
|
||||||
if (eval.e[i] < 0.0) {
|
if (eval.e[i] < 0.0f) {
|
||||||
z = w;
|
z = w;
|
||||||
s = r;
|
s = r;
|
||||||
} else {
|
} else {
|
||||||
l = i;
|
l = i;
|
||||||
if (eval.e[i] == 0.0) {
|
if (eval.e[i] == 0.0f) {
|
||||||
if (w != 0.0) {
|
if (w != 0.0f) {
|
||||||
eval.H[i][n] = -r / w;
|
eval.H[i][n] = -r / w;
|
||||||
} else {
|
} else {
|
||||||
eval.H[i][n] = -r / (eps * norm);
|
eval.H[i][n] = -r / (eps * norm);
|
||||||
@@ -735,7 +735,7 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
}
|
}
|
||||||
w = eval.H[i][i] - p;
|
w = eval.H[i][i] - p;
|
||||||
|
|
||||||
if (eval.e[i] < 0.0) {
|
if (eval.e[i] < 0.0f) {
|
||||||
z = w;
|
z = w;
|
||||||
r = ra;
|
r = ra;
|
||||||
s = sa;
|
s = sa;
|
||||||
@@ -752,8 +752,8 @@ void Eigenhqr2 (sEigenvalue& eval) {
|
|||||||
x = eval.H[i][i+1];
|
x = eval.H[i][i+1];
|
||||||
y = eval.H[i+1][i];
|
y = eval.H[i+1][i];
|
||||||
vr = (eval.d[i] - p) * (eval.d[i] - p) + eval.e[i] * eval.e[i] - q * q;
|
vr = (eval.d[i] - p) * (eval.d[i] - p) + eval.e[i] * eval.e[i] - q * q;
|
||||||
vi = (eval.d[i] - p) * 2.0 * q;
|
vi = (eval.d[i] - p) * 2.0f * q;
|
||||||
if ((vr == 0.0) && (vi == 0.0)) {
|
if ((vr == 0.0f) && (vi == 0.0f)) {
|
||||||
vr = eps * norm * (fabs(w) + fabs(q) +
|
vr = eps * norm * (fabs(w) + fabs(q) +
|
||||||
fabs(x) + fabs(y) + fabs(z));
|
fabs(x) + fabs(y) + fabs(z));
|
||||||
}
|
}
|
||||||
|
|||||||
@@ -574,17 +574,17 @@ void FLUID_3D::artificialDampingSL(int zBegin, int zEnd) {
|
|||||||
for (int y = 1; y < _res[1]-1; y++)
|
for (int y = 1; y < _res[1]-1; y++)
|
||||||
for (int x = 1+(y+z)%2; x < _res[0]-1; x+=2) {
|
for (int x = 1+(y+z)%2; x < _res[0]-1; x+=2) {
|
||||||
const int index = x + y*_res[0] + z * _slabSize;
|
const int index = x + y*_res[0] + z * _slabSize;
|
||||||
_xForce[index] = (1-w)*_xVelocityTemp[index] + 1./6. * w*(
|
_xForce[index] = (1-w)*_xVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_xVelocityTemp[index+1] + _xVelocityTemp[index-1] +
|
_xVelocityTemp[index+1] + _xVelocityTemp[index-1] +
|
||||||
_xVelocityTemp[index+_res[0]] + _xVelocityTemp[index-_res[0]] +
|
_xVelocityTemp[index+_res[0]] + _xVelocityTemp[index-_res[0]] +
|
||||||
_xVelocityTemp[index+_slabSize] + _xVelocityTemp[index-_slabSize] );
|
_xVelocityTemp[index+_slabSize] + _xVelocityTemp[index-_slabSize] );
|
||||||
|
|
||||||
_yForce[index] = (1-w)*_yVelocityTemp[index] + 1./6. * w*(
|
_yForce[index] = (1-w)*_yVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_yVelocityTemp[index+1] + _yVelocityTemp[index-1] +
|
_yVelocityTemp[index+1] + _yVelocityTemp[index-1] +
|
||||||
_yVelocityTemp[index+_res[0]] + _yVelocityTemp[index-_res[0]] +
|
_yVelocityTemp[index+_res[0]] + _yVelocityTemp[index-_res[0]] +
|
||||||
_yVelocityTemp[index+_slabSize] + _yVelocityTemp[index-_slabSize] );
|
_yVelocityTemp[index+_slabSize] + _yVelocityTemp[index-_slabSize] );
|
||||||
|
|
||||||
_zForce[index] = (1-w)*_zVelocityTemp[index] + 1./6. * w*(
|
_zForce[index] = (1-w)*_zVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_zVelocityTemp[index+1] + _zVelocityTemp[index-1] +
|
_zVelocityTemp[index+1] + _zVelocityTemp[index-1] +
|
||||||
_zVelocityTemp[index+_res[0]] + _zVelocityTemp[index-_res[0]] +
|
_zVelocityTemp[index+_res[0]] + _zVelocityTemp[index-_res[0]] +
|
||||||
_zVelocityTemp[index+_slabSize] + _zVelocityTemp[index-_slabSize] );
|
_zVelocityTemp[index+_slabSize] + _zVelocityTemp[index-_slabSize] );
|
||||||
@@ -596,17 +596,17 @@ void FLUID_3D::artificialDampingSL(int zBegin, int zEnd) {
|
|||||||
for (int y = 1; y < _res[1]-1; y++)
|
for (int y = 1; y < _res[1]-1; y++)
|
||||||
for (int x = 1+(y+z+1)%2; x < _res[0]-1; x+=2) {
|
for (int x = 1+(y+z+1)%2; x < _res[0]-1; x+=2) {
|
||||||
const int index = x + y*_res[0] + z * _slabSize;
|
const int index = x + y*_res[0] + z * _slabSize;
|
||||||
_xForce[index] = (1-w)*_xVelocityTemp[index] + 1./6. * w*(
|
_xForce[index] = (1-w)*_xVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_xVelocityTemp[index+1] + _xVelocityTemp[index-1] +
|
_xVelocityTemp[index+1] + _xVelocityTemp[index-1] +
|
||||||
_xVelocityTemp[index+_res[0]] + _xVelocityTemp[index-_res[0]] +
|
_xVelocityTemp[index+_res[0]] + _xVelocityTemp[index-_res[0]] +
|
||||||
_xVelocityTemp[index+_slabSize] + _xVelocityTemp[index-_slabSize] );
|
_xVelocityTemp[index+_slabSize] + _xVelocityTemp[index-_slabSize] );
|
||||||
|
|
||||||
_yForce[index] = (1-w)*_yVelocityTemp[index] + 1./6. * w*(
|
_yForce[index] = (1-w)*_yVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_yVelocityTemp[index+1] + _yVelocityTemp[index-1] +
|
_yVelocityTemp[index+1] + _yVelocityTemp[index-1] +
|
||||||
_yVelocityTemp[index+_res[0]] + _yVelocityTemp[index-_res[0]] +
|
_yVelocityTemp[index+_res[0]] + _yVelocityTemp[index-_res[0]] +
|
||||||
_yVelocityTemp[index+_slabSize] + _yVelocityTemp[index-_slabSize] );
|
_yVelocityTemp[index+_slabSize] + _yVelocityTemp[index-_slabSize] );
|
||||||
|
|
||||||
_zForce[index] = (1-w)*_zVelocityTemp[index] + 1./6. * w*(
|
_zForce[index] = (1-w)*_zVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_zVelocityTemp[index+1] + _zVelocityTemp[index-1] +
|
_zVelocityTemp[index+1] + _zVelocityTemp[index-1] +
|
||||||
_zVelocityTemp[index+_res[0]] + _zVelocityTemp[index-_res[0]] +
|
_zVelocityTemp[index+_res[0]] + _zVelocityTemp[index-_res[0]] +
|
||||||
_zVelocityTemp[index+_slabSize] + _zVelocityTemp[index-_slabSize] );
|
_zVelocityTemp[index+_slabSize] + _zVelocityTemp[index-_slabSize] );
|
||||||
@@ -636,17 +636,17 @@ void FLUID_3D::artificialDampingExactSL(int pos) {
|
|||||||
* Uses xForce as temporary storage to allow other threads to read
|
* Uses xForce as temporary storage to allow other threads to read
|
||||||
* old values from xVelocityTemp
|
* old values from xVelocityTemp
|
||||||
*/
|
*/
|
||||||
_xForce[index] = (1-w)*_xVelocityTemp[index] + 1./6. * w*(
|
_xForce[index] = (1-w)*_xVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_xVelocityTemp[index+1] + _xVelocityTemp[index-1] +
|
_xVelocityTemp[index+1] + _xVelocityTemp[index-1] +
|
||||||
_xVelocityTemp[index+_res[0]] + _xVelocityTemp[index-_res[0]] +
|
_xVelocityTemp[index+_res[0]] + _xVelocityTemp[index-_res[0]] +
|
||||||
_xVelocityTemp[index+_slabSize] + _xVelocityTemp[index-_slabSize] );
|
_xVelocityTemp[index+_slabSize] + _xVelocityTemp[index-_slabSize] );
|
||||||
|
|
||||||
_yForce[index] = (1-w)*_yVelocityTemp[index] + 1./6. * w*(
|
_yForce[index] = (1-w)*_yVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_yVelocityTemp[index+1] + _yVelocityTemp[index-1] +
|
_yVelocityTemp[index+1] + _yVelocityTemp[index-1] +
|
||||||
_yVelocityTemp[index+_res[0]] + _yVelocityTemp[index-_res[0]] +
|
_yVelocityTemp[index+_res[0]] + _yVelocityTemp[index-_res[0]] +
|
||||||
_yVelocityTemp[index+_slabSize] + _yVelocityTemp[index-_slabSize] );
|
_yVelocityTemp[index+_slabSize] + _yVelocityTemp[index-_slabSize] );
|
||||||
|
|
||||||
_zForce[index] = (1-w)*_zVelocityTemp[index] + 1./6. * w*(
|
_zForce[index] = (1-w)*_zVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_zVelocityTemp[index+1] + _zVelocityTemp[index-1] +
|
_zVelocityTemp[index+1] + _zVelocityTemp[index-1] +
|
||||||
_zVelocityTemp[index+_res[0]] + _zVelocityTemp[index-_res[0]] +
|
_zVelocityTemp[index+_res[0]] + _zVelocityTemp[index-_res[0]] +
|
||||||
_zVelocityTemp[index+_slabSize] + _zVelocityTemp[index-_slabSize] );
|
_zVelocityTemp[index+_slabSize] + _zVelocityTemp[index-_slabSize] );
|
||||||
@@ -663,17 +663,17 @@ void FLUID_3D::artificialDampingExactSL(int pos) {
|
|||||||
* Uses xForce as temporary storage to allow other threads to read
|
* Uses xForce as temporary storage to allow other threads to read
|
||||||
* old values from xVelocityTemp
|
* old values from xVelocityTemp
|
||||||
*/
|
*/
|
||||||
_xForce[index] = (1-w)*_xVelocityTemp[index] + 1./6. * w*(
|
_xForce[index] = (1-w)*_xVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_xVelocityTemp[index+1] + _xVelocityTemp[index-1] +
|
_xVelocityTemp[index+1] + _xVelocityTemp[index-1] +
|
||||||
_xVelocityTemp[index+_res[0]] + _xVelocityTemp[index-_res[0]] +
|
_xVelocityTemp[index+_res[0]] + _xVelocityTemp[index-_res[0]] +
|
||||||
_xVelocityTemp[index+_slabSize] + _xVelocityTemp[index-_slabSize] );
|
_xVelocityTemp[index+_slabSize] + _xVelocityTemp[index-_slabSize] );
|
||||||
|
|
||||||
_yForce[index] = (1-w)*_yVelocityTemp[index] + 1./6. * w*(
|
_yForce[index] = (1-w)*_yVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_yVelocityTemp[index+1] + _yVelocityTemp[index-1] +
|
_yVelocityTemp[index+1] + _yVelocityTemp[index-1] +
|
||||||
_yVelocityTemp[index+_res[0]] + _yVelocityTemp[index-_res[0]] +
|
_yVelocityTemp[index+_res[0]] + _yVelocityTemp[index-_res[0]] +
|
||||||
_yVelocityTemp[index+_slabSize] + _yVelocityTemp[index-_slabSize] );
|
_yVelocityTemp[index+_slabSize] + _yVelocityTemp[index-_slabSize] );
|
||||||
|
|
||||||
_zForce[index] = (1-w)*_zVelocityTemp[index] + 1./6. * w*(
|
_zForce[index] = (1-w)*_zVelocityTemp[index] + 1.0f/6.0f * w*(
|
||||||
_zVelocityTemp[index+1] + _zVelocityTemp[index-1] +
|
_zVelocityTemp[index+1] + _zVelocityTemp[index-1] +
|
||||||
_zVelocityTemp[index+_res[0]] + _zVelocityTemp[index-_res[0]] +
|
_zVelocityTemp[index+_res[0]] + _zVelocityTemp[index-_res[0]] +
|
||||||
_zVelocityTemp[index+_slabSize] + _zVelocityTemp[index-_slabSize] );
|
_zVelocityTemp[index+_slabSize] + _zVelocityTemp[index-_slabSize] );
|
||||||
@@ -1244,27 +1244,27 @@ void FLUID_3D::setObstaclePressure(float *_pressure, int zBegin, int zEnd)
|
|||||||
float pcnt = 0.;
|
float pcnt = 0.;
|
||||||
if (left && !right) {
|
if (left && !right) {
|
||||||
_pressure[index] += _pressure[index + 1];
|
_pressure[index] += _pressure[index + 1];
|
||||||
pcnt += 1.;
|
pcnt += 1.0f;
|
||||||
}
|
}
|
||||||
if (!left && right) {
|
if (!left && right) {
|
||||||
_pressure[index] += _pressure[index - 1];
|
_pressure[index] += _pressure[index - 1];
|
||||||
pcnt += 1.;
|
pcnt += 1.0f;
|
||||||
}
|
}
|
||||||
if (up && !down) {
|
if (up && !down) {
|
||||||
_pressure[index] += _pressure[index - _xRes];
|
_pressure[index] += _pressure[index - _xRes];
|
||||||
pcnt += 1.;
|
pcnt += 1.0f;
|
||||||
}
|
}
|
||||||
if (!up && down) {
|
if (!up && down) {
|
||||||
_pressure[index] += _pressure[index + _xRes];
|
_pressure[index] += _pressure[index + _xRes];
|
||||||
pcnt += 1.;
|
pcnt += 1.0f;
|
||||||
}
|
}
|
||||||
if (top && !bottom) {
|
if (top && !bottom) {
|
||||||
_pressure[index] += _pressure[index - _slabSize];
|
_pressure[index] += _pressure[index - _slabSize];
|
||||||
pcnt += 1.;
|
pcnt += 1.0f;
|
||||||
}
|
}
|
||||||
if (!top && bottom) {
|
if (!top && bottom) {
|
||||||
_pressure[index] += _pressure[index + _slabSize];
|
_pressure[index] += _pressure[index + _slabSize];
|
||||||
pcnt += 1.;
|
pcnt += 1.0f;
|
||||||
}
|
}
|
||||||
|
|
||||||
if(pcnt > 0.000001f)
|
if(pcnt > 0.000001f)
|
||||||
@@ -1369,7 +1369,7 @@ void FLUID_3D::addVorticity(int zBegin, int zEnd)
|
|||||||
// set flame vorticity from RNA value
|
// set flame vorticity from RNA value
|
||||||
float flame_vorticity = (*_flame_vorticity)/_constantScaling;
|
float flame_vorticity = (*_flame_vorticity)/_constantScaling;
|
||||||
//int x,y,z,index;
|
//int x,y,z,index;
|
||||||
if(_vorticityEps+flame_vorticity<=0.) return;
|
if(_vorticityEps+flame_vorticity<=0.0f) return;
|
||||||
|
|
||||||
int _blockSize=zEnd-zBegin;
|
int _blockSize=zEnd-zBegin;
|
||||||
int _blockTotalCells = _slabSize * (_blockSize+2);
|
int _blockTotalCells = _slabSize * (_blockSize+2);
|
||||||
|
|||||||
@@ -76,12 +76,12 @@ void FLUID_3D::solveHeat(float* field, float* b, unsigned char* skip)
|
|||||||
if (!skip[index - _slabSize]) _Acenter[index] += heatConst;
|
if (!skip[index - _slabSize]) _Acenter[index] += heatConst;
|
||||||
|
|
||||||
_residual[index] = b[index] - (_Acenter[index] * field[index] +
|
_residual[index] = b[index] - (_Acenter[index] * field[index] +
|
||||||
field[index - 1] * (skip[index - 1] ? 0.0 : -heatConst) +
|
field[index - 1] * (skip[index - 1] ? 0.0f : -heatConst) +
|
||||||
field[index + 1] * (skip[index + 1] ? 0.0 : -heatConst) +
|
field[index + 1] * (skip[index + 1] ? 0.0f : -heatConst) +
|
||||||
field[index - _xRes] * (skip[index - _xRes] ? 0.0 : -heatConst) +
|
field[index - _xRes] * (skip[index - _xRes] ? 0.0f : -heatConst) +
|
||||||
field[index + _xRes] * (skip[index + _xRes] ? 0.0 : -heatConst) +
|
field[index + _xRes] * (skip[index + _xRes] ? 0.0f : -heatConst) +
|
||||||
field[index - _slabSize] * (skip[index - _slabSize] ? 0.0 : -heatConst) +
|
field[index - _slabSize] * (skip[index - _slabSize] ? 0.0f : -heatConst) +
|
||||||
field[index + _slabSize] * (skip[index + _slabSize] ? 0.0 : -heatConst));
|
field[index + _slabSize] * (skip[index + _slabSize] ? 0.0f : -heatConst));
|
||||||
}
|
}
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
@@ -111,12 +111,12 @@ void FLUID_3D::solveHeat(float* field, float* b, unsigned char* skip)
|
|||||||
{
|
{
|
||||||
|
|
||||||
_q[index] = (_Acenter[index] * _direction[index] +
|
_q[index] = (_Acenter[index] * _direction[index] +
|
||||||
_direction[index - 1] * (skip[index - 1] ? 0.0 : -heatConst) +
|
_direction[index - 1] * (skip[index - 1] ? 0.0f : -heatConst) +
|
||||||
_direction[index + 1] * (skip[index + 1] ? 0.0 : -heatConst) +
|
_direction[index + 1] * (skip[index + 1] ? 0.0f : -heatConst) +
|
||||||
_direction[index - _xRes] * (skip[index - _xRes] ? 0.0 : -heatConst) +
|
_direction[index - _xRes] * (skip[index - _xRes] ? 0.0f : -heatConst) +
|
||||||
_direction[index + _xRes] * (skip[index + _xRes] ? 0.0 : -heatConst) +
|
_direction[index + _xRes] * (skip[index + _xRes] ? 0.0f : -heatConst) +
|
||||||
_direction[index - _slabSize] * (skip[index - _slabSize] ? 0.0 : -heatConst) +
|
_direction[index - _slabSize] * (skip[index - _slabSize] ? 0.0f : -heatConst) +
|
||||||
_direction[index + _slabSize] * (skip[index + _slabSize] ? 0.0 : -heatConst));
|
_direction[index + _slabSize] * (skip[index + _slabSize] ? 0.0f : -heatConst));
|
||||||
}
|
}
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
@@ -199,20 +199,20 @@ void FLUID_3D::solvePressurePre(float* field, float* b, unsigned char* skip)
|
|||||||
if (!skip[index])
|
if (!skip[index])
|
||||||
{
|
{
|
||||||
// set the matrix to the Poisson stencil in order
|
// set the matrix to the Poisson stencil in order
|
||||||
if (!skip[index + 1]) Acenter += 1.;
|
if (!skip[index + 1]) Acenter += 1.0f;
|
||||||
if (!skip[index - 1]) Acenter += 1.;
|
if (!skip[index - 1]) Acenter += 1.0f;
|
||||||
if (!skip[index + _xRes]) Acenter += 1.;
|
if (!skip[index + _xRes]) Acenter += 1.0f;
|
||||||
if (!skip[index - _xRes]) Acenter += 1.;
|
if (!skip[index - _xRes]) Acenter += 1.0f;
|
||||||
if (!skip[index + _slabSize]) Acenter += 1.;
|
if (!skip[index + _slabSize]) Acenter += 1.0f;
|
||||||
if (!skip[index - _slabSize]) Acenter += 1.;
|
if (!skip[index - _slabSize]) Acenter += 1.0f;
|
||||||
|
|
||||||
_residual[index] = b[index] - (Acenter * field[index] +
|
_residual[index] = b[index] - (Acenter * field[index] +
|
||||||
field[index - 1] * (skip[index - 1] ? 0.0 : -1.0f)+
|
field[index - 1] * (skip[index - 1] ? 0.0f : -1.0f) +
|
||||||
field[index + 1] * (skip[index + 1] ? 0.0 : -1.0f)+
|
field[index + 1] * (skip[index + 1] ? 0.0f : -1.0f) +
|
||||||
field[index - _xRes] * (skip[index - _xRes] ? 0.0 : -1.0f)+
|
field[index - _xRes] * (skip[index - _xRes] ? 0.0f : -1.0f)+
|
||||||
field[index + _xRes] * (skip[index + _xRes] ? 0.0 : -1.0f)+
|
field[index + _xRes] * (skip[index + _xRes] ? 0.0f : -1.0f)+
|
||||||
field[index - _slabSize] * (skip[index - _slabSize] ? 0.0 : -1.0f)+
|
field[index - _slabSize] * (skip[index - _slabSize] ? 0.0f : -1.0f)+
|
||||||
field[index + _slabSize] * (skip[index + _slabSize] ? 0.0 : -1.0f) );
|
field[index + _slabSize] * (skip[index + _slabSize] ? 0.0f : -1.0f) );
|
||||||
}
|
}
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
@@ -220,10 +220,10 @@ void FLUID_3D::solvePressurePre(float* field, float* b, unsigned char* skip)
|
|||||||
}
|
}
|
||||||
|
|
||||||
// P^-1
|
// P^-1
|
||||||
if(Acenter < 1.0)
|
if(Acenter < 1.0f)
|
||||||
_Precond[index] = 0.0;
|
_Precond[index] = 0.0;
|
||||||
else
|
else
|
||||||
_Precond[index] = 1.0 / Acenter;
|
_Precond[index] = 1.0f / Acenter;
|
||||||
|
|
||||||
// p = P^-1 * r
|
// p = P^-1 * r
|
||||||
_direction[index] = _residual[index] * _Precond[index];
|
_direction[index] = _residual[index] * _Precond[index];
|
||||||
@@ -237,7 +237,7 @@ void FLUID_3D::solvePressurePre(float* field, float* b, unsigned char* skip)
|
|||||||
//while ((i < _iterations) && (deltaNew > eps*delta0))
|
//while ((i < _iterations) && (deltaNew > eps*delta0))
|
||||||
float maxR = 2.0f * eps;
|
float maxR = 2.0f * eps;
|
||||||
// while (i < _iterations)
|
// while (i < _iterations)
|
||||||
while ((i < _iterations) && (maxR > 0.001*eps))
|
while ((i < _iterations) && (maxR > 0.001f * eps))
|
||||||
{
|
{
|
||||||
|
|
||||||
float alpha = 0.0f;
|
float alpha = 0.0f;
|
||||||
@@ -252,20 +252,20 @@ void FLUID_3D::solvePressurePre(float* field, float* b, unsigned char* skip)
|
|||||||
if (!skip[index])
|
if (!skip[index])
|
||||||
{
|
{
|
||||||
// set the matrix to the Poisson stencil in order
|
// set the matrix to the Poisson stencil in order
|
||||||
if (!skip[index + 1]) Acenter += 1.;
|
if (!skip[index + 1]) Acenter += 1.0f;
|
||||||
if (!skip[index - 1]) Acenter += 1.;
|
if (!skip[index - 1]) Acenter += 1.0f;
|
||||||
if (!skip[index + _xRes]) Acenter += 1.;
|
if (!skip[index + _xRes]) Acenter += 1.0f;
|
||||||
if (!skip[index - _xRes]) Acenter += 1.;
|
if (!skip[index - _xRes]) Acenter += 1.0f;
|
||||||
if (!skip[index + _slabSize]) Acenter += 1.;
|
if (!skip[index + _slabSize]) Acenter += 1.0f;
|
||||||
if (!skip[index - _slabSize]) Acenter += 1.;
|
if (!skip[index - _slabSize]) Acenter += 1.0f;
|
||||||
|
|
||||||
_q[index] = Acenter * _direction[index] +
|
_q[index] = Acenter * _direction[index] +
|
||||||
_direction[index - 1] * (skip[index - 1] ? 0.0 : -1.0f) +
|
_direction[index - 1] * (skip[index - 1] ? 0.0f : -1.0f) +
|
||||||
_direction[index + 1] * (skip[index + 1] ? 0.0 : -1.0f) +
|
_direction[index + 1] * (skip[index + 1] ? 0.0f : -1.0f) +
|
||||||
_direction[index - _xRes] * (skip[index - _xRes] ? 0.0 : -1.0f) +
|
_direction[index - _xRes] * (skip[index - _xRes] ? 0.0f : -1.0f) +
|
||||||
_direction[index + _xRes] * (skip[index + _xRes] ? 0.0 : -1.0f)+
|
_direction[index + _xRes] * (skip[index + _xRes] ? 0.0f : -1.0f)+
|
||||||
_direction[index - _slabSize] * (skip[index - _slabSize] ? 0.0 : -1.0f) +
|
_direction[index - _slabSize] * (skip[index - _slabSize] ? 0.0f : -1.0f) +
|
||||||
_direction[index + _slabSize] * (skip[index + _slabSize] ? 0.0 : -1.0f);
|
_direction[index + _slabSize] * (skip[index + _slabSize] ? 0.0f : -1.0f);
|
||||||
}
|
}
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
|
|||||||
@@ -317,12 +317,12 @@ void FLUID_3D::advectFieldSemiLagrange(const float dt, const float* velx, const
|
|||||||
float zTrace = z - dt * velz[index];
|
float zTrace = z - dt * velz[index];
|
||||||
|
|
||||||
// clamp backtrace to grid boundaries
|
// clamp backtrace to grid boundaries
|
||||||
if (xTrace < 0.5) xTrace = 0.5;
|
if (xTrace < 0.5f) xTrace = 0.5f;
|
||||||
if (xTrace > xres - 1.5) xTrace = xres - 1.5;
|
if (xTrace > xres - 1.5f) xTrace = xres - 1.5f;
|
||||||
if (yTrace < 0.5) yTrace = 0.5;
|
if (yTrace < 0.5f) yTrace = 0.5f;
|
||||||
if (yTrace > yres - 1.5) yTrace = yres - 1.5;
|
if (yTrace > yres - 1.5f) yTrace = yres - 1.5f;
|
||||||
if (zTrace < 0.5) zTrace = 0.5;
|
if (zTrace < 0.5f) zTrace = 0.5f;
|
||||||
if (zTrace > zres - 1.5) zTrace = zres - 1.5;
|
if (zTrace > zres - 1.5f) zTrace = zres - 1.5f;
|
||||||
|
|
||||||
// locate neighbors to interpolate
|
// locate neighbors to interpolate
|
||||||
const int x0 = (int)xTrace;
|
const int x0 = (int)xTrace;
|
||||||
@@ -403,7 +403,7 @@ void FLUID_3D::advectFieldMacCormack2(const float dt, const float* xVelocity, co
|
|||||||
|
|
||||||
|
|
||||||
// phiHatN = A^R(phiHatN1)
|
// phiHatN = A^R(phiHatN1)
|
||||||
advectFieldSemiLagrange( -1.0*dt, xVelocity, yVelocity, zVelocity, phiHatN, t1, res, zBegin, zEnd); // uses wide data from old field and velocities (both are whole)
|
advectFieldSemiLagrange( -1.0f*dt, xVelocity, yVelocity, zVelocity, phiHatN, t1, res, zBegin, zEnd); // uses wide data from old field and velocities (both are whole)
|
||||||
|
|
||||||
// phiN1 = phiHatN1 + (phiN - phiHatN) / 2
|
// phiN1 = phiHatN1 + (phiN - phiHatN) / 2
|
||||||
const int border = 0;
|
const int border = 0;
|
||||||
@@ -456,12 +456,12 @@ void FLUID_3D::clampExtrema(const float dt, const float* velx, const float* vely
|
|||||||
float zTrace = z - dt * velz[index];
|
float zTrace = z - dt * velz[index];
|
||||||
|
|
||||||
// clamp backtrace to grid boundaries
|
// clamp backtrace to grid boundaries
|
||||||
if (xTrace < 0.5) xTrace = 0.5;
|
if (xTrace < 0.5f) xTrace = 0.5f;
|
||||||
if (xTrace > xres - 1.5) xTrace = xres - 1.5;
|
if (xTrace > xres - 1.5f) xTrace = xres - 1.5f;
|
||||||
if (yTrace < 0.5) yTrace = 0.5;
|
if (yTrace < 0.5f) yTrace = 0.5f;
|
||||||
if (yTrace > yres - 1.5) yTrace = yres - 1.5;
|
if (yTrace > yres - 1.5f) yTrace = yres - 1.5f;
|
||||||
if (zTrace < 0.5) zTrace = 0.5;
|
if (zTrace < 0.5f) zTrace = 0.5f;
|
||||||
if (zTrace > zres - 1.5) zTrace = zres - 1.5;
|
if (zTrace > zres - 1.5f) zTrace = zres - 1.5f;
|
||||||
|
|
||||||
// locate neighbors to interpolate
|
// locate neighbors to interpolate
|
||||||
const int x0 = (int)xTrace;
|
const int x0 = (int)xTrace;
|
||||||
|
|||||||
@@ -62,12 +62,12 @@ static inline float lerp(float* field, float x, float y, int res) {
|
|||||||
//////////////////////////////////////////////////////////////////////////////////////////
|
//////////////////////////////////////////////////////////////////////////////////////////
|
||||||
static inline float lerp3d(float* field, float x, float y, float z, int xres, int yres, int zres) {
|
static inline float lerp3d(float* field, float x, float y, float z, int xres, int yres, int zres) {
|
||||||
// clamp pos to grid boundaries
|
// clamp pos to grid boundaries
|
||||||
if (x < 0.5) x = 0.5;
|
if (x < 0.5f) x = 0.5f;
|
||||||
if (x > xres - 1.5) x = xres - 1.5;
|
if (x > xres - 1.5f) x = xres - 1.5f;
|
||||||
if (y < 0.5) y = 0.5;
|
if (y < 0.5f) y = 0.5f;
|
||||||
if (y > yres - 1.5) y = yres - 1.5;
|
if (y > yres - 1.5f) y = yres - 1.5f;
|
||||||
if (z < 0.5) z = 0.5;
|
if (z < 0.5f) z = 0.5f;
|
||||||
if (z > zres - 1.5) z = zres - 1.5;
|
if (z > zres - 1.5f) z = zres - 1.5f;
|
||||||
|
|
||||||
// locate neighbors to interpolate
|
// locate neighbors to interpolate
|
||||||
const int x0 = (int)x;
|
const int x0 = (int)x;
|
||||||
@@ -113,12 +113,12 @@ template <class T>
|
|||||||
static inline float lerp3dToFloat(T* field1,
|
static inline float lerp3dToFloat(T* field1,
|
||||||
float x, float y, float z, int xres, int yres, int zres) {
|
float x, float y, float z, int xres, int yres, int zres) {
|
||||||
// clamp pos to grid boundaries
|
// clamp pos to grid boundaries
|
||||||
if (x < 0.5) x = 0.5;
|
if (x < 0.5f) x = 0.5f;
|
||||||
if (x > xres - 1.5) x = xres - 1.5;
|
if (x > xres - 1.5f) x = xres - 1.5f;
|
||||||
if (y < 0.5) y = 0.5;
|
if (y < 0.5f) y = 0.5f;
|
||||||
if (y > yres - 1.5) y = yres - 1.5;
|
if (y > yres - 1.5f) y = yres - 1.5f;
|
||||||
if (z < 0.5) z = 0.5;
|
if (z < 0.5f) z = 0.5f;
|
||||||
if (z > zres - 1.5) z = zres - 1.5;
|
if (z > zres - 1.5f) z = zres - 1.5f;
|
||||||
|
|
||||||
// locate neighbors to interpolate
|
// locate neighbors to interpolate
|
||||||
const int x0 = (int)x;
|
const int x0 = (int)x;
|
||||||
@@ -164,12 +164,12 @@ static inline float lerp3dToFloat(T* field1,
|
|||||||
static inline Vec3 lerp3dVec(float* field1, float* field2, float* field3,
|
static inline Vec3 lerp3dVec(float* field1, float* field2, float* field3,
|
||||||
float x, float y, float z, int xres, int yres, int zres) {
|
float x, float y, float z, int xres, int yres, int zres) {
|
||||||
// clamp pos to grid boundaries
|
// clamp pos to grid boundaries
|
||||||
if (x < 0.5) x = 0.5;
|
if (x < 0.5f) x = 0.5f;
|
||||||
if (x > xres - 1.5) x = xres - 1.5;
|
if (x > xres - 1.5f) x = xres - 1.5f;
|
||||||
if (y < 0.5) y = 0.5;
|
if (y < 0.5f) y = 0.5f;
|
||||||
if (y > yres - 1.5) y = yres - 1.5;
|
if (y > yres - 1.5f) y = yres - 1.5f;
|
||||||
if (z < 0.5) z = 0.5;
|
if (z < 0.5f) z = 0.5f;
|
||||||
if (z > zres - 1.5) z = zres - 1.5;
|
if (z > zres - 1.5f) z = zres - 1.5f;
|
||||||
|
|
||||||
// locate neighbors to interpolate
|
// locate neighbors to interpolate
|
||||||
const int x0 = (int)x;
|
const int x0 = (int)x;
|
||||||
|
|||||||
@@ -51,10 +51,10 @@ sLU computeLU( float a[3][3])
|
|||||||
int kmax = min(i,j);
|
int kmax = min(i,j);
|
||||||
double s = 0.0;
|
double s = 0.0;
|
||||||
for (int k = 0; k < kmax; k++) {
|
for (int k = 0; k < kmax; k++) {
|
||||||
s += result.values[i][k]*LUcolj[k];
|
s += (double)(result.values[i][k]*LUcolj[k]);
|
||||||
}
|
}
|
||||||
|
|
||||||
result.values[i][j] = LUcolj[i] -= s;
|
result.values[i][j] = LUcolj[i] -= (float)s;
|
||||||
}
|
}
|
||||||
|
|
||||||
// Find pivot and exchange if necessary.
|
// Find pivot and exchange if necessary.
|
||||||
@@ -80,7 +80,7 @@ sLU computeLU( float a[3][3])
|
|||||||
|
|
||||||
// Compute multipliers.
|
// Compute multipliers.
|
||||||
|
|
||||||
if ((j < m) && (result.values[j][j] != 0.0)) {
|
if ((j < m) && (result.values[j][j] != 0.0f)) {
|
||||||
for (int i = j+1; i < m; i++) {
|
for (int i = j+1; i < m; i++) {
|
||||||
result.values[i][j] /= result.values[j][j];
|
result.values[i][j] /= result.values[j][j];
|
||||||
}
|
}
|
||||||
|
|||||||
@@ -730,7 +730,7 @@ inline Real normHelper(const Vector3Dim<Real> &v) {
|
|||||||
return norm(v);
|
return norm(v);
|
||||||
}
|
}
|
||||||
inline Real normHelper(const Real &v) {
|
inline Real normHelper(const Real &v) {
|
||||||
return (0. < v) ? v : -v ;
|
return (0.0f < v) ? v : -v ;
|
||||||
}
|
}
|
||||||
inline Real normHelper(const int &v) {
|
inline Real normHelper(const int &v) {
|
||||||
return (0 < v) ? (Real)(v) : (Real)(-v) ;
|
return (0 < v) ? (Real)(v) : (Real)(-v) ;
|
||||||
|
|||||||
@@ -395,20 +395,20 @@ static inline float WNoiseDx(Vec3 p, float* data) {
|
|||||||
int c[3], mid[3], n = noiseTileSize;
|
int c[3], mid[3], n = noiseTileSize;
|
||||||
float w[3][3], t, result = 0;
|
float w[3][3], t, result = 0;
|
||||||
|
|
||||||
mid[0] = (int)ceil(p[0] - 0.5);
|
mid[0] = (int)ceil(p[0] - 0.5f);
|
||||||
t = mid[0] - (p[0] - 0.5);
|
t = mid[0] - (p[0] - 0.5f);
|
||||||
w[0][0] = -t;
|
w[0][0] = -t;
|
||||||
w[0][2] = (1.f - t);
|
w[0][2] = (1.f - t);
|
||||||
w[0][1] = 2.0f * t - 1.0f;
|
w[0][1] = 2.0f * t - 1.0f;
|
||||||
|
|
||||||
mid[1] = (int)ceil(p[1] - 0.5);
|
mid[1] = (int)ceil(p[1] - 0.5f);
|
||||||
t = mid[1] - (p[1] - 0.5);
|
t = mid[1] - (p[1] - 0.5f);
|
||||||
w[1][0] = t * t / 2;
|
w[1][0] = t * t / 2;
|
||||||
w[1][2] = (1 - t) * (1 - t) / 2;
|
w[1][2] = (1 - t) * (1 - t) / 2;
|
||||||
w[1][1] = 1 - w[1][0] - w[1][2];
|
w[1][1] = 1 - w[1][0] - w[1][2];
|
||||||
|
|
||||||
mid[2] = (int)ceil(p[2] - 0.5);
|
mid[2] = (int)ceil(p[2] - 0.5f);
|
||||||
t = mid[2] - (p[2] - 0.5);
|
t = mid[2] - (p[2] - 0.5f);
|
||||||
w[2][0] = t * t / 2;
|
w[2][0] = t * t / 2;
|
||||||
w[2][2] = (1 - t) * (1 - t)/2;
|
w[2][2] = (1 - t) * (1 - t)/2;
|
||||||
w[2][1] = 1 - w[2][0] - w[2][2];
|
w[2][1] = 1 - w[2][0] - w[2][2];
|
||||||
@@ -437,20 +437,20 @@ static inline float WNoiseDy(Vec3 p, float* data) {
|
|||||||
int c[3], mid[3], n=noiseTileSize;
|
int c[3], mid[3], n=noiseTileSize;
|
||||||
float w[3][3], t, result =0;
|
float w[3][3], t, result =0;
|
||||||
|
|
||||||
mid[0] = (int)ceil(p[0] - 0.5);
|
mid[0] = (int)ceil(p[0] - 0.5f);
|
||||||
t = mid[0]-(p[0] - 0.5);
|
t = mid[0]-(p[0] - 0.5f);
|
||||||
w[0][0] = t * t / 2;
|
w[0][0] = t * t / 2;
|
||||||
w[0][2] = (1 - t) * (1 - t) / 2;
|
w[0][2] = (1 - t) * (1 - t) / 2;
|
||||||
w[0][1] = 1 - w[0][0] - w[0][2];
|
w[0][1] = 1 - w[0][0] - w[0][2];
|
||||||
|
|
||||||
mid[1] = (int)ceil(p[1] - 0.5);
|
mid[1] = (int)ceil(p[1] - 0.5f);
|
||||||
t = mid[1]-(p[1] - 0.5);
|
t = mid[1]-(p[1] - 0.5f);
|
||||||
w[1][0] = -t;
|
w[1][0] = -t;
|
||||||
w[1][2] = (1.f - t);
|
w[1][2] = (1.f - t);
|
||||||
w[1][1] = 2.0f * t - 1.0f;
|
w[1][1] = 2.0f * t - 1.0f;
|
||||||
|
|
||||||
mid[2] = (int)ceil(p[2] - 0.5);
|
mid[2] = (int)ceil(p[2] - 0.5f);
|
||||||
t = mid[2] - (p[2] - 0.5);
|
t = mid[2] - (p[2] - 0.5f);
|
||||||
w[2][0] = t * t / 2;
|
w[2][0] = t * t / 2;
|
||||||
w[2][2] = (1 - t) * (1 - t)/2;
|
w[2][2] = (1 - t) * (1 - t)/2;
|
||||||
w[2][1] = 1 - w[2][0] - w[2][2];
|
w[2][1] = 1 - w[2][0] - w[2][2];
|
||||||
@@ -480,20 +480,20 @@ static inline float WNoiseDz(Vec3 p, float* data) {
|
|||||||
int c[3], mid[3], n=noiseTileSize;
|
int c[3], mid[3], n=noiseTileSize;
|
||||||
float w[3][3], t, result =0;
|
float w[3][3], t, result =0;
|
||||||
|
|
||||||
mid[0] = (int)ceil(p[0] - 0.5);
|
mid[0] = (int)ceil(p[0] - 0.5f);
|
||||||
t = mid[0]-(p[0] - 0.5);
|
t = mid[0]-(p[0] - 0.5f);
|
||||||
w[0][0] = t * t / 2;
|
w[0][0] = t * t / 2;
|
||||||
w[0][2] = (1 - t) * (1 - t) / 2;
|
w[0][2] = (1 - t) * (1 - t) / 2;
|
||||||
w[0][1] = 1 - w[0][0] - w[0][2];
|
w[0][1] = 1 - w[0][0] - w[0][2];
|
||||||
|
|
||||||
mid[1] = (int)ceil(p[1] - 0.5);
|
mid[1] = (int)ceil(p[1] - 0.5f);
|
||||||
t = mid[1]-(p[1] - 0.5);
|
t = mid[1]-(p[1] - 0.5f);
|
||||||
w[1][0] = t * t / 2;
|
w[1][0] = t * t / 2;
|
||||||
w[1][2] = (1 - t) * (1 - t) / 2;
|
w[1][2] = (1 - t) * (1 - t) / 2;
|
||||||
w[1][1] = 1 - w[1][0] - w[1][2];
|
w[1][1] = 1 - w[1][0] - w[1][2];
|
||||||
|
|
||||||
mid[2] = (int)ceil(p[2] - 0.5);
|
mid[2] = (int)ceil(p[2] - 0.5f);
|
||||||
t = mid[2] - (p[2] - 0.5);
|
t = mid[2] - (p[2] - 0.5f);
|
||||||
w[2][0] = -t;
|
w[2][0] = -t;
|
||||||
w[2][2] = (1.f - t);
|
w[2][2] = (1.f - t);
|
||||||
w[2][1] = 2.0f * t - 1.0f;
|
w[2][1] = 2.0f * t - 1.0f;
|
||||||
|
|||||||
@@ -64,14 +64,14 @@ WTURBULENCE::WTURBULENCE(int xResSm, int yResSm, int zResSm, int amplify, int no
|
|||||||
// DG - RNA-fied _strength = 2.;
|
// DG - RNA-fied _strength = 2.;
|
||||||
|
|
||||||
// add the corresponding octaves of noise
|
// add the corresponding octaves of noise
|
||||||
_octaves = (int)(log((float)_amplify) / log(2.0f) + 0.5); // XXX DEBUG/ TODO: int casting correct? - dg
|
_octaves = (int)(log((float)_amplify) / log(2.0f) + 0.5f); // XXX DEBUG/ TODO: int casting correct? - dg
|
||||||
|
|
||||||
// noise resolution
|
// noise resolution
|
||||||
_xResBig = _amplify * xResSm;
|
_xResBig = _amplify * xResSm;
|
||||||
_yResBig = _amplify * yResSm;
|
_yResBig = _amplify * yResSm;
|
||||||
_zResBig = _amplify * zResSm;
|
_zResBig = _amplify * zResSm;
|
||||||
_resBig = Vec3Int(_xResBig, _yResBig, _zResBig);
|
_resBig = Vec3Int(_xResBig, _yResBig, _zResBig);
|
||||||
_invResBig = Vec3(1./(float)_resBig[0], 1./(float)_resBig[1], 1./(float)_resBig[2]);
|
_invResBig = Vec3(1.0f/(float)_resBig[0], 1.0f/(float)_resBig[1], 1.0f/(float)_resBig[2]);
|
||||||
_slabSizeBig = _xResBig*_yResBig;
|
_slabSizeBig = _xResBig*_yResBig;
|
||||||
_totalCellsBig = _slabSizeBig * _zResBig;
|
_totalCellsBig = _slabSizeBig * _zResBig;
|
||||||
|
|
||||||
@@ -80,7 +80,7 @@ WTURBULENCE::WTURBULENCE(int xResSm, int yResSm, int zResSm, int amplify, int no
|
|||||||
_yResSm = yResSm;
|
_yResSm = yResSm;
|
||||||
_zResSm = zResSm;
|
_zResSm = zResSm;
|
||||||
_resSm = Vec3Int(xResSm, yResSm, zResSm);
|
_resSm = Vec3Int(xResSm, yResSm, zResSm);
|
||||||
_invResSm = Vec3(1./(float)_resSm[0], 1./(float)_resSm[1], 1./(float)_resSm[2] );
|
_invResSm = Vec3(1.0f/(float)_resSm[0], 1.0f/(float)_resSm[1], 1.0f/(float)_resSm[2] );
|
||||||
_slabSizeSm = _xResSm*_yResSm;
|
_slabSizeSm = _xResSm*_yResSm;
|
||||||
_totalCellsSm = _slabSizeSm * _zResSm;
|
_totalCellsSm = _slabSizeSm * _zResSm;
|
||||||
|
|
||||||
@@ -115,9 +115,9 @@ WTURBULENCE::WTURBULENCE(int xResSm, int yResSm, int zResSm, int amplify, int no
|
|||||||
_tcTemp = new float[_totalCellsSm];
|
_tcTemp = new float[_totalCellsSm];
|
||||||
|
|
||||||
// map all
|
// map all
|
||||||
const float dx = 1./(float)(_resSm[0]);
|
const float dx = 1.0f/(float)(_resSm[0]);
|
||||||
const float dy = 1./(float)(_resSm[1]);
|
const float dy = 1.0f/(float)(_resSm[1]);
|
||||||
const float dz = 1./(float)(_resSm[2]);
|
const float dz = 1.0f/(float)(_resSm[2]);
|
||||||
int index = 0;
|
int index = 0;
|
||||||
for (int z = 0; z < _zResSm; z++)
|
for (int z = 0; z < _zResSm; z++)
|
||||||
for (int y = 0; y < _yResSm; y++)
|
for (int y = 0; y < _yResSm; y++)
|
||||||
@@ -438,13 +438,13 @@ void WTURBULENCE::resetTextureCoordinates(float *_eigMin, float *_eigMax)
|
|||||||
{
|
{
|
||||||
// allowed deformation of the textures
|
// allowed deformation of the textures
|
||||||
const float limit = 2.f;
|
const float limit = 2.f;
|
||||||
const float limitInv = 1./limit;
|
const float limitInv = 1.0f/limit;
|
||||||
|
|
||||||
// standard reset
|
// standard reset
|
||||||
int resets = 0;
|
int resets = 0;
|
||||||
const float dx = 1./(float)(_resSm[0]);
|
const float dx = 1.0f/(float)(_resSm[0]);
|
||||||
const float dy = 1./(float)(_resSm[1]);
|
const float dy = 1.0f/(float)(_resSm[1]);
|
||||||
const float dz = 1./(float)(_resSm[2]);
|
const float dz = 1.0f/(float)(_resSm[2]);
|
||||||
|
|
||||||
for (int z = 1; z < _zResSm-1; z++)
|
for (int z = 1; z < _zResSm-1; z++)
|
||||||
for (int y = 1; y < _yResSm-1; y++)
|
for (int y = 1; y < _yResSm-1; y++)
|
||||||
@@ -988,13 +988,13 @@ void WTURBULENCE::stepTurbulenceFull(float dtOrg, float* xvel, float* yvel, floa
|
|||||||
|
|
||||||
// base amplitude for octave 0
|
// base amplitude for octave 0
|
||||||
float coefficient = sqrtf(2.0f * fabs(energy));
|
float coefficient = sqrtf(2.0f * fabs(energy));
|
||||||
const float amplitude = *_strength * fabs(0.5 * coefficient) * persistence;
|
const float amplitude = *_strength * fabs(0.5f * coefficient) * persistence;
|
||||||
|
|
||||||
// add noise to velocity, but only if the turbulence is
|
// add noise to velocity, but only if the turbulence is
|
||||||
// sufficiently undeformed, and the energy is large enough
|
// sufficiently undeformed, and the energy is large enough
|
||||||
// to make a difference
|
// to make a difference
|
||||||
const bool addNoise = eigMax[indexSmall] < 2. &&
|
const bool addNoise = eigMax[indexSmall] < 2.0f &&
|
||||||
eigMin[indexSmall] > 0.5;
|
eigMin[indexSmall] > 0.5f;
|
||||||
if (addNoise && amplitude > _cullingThreshold) {
|
if (addNoise && amplitude > _cullingThreshold) {
|
||||||
// base amplitude for octave 0
|
// base amplitude for octave 0
|
||||||
float amplitudeScaled = amplitude;
|
float amplitudeScaled = amplitude;
|
||||||
@@ -1029,7 +1029,7 @@ void WTURBULENCE::stepTurbulenceFull(float dtOrg, float* xvel, float* yvel, floa
|
|||||||
// zero out velocity inside obstacles
|
// zero out velocity inside obstacles
|
||||||
float obsCheck = INTERPOLATE::lerp3dToFloat(
|
float obsCheck = INTERPOLATE::lerp3dToFloat(
|
||||||
obstacles, posSm[0], posSm[1], posSm[2], _xResSm, _yResSm, _zResSm);
|
obstacles, posSm[0], posSm[1], posSm[2], _xResSm, _yResSm, _zResSm);
|
||||||
if (obsCheck > 0.95)
|
if (obsCheck > 0.95f)
|
||||||
bigUx[index] = bigUy[index] = bigUz[index] = 0.;
|
bigUx[index] = bigUy[index] = bigUz[index] = 0.;
|
||||||
} // xyz*/
|
} // xyz*/
|
||||||
|
|
||||||
|
|||||||
@@ -969,7 +969,7 @@ static void emit_from_particles(Object *flow_ob, SmokeDomainSettings *sds, Smoke
|
|||||||
sim.scene = scene;
|
sim.scene = scene;
|
||||||
sim.ob = flow_ob;
|
sim.ob = flow_ob;
|
||||||
sim.psys = psys;
|
sim.psys = psys;
|
||||||
sim.rng = BLI_rng_new(sim.rng);
|
sim.rng = BLI_rng_new(psys->seed);
|
||||||
|
|
||||||
if (psys->part->type == PART_HAIR)
|
if (psys->part->type == PART_HAIR)
|
||||||
{
|
{
|
||||||
|
|||||||
@@ -340,7 +340,7 @@ int mathutils_deepcopy_args_check(PyObject *args)
|
|||||||
/* Mathutils Callbacks */
|
/* Mathutils Callbacks */
|
||||||
|
|
||||||
/* for mathutils internal use only, eventually should re-alloc but to start with we only have a few users */
|
/* for mathutils internal use only, eventually should re-alloc but to start with we only have a few users */
|
||||||
#define MATHUTILS_TOT_CB 13
|
#define MATHUTILS_TOT_CB 15
|
||||||
static Mathutils_Callback *mathutils_callbacks[MATHUTILS_TOT_CB] = {NULL};
|
static Mathutils_Callback *mathutils_callbacks[MATHUTILS_TOT_CB] = {NULL};
|
||||||
|
|
||||||
unsigned char Mathutils_RegisterCallback(Mathutils_Callback *cb)
|
unsigned char Mathutils_RegisterCallback(Mathutils_Callback *cb)
|
||||||
|
|||||||
Reference in New Issue
Block a user