Вихревой шум (curl noise)

Довольно часто возникает задача по построении движении жидкости или газа (воды, дыма и т.п.). При этом использование классических численных методов (т.е. численного решения уравнения Навье-Стокса) оказывается часто слишком сложным и дорогим. И тут на помощь нам приходят процедурные методы - они достаточно дешевы и просты и способны создавать правдоподобно выглядящую анимацию.

Фактически нам нужно задать поле скорости, т.е. вектор скорости для каждой точки. Самый простой способ - создание поля скоростей как значения шумовой функции noise. Тогда для построения анимации системы частиц мы на каждом шаге для каждой частицы по ее положению находим соответствующее значение скорости и сдвигаем частицу вдоль него.

К сожалению этот подход хотя и очень прост, но приводит к неправдоподобному движению - у нас возникают точки сгущения частиц и точки, которые частицы избегают. В результате итоговое движение оказывается не правдоподобным.

И здесь мы можем воспользоваться технологией, известной как вихревой шум (curl noise). Мы будем также использовать шумовую функцию noise, но на этот раз мы будем задавать значение скорости в точке не как значение шумовой функции, а выражать его через ее первые производные.

В двухмерном случае мы определим значение скорости v в произвольной точке как (n'y, -n'x). Поскольку чаще всего шумовая функция задается при помощи вух- или трехмерной текстуры, то мы будем считать требуемые производные численно, используя симметричную разностную схему:

Тогда процесс анимации двухмерной системы частиц при помощи описанного приема можно реализовать при помощи следующего вычислительного шейдера:

-- compute
#version 430

layout( local_size_x = 64 ) in;

uniform float deltaT;
uniform float maxDist       = 45.0;

uniform sampler2D   noiseMap;

layout(std430, binding = 0) buffer Pos
{
    vec4 position [];
};

const float eps = 1.0 / 256.0;
const vec2 dx = vec2 ( eps, 0 );
const vec2 dy = vec2 ( 0, eps );

void main()
{
    uint  idx  = gl_GlobalInvocationID.x;
    vec2  p = position [idx].xy;
    float a = textureLod ( noiseMap, 0.25*p + dy, 0 ).r - textureLod ( noiseMap, 0.25*p - dy, 0 ).r;
    float b = textureLod ( noiseMap, 0.25*p + dx, 0 ).r - textureLod ( noiseMap, 0.25*p - dx, 0 ).r;
    vec2  curl = vec2 ( a, -b ) / eps;
     
    position [idx] = vec4 ( p + curl*0.00001, 0.0, 1.0 );
}

В трехмерном случае мы также будем выражать вектор скорости через первые производные шумовой функции, только формула для скорости станет чуть более сложной:

Производные мы будем считать через ту же самую разностную схему, что и ранее. В результате для трехмерного случая мы приходим к следующему вычислительному шейдеру

-- compute
#version 430

layout( local_size_x = 64 ) in;

uniform float deltaT;
uniform float maxDist = 45.0;

uniform sampler3D   noiseMap;

layout(std430, binding = 0) buffer Pos
{
    vec4 position [];
};

const float eps = 1.0 / 64.0;
const vec3 dx   = vec3 ( eps, 0, 0 );
const vec3 dy   = vec3 ( 0, eps, 0 );
const vec3 dz   = vec3 ( 0, 0, eps );
const float scale = 0.2;

void main()
{
    uint  idx  = gl_GlobalInvocationID.x;
    vec3  p = position [idx].xyz;
    float nx   = textureLod ( noiseMap, scale * p + dx, 0 ).r - textureLod ( noiseMap, scale * p - dx, 0 ).r;
    float ny   = textureLod ( noiseMap, scale * p + dy, 0 ).r - textureLod ( noiseMap, scale * p - dy, 0 ).r;
    float nz   = textureLod ( noiseMap, scale * p + dz, 0 ).r - textureLod ( noiseMap, scale * p - dz, 0 ).r;
    vec3  curl = vec3 ( ny - nz, nz - nx,  nx - ny ) / eps;
     
    position [idx] = vec4 ( p + curl*0.0003, 1.0 );
}

Рис. Получаемая анимация трехмерной системы частиц

Можно заметить, что формула для скорости в трехмерном случае это просто применение к шумовой функции такого дифференциального оператора как ротор (curl).

Важным свойством получаемого таким образом векторного поля скоростей является то, для него ( и в 2D и в 3D) его дивергенция (дифференциальный оператор div) тождественно равен нулю. Тем самым данное поле скоростей описывает движение несжимаемой жидкости, что и дает реалистичность получаемой анимации.

В каждой точке вектор скорости будет перпендикулярен градиенту шумовой функции, т.е. направлению ее скорейшего роста. Фактически получаемое таким образом поле скоростей задает движение вдоль линии (поверхности в 3D) уровня шумовой функции.

Вместо "чистой" шумовой функции можно использовать турбулентность или функцию fBm для получения более турбулентного движения.

Примеры кода к данной статье можно скачать из репозитория https://github.com/steps3d/graphics-book.git - они находятся в папке Code/python-code.

Ссылки

https://www.cs.ubc.ca/~rbridson/docs/bridson-siggraph2007-curlnoise.pdf

https://raw.githubusercontent.com/petewerner/misc/master/Curl%20Noise%20Slides.pdf

http://petewerner.blogspot.com/2015/02/intro-to-curl-noise.html

https://www.bit-101.com/blog/2021/07/curl-noise-demystified/