Question about simplex-noise

316cf5620d41dfc55bf42e879a52f457
0
tuohirv 101 Aug 05, 2010 at 07:01

Hi,

I’ve been working with a project requiring a good noise-function once again. I implemented the basic simplex-noise instead of perlin-noise this time. I’m attracted by a fact that there should be no visible artifacts (The main problem with every perlin-noise i’ve written)… I checked the reference implementation at:http://staffwww.itn.liu.se/\~stegu/simplexnoise/simplexnoise.pdf and used it as a base for my implementation. I removed floating points and static, constant gradients and replaced them with full-random - gradients…

When i was reimplementing the noise-function. After i transformed the original coordinates into “hypercube”-space… i spent some time trying to do the interpolation in this space directly. Instead of transforming the corners back to the real space. I managed to implement the noise with linear interpolation….. The problem is, no matter what i tried. I couldn’t make it look any better (than linear). I have a strong feeling that the whole function could be implemented without backward-transform (skew)… I don’t have any mathematical base to say so though. If this is possible, the simplex-noise implementations would increase significantly in speed.

I was wondering that maybe someone have been fighting with the same issue? It would be enough for me to know if it can, or if it cannot be done.

Thanks, Tuomo

5 Replies

Please log in or register to post a reply.

5225bc0c3bf66f4c275c332de6388d1f
0
SyntaxError 101 Nov 16, 2010 at 04:31

I’m not a math wiz but the skew seems pretty cheap to me. You really think it would save that much time not to do it? I would guess that it might take more time to do it without the skew because the contribution calculation is now harder and not the same for every corner. I also have nothing to back this up but for instance in 2D your tidy little circles become more like ellipses for two of the three corners. However if you figure out anything cool I’d be happy to hear about it :-) I use simplex noise for everything.

316cf5620d41dfc55bf42e879a52f457
0
tuohirv 101 Dec 07, 2010 at 11:01

Well, all of the separate operations seems cheap I agree. However, when ran with CPU only (as in my case) I would like the algorithm to be much faster. Its true that the skew/unskew is fast way to transform the space into more usable one for this case, but if it can be done without them, directly, it could speedup the algorithm.. Any operations on the lowest level can be significant to the speed. I’m well aware that the algorithm has been studied and widely used a lot but still I wouldnt want to give up hope it can be made even better.. Its difficult to just “accept” there are algorithms that are already as good as they possible can be :)

5225bc0c3bf66f4c275c332de6388d1f
0
SyntaxError 101 Dec 07, 2010 at 18:56

Funny thing is, after I posted a reply I stated working on something to speed SN up for myself because it turns out things were too slow for me too :) (At least on the GPU side). What I’m trying is pre-calcuating cubes on the CPU that do most of the claculation and then using trilinear interpoation on those cubes to get the results on the GPU. This is probably only useful on a GPU but if that’s where you need it, I’ll post the results once (or rather if) I get it working properly.

316cf5620d41dfc55bf42e879a52f457
0
tuohirv 101 Jan 19, 2011 at 08:10

for reference:

Here is my simplex-noise function:

inline int dot(const unsigned short packed_grad, const int x, const int y ) {
   return (((packed_grad&255)-127)*x + (((packed_grad>>8)&255)-127)*y)>>7;
};

int CTNoise::simplexNoise2D(int x, int y) {

            // convert multipliers from 16bit fixedpoint to tnoise-fixedpoint
    int skew_mul = (23987*TNOISE_VAL)>>16;
    int unskew_mul = (13849*TNOISE_VAL)>>16;


    int temp = (((x+y) * skew_mul)>>TNOISE_BITS);
    // skewed coordinates.
    int mx = ((x+temp)>>TNOISE_BITS);
    int my = ((y+temp)>>TNOISE_BITS);


    int point0_x = mx<<TNOISE_BITS;
    int point0_y = my<<TNOISE_BITS;
    temp = ((point0_x + point0_y) * unskew_mul)>>TNOISE_BITS;
    point0_x = x - (point0_x-temp);
    point0_y = y - (point0_y-temp);



    int rval = 0;

           // FIRST
    temp = (TNOISE_VAL>>1) - ((point0_x*point0_x + point0_y*point0_y)>>TNOISE_BITS);
    if (temp>0) {
        temp = ((temp*temp*temp)>>15);  // gives maximum multiplier of 4096
        rval+= dot( m_rtable[(my + m_rtable[(mx)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND], point0_x, point0_y ) * temp;
    };

        // LAST
    int point1_x = point0_x - TNOISE_VAL + 2*unskew_mul;
    int point1_y = point0_y - TNOISE_VAL + 2*unskew_mul;
    temp = (TNOISE_VAL>>1) - ((point1_x*point1_x + point1_y*point1_y)>>TNOISE_BITS);
    if (temp>0) {
        temp = ((temp*temp*temp)>>15);
        rval+= dot( m_rtable[(my + 1 + m_rtable[(mx+1)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND], point1_x, point1_y ) * temp;
    };


    point1_x = point0_x + unskew_mul;
    point1_y = point0_y + unskew_mul;

        // MIDDLE point
    if (point0_x > point0_y) {          // lower triangle
        point1_x -= TNOISE_VAL;

        temp = (TNOISE_VAL>>1) - ((point1_x*point1_x + point1_y*point1_y)>>TNOISE_BITS);
        if (temp>0) {
            temp = ((temp*temp*temp)>>15);
            rval+= dot( m_rtable[(my + m_rtable[(mx+1)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND], point1_x, point1_y ) * temp;
        };


    } else {                            // upper triangle
        point1_y -= TNOISE_VAL;

        temp = (TNOISE_VAL>>1) - ((point1_x*point1_x + point1_y*point1_y)>>TNOISE_BITS);
        if (temp>0) {
            temp = ((temp*temp*temp)>>15);
            rval+= dot( m_rtable[(my + 1 + m_rtable[(mx)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND], point1_x, point1_y ) * temp;
        };
    }
    return ((rval * 13 )>>8);
};

And here is a version with linear interpolation (only).. The interpolation is done in hypercube-space,…

int CTNoise::simplexNoise2D_linear(int x, int y) {
    int temp = (((x+y)*94)>>8); // unskew multiply by 23/64
    x+=temp;
    y+=temp;
    int mx = (x>>TNOISE_BITS);
    int my = (y>>TNOISE_BITS);
    x&=TNOISE_AND;
    y&=TNOISE_AND;

    if (x>y) {
        x = TNOISE_VAL-x;
        temp = (TNOISE_VAL - (x+y));
        if (temp>0) temp = ((int)m_rtable[((my) + m_rtable[(mx+1)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND] * temp); else temp = 0;
        temp += ((int)m_rtable[((my) + m_rtable[(mx)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND]*x);
        temp += ((int)m_rtable[((my+1) + m_rtable[(mx+1)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND]*y);
    } else {
        y = TNOISE_VAL-y;
        temp = TNOISE_VAL - (x+y);
        if (temp>0) temp = ((int)m_rtable[((my+1) + m_rtable[(mx)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND] * temp); else temp = 0;
        temp += ((int)m_rtable[((my) + m_rtable[(mx)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND] * y);
        temp += ((int)m_rtable[((my+1) + m_rtable[(mx+1)&TNOISE_RTABLE_AND])&TNOISE_RTABLE_AND] * x);
    }
    return -65535 + (temp>>(TNOISE_BITS-1));
};

As you know, the difference of complexity (and performance) increases exponentially when dimensions are added.

6eaf0e08fe36b2c23ca096562dd7a8b7
0
__________Smile_ 101 Jan 25, 2011 at 02:35

Interpolation can be effective done in transformed space if gradient table has entries in the same space.

Squared distance in such transformed space can be computed with following formula (3D case):
r2 = Vx\^2 + Vy\^2 + Vz\^2 - g*(Vx + Vy + Vz)\^2
g = (1 - a\^2)/n, a – scaling factor along main diagonal, n = 3

Full noise can be computed effectively with following algorithm.

Consider B[0], B[1], B[2] – fractional part of the transformed coordinates and K[0], K[1], K[2] – corresponding indexes.

First, determine order of the values in B[], say I[0] – index of the smallest, I[1] – middle and I[2] – largest value from array.
This step effectively is choosing correct simplex from hypercube.

Second, calculate start squared distance
R2 = B[0]\^2 + B[1]\^2 + B[1]\^2 - g*(B[0] + B[1] + B[2])\^2,
and accumulator
C = 2*g*(B[0] + B[1] + B[2]) - g + 1.
Use this values for the first contribution to final result
Z = F(R2) * (grad[K] * B).

Then iterate, i = 0..n-1.
R2 += C - 2*B[I], C -= 2*g,
K[I]++, B[I] -= 1,
Z += F(R2) * (grad[K] * B).

I think this implementation can be quite fast even in large dimension case.