/*
* __ .__ .__ ._____.
* _/ |_ _______ __|__| ____ | | |__\_ |__ ______
* \ __\/ _ \ \/ / |/ ___\| | | || __ \ / ___/
* | | ( <_> > <| \ \___| |_| || \_\ \\___ \
* |__| \____/__/\_ \__|\___ >____/__||___ /____ >
* \/ \/ \/ \/
*
* Copyright (c) 2006-2011 Karsten Schmidt
*
* This library is free software; you can redistribute it and/or
* modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation; either
* version 2.1 of the License, or (at your option) any later version.
*
* http://creativecommons.org/licenses/LGPL/2.1/
*
* This library is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with this library; if not, write to the Free Software
* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA
*/
package toxi.sim.grayscott;
import toxi.geom.Rect;
import toxi.math.MathUtils;
/**
* Implementation of the Gray-Scott reaction diffusion model described in detail
* on the links below:
* <ul>
* <li>http://groups.csail.mit.edu/mac/projects/amorphous/GrayScott/</li>
* <li>http://ix.cs.uoregon.edu/~jlidbeck/java/rd/</li>
* <li>http://www.mrob.com/pub/comp/xmorphia/</li>
* </ul>
*/
public class GrayScott {
public float[] u,
/**
*
*/
/**
*
*/
v;
protected float[] uu,
/**
*
*/
/**
*
*/
vv;
protected int width,
/**
*
*/
/**
*
*/
height;
protected float f,
/**
*
*/
/**
*
*/
k;
protected float dU,
/**
*
*/
/**
*
*/
dV;
/**
*
*/
protected boolean isTiling;
/**
*
* @param width
* @param height
* @param wrap
*/
public GrayScott(int width, int height, boolean wrap) {
this.width = width;
this.height = height;
this.u = new float[width * height];
this.v = new float[u.length];
this.uu = new float[u.length];
this.vv = new float[u.length];
this.isTiling = wrap;
reset();
// default config
setCoefficients(0.023f, 0.077f, 0.16f, 0.08f);
}
/**
* Convenience method to access U array using 2D coordinates.
*
* @param x
* @param y
* @return current u value at the given position
*/
public float getCurrentUAt(int x, int y) {
if (y >= 0 && y < height && x >= 0 && x < width) {
return u[y * width + x];
}
return 0;
}
/**
* Convenience method to access V array using 2D coordinates.
*
* @param x
* @param y
* @return current v value at the given position
*/
public float getCurrentVAt(int x, int y) {
if (y >= 0 && y < height && x >= 0 && x < width) {
return v[y * width + x];
}
return 0;
}
/**
* @return the diffuse U coefficient
*/
public float getDiffuseU() {
return dU;
}
/**
* @return the diffuse V coefficient
*/
public float getDiffuseV() {
return dV;
}
/**
* @return the F coefficient
*/
public float getF() {
return f;
}
/**
* Extension point for subclasses to modulate the F coefficient of the
* reaction diffusion, based on spatial (or other) parameters. This method
* is called for every cell/pixel of the simulation space from the main
* {@link #update(float)} cycle and can be used to create parameter
* gradients, animations and other spatial or temporal modulations.
*
* @param x
* @param y
* @return the active F coefficient at the given position
*/
public float getFCoeffAt(int x, int y) {
return f;
}
/**
* @return the K coefficient
*/
public float getK() {
return k;
}
/**
* Extension point for subclasses to modulate the K coefficient of the
* reaction diffusion, based on spatial (or other) parameters. This method
* is called for every cell/pixel of the simulation space and can be used to
* create parameter gradients, animations and other spatial or temporal
* modulations.
*
* @param x
* @param y
* @return the active K coefficient at the given position
*/
public float getKCoeffAt(int x, int y) {
return k;
}
/**
* @return the isTiling
*/
public boolean isTiling() {
return isTiling;
}
/**
* Resets the simulation matrix to an initial, clean state.
*/
public final void reset() {
for (int i = 0; i < uu.length; i++) {
uu[i] = 1.0f;
vv[i] = 0.0f;
}
}
/**
*
* @param pixels
* @param imgWidth
* @param imgHeight
*/
public void seedImage(int[] pixels, int imgWidth, int imgHeight) {
int xo = MathUtils.clip((width - imgWidth) / 2, 0, width - 1);
int yo = MathUtils.clip((height - imgHeight) / 2, 0, height - 1);
imgWidth = MathUtils.min(imgWidth, width);
imgHeight = MathUtils.min(imgHeight, height);
for (int y = 0; y < imgHeight; y++) {
int i = y * imgWidth;
for (int x = 0; x < imgWidth; x++) {
if (0 < (pixels[i + x] & 0xff)) {
int idx = (yo + y) * width + xo + x;
uu[idx] = 0.5f;
vv[idx] = 0.25f;
}
}
}
}
/**
*
* @param f
* @param k
* @param dU
* @param dV
*/
public final void setCoefficients(float f, float k, float dU, float dV) {
this.f = f;
this.k = k;
this.dU = dU;
this.dV = dV;
}
/**
* @param dU
* the diffuse U coefficient to set
*/
public void setDiffuseU(float dU) {
this.dU = dU;
}
/**
* @param dV
* the diffuse V coefficient to set
*/
public void setDiffuseV(float dV) {
this.dV = dV;
}
/**
* @param f
* the F coefficient to set
*/
public void setF(float f) {
this.f = f;
}
/**
* @param k
* the K coefficient to set
*/
public void setK(float k) {
this.k = k;
}
/**
* @param x
* @param y
* @param w
* @param h
*/
public void setRect(int x, int y, int w, int h) {
int mix = MathUtils.clip(x - w / 2, 0, width);
int max = MathUtils.clip(x + w / 2, 0, width);
int miy = MathUtils.clip(y - h / 2, 0, height);
int may = MathUtils.clip(y + h / 2, 0, height);
for (int yy = miy; yy < may; yy++) {
for (int xx = mix; xx < max; xx++) {
int idx = yy * width + xx;
uu[idx] = 0.5f;
vv[idx] = 0.25f;
}
}
}
/**
*
* @param r
*/
public void setRect(Rect r) {
setRect((int) r.x, (int) r.y, (int) r.width, (int) r.height);
}
/**
* @param isTiling
* the isTiling to set
*/
public void setTiling(boolean isTiling) {
this.isTiling = isTiling;
}
/**
*
* @param t
*/
public void update(float t) {
t = MathUtils.clip(t, 0, 1f);
int w1 = width - 1;
int h1 = height - 1;
for (int y = 1; y < h1; y++) {
for (int x = 1; x < w1; x++) {
int idx = y * width + x;
int top = idx - width;
int bottom = idx + width;
int left = idx - 1;
int right = idx + 1;
float currF = getFCoeffAt(x, y);
float currK = getKCoeffAt(x, y);
float currU = uu[idx];
float currV = vv[idx];
float d2 = currU * currV * currV;
u[idx] = MathUtils
.max(0,
currU
+ t
* ((dU
* ((uu[right] + uu[left]
+ uu[bottom] + uu[top]) - 4 * currU) - d2) + currF
* (1.0f - currU)));
v[idx] = MathUtils
.max(0,
currV
+ t
* ((dV
* ((vv[right] + vv[left]
+ vv[bottom] + vv[top]) - 4 * currV) + d2) - currK
* currV));
}
}
if (isTiling) {
int w2 = w1 - 1;
int idxH1 = h1 * width;
int idxH2 = (h1 - 1) * width;
for (int x = 0; x < width; x++) {
int left = (x == 0 ? w1 : x - 1);
int right = (x == w1 ? 0 : x + 1);
int idx = idxH1 + x;
float cu = uu[x];
float cv = vv[x];
float cui = uu[idx];
float cvi = vv[idx];
float d = cu * cv * cv;
u[x] = MathUtils
.max(0,
cu
+ t
* ((dU
* ((uu[right] + uu[left]
+ uu[width + x] + cui) - 4 * cu) - d) + f
* (1.0f - cu)));
v[x] = MathUtils
.max(0,
cv
+ t
* ((dV
* ((vv[right] + vv[left]
+ vv[width + x] + cvi) - 4 * cv) + d) - k
* cv));
d = cui * cvi * cvi;
u[idx] = MathUtils
.max(0,
cui
+ t
* ((dU
* ((uu[idxH1 + right]
+ uu[idxH1 + left] + cu + uu[idxH2
+ x]) - 4 * cui) - d) + f
* (1.0f - cui)));
v[idx] = MathUtils
.max(0,
cvi
+ t
* ((dU
* ((vv[idxH1 + right]
+ vv[idxH1 + left] + cv + vv[idxH2
+ x]) - 4 * cvi) + d) - k
* cvi));
}
for (int y = 0; y < height; y++) {
int idx = y * width;
int idxW1 = idx + w1;
int idxW2 = idx + w2;
float cu = uu[idx];
float cv = vv[idx];
float cui = uu[idxW1];
float cvi = vv[idxW1];
float d = cu * cv * cv;
int up = (y == 0 ? h1 : y - 1) * width;
int down = (y == h1 ? 0 : y + 1) * width;
u[idx] = MathUtils
.max(0,
cu
+ t
* ((dU
* ((uu[idx + 1] + cui
+ uu[down] + uu[up]) - 4 * cu) - d) + f
* (1.0f - cu)));
v[idx] = MathUtils
.max(0,
cv
+ t
* ((dV
* ((vv[idx + 1] + cvi
+ vv[down] + vv[up]) - 4 * cv) + d) - k
* cv));
d = cui * cvi * cvi;
u[idxW1] = MathUtils
.max(0,
cui
+ t
* ((dU
* ((cu + uu[idxW2]
+ uu[down + w1] + uu[up
+ w1]) - 4 * cui) - d) + f
* (1.0f - cui)));
v[idxW1] = MathUtils
.max(0,
cvi
+ t
* ((dV
* ((cv + vv[idxW2]
+ vv[down + w1] + vv[up
+ w1]) - 4 * cvi) + d) - k
* cvi));
}
}
System.arraycopy(u, 0, uu, 0, u.length);
System.arraycopy(v, 0, vv, 0, v.length);
}
}