GollyGang / ready

A cross-platform implementation of various reaction-diffusion systems and PDEs.
GNU General Public License v3.0
780 stars 60 forks source link

Support block size in formula image rules #125

Closed timhutton closed 3 years ago

timhutton commented 3 years ago

This closes #124

Supports 4x1x1 (float4 or double4) and 1x1x1 (float or double).

@danwills I hope this is useful for your Houdini efforts. Now you should be able to call system->SetBlockSizeX(1) and then system->GetKernel() will return: (e.g. for advection.vti)

__kernel void rd_compute(__global float *a_in,__global float *a_out)
{
    // parameters:
    const float dx = 0.10000000f;
    const float timestep = 0.00010000f;

    // keywords needed for the formula:
    const int index_x = get_global_id(0);
    const int index_y = get_global_id(1);
    const int index_z = get_global_id(2);
    const int X = get_global_size(0);
    const int Y = get_global_size(1);
    const int Z = get_global_size(2);
    const int index_here = X*(Y*index_z + index_y) + index_x;
    const float a_w = a_in[X* (Y * index_z + index_y) + ((index_x-1 + X) & (X - 1))];
    float a = a_in[index_here];
    const float a_e = a_in[X* (Y * index_z + index_y) + ((index_x+1 + X) & (X - 1))];
    const float x_gradient_a = (-1 * a_w + a_e) / (2 * dx);
    float delta_a = 0.0f;

    // the formula:
    delta_a = -x_gradient_a;

    // forward-Euler update step:
    a_out[index_here] = a + timestep * delta_a;
}

Let me know how else to make the interop with Houdini easier. Would it be useful to have a system->GetFormulaKeywords() that returned this part:

    const float x_gradient_a = (-1 * a_w + a_e) / (2 * dx);

Or something like that?

danwills commented 3 years ago

That looks absolutely perfect Tim! I'll merge your pull request and check it out soon. One thing that I'm currently doing in the plugin is stripping out all the const's.. but maybe they will actually be ok once it's not float4. I am intending to try to use the whole kernel that Ready builds, so there's no need for a single method that returns the keywords section any more. Reading the code it did seem like it might be a good thing to do to it anyway though? (even just for organisational/modularity reasons?)