d_reversed = reversed >> (32u - bits);\\n let upper_bits = (num >> bits) << bits;\\n return shifted_reversed + upper_bits;\\n}\\n\\nfn getIndexPair(i: u32, it: u32) -> vec2 {\\n let k1 = reverseLowestBits(2u * i, it + 1u);\\n let k2 = k1 + (1u << it);\\n return vec2u(k1, k2);\\n}\\n\\nfn expi(angle: float) -> vec2f {\\n return vec2f(cos(angle), sin(angle));\\n}\\n\\nfn cmul(a: vec2f, b: vec2f) -> vec2f {\\n return vec2f(a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x);\\n}\\n\\nfn unityRoot(k: u32, N: u32, inverse: bool) -> vec2f {\\n let d = select(-1.0, 1.0, inverse);\\n return expi(2.0 * d * PI * float(k) / float(N));\\n}\\n\\nfn fftshift(index: u32) -> u32 \\n{\\n return (index + SIZE / 2) % SIZE;\\n}\\n\\nfn ifftshift(index: u32) -> u32 \\n{\\n return (index + (SIZE + 1) / 2) % SIZE;\\n}\\n\\nvar TEMP: array;\\n\\nfn fft(index: u32, group: u32, axis: u32, inverse: bool) {\\n for (var iteration = 0u; iteration <= u32(IT_NUM); iteration++)\\n {\\n var ids = getIndexPair(index, select(iteration, IT_NUM, iteration == 0u));\\n let v1 = select(TEMP[ids.x], sim[getAxisIndex(ids.x, group, axis)], iteration == 0u);\\n let v2 = select(TEMP[ids.y], sim[getAxisIndex(ids.y, group, axis)], iteration == 0u);\\n\\n ids = getIndexPair(index, iteration);\\n let rootIndex = (ids.x & ((1u << iteration) - 1u)) << (u32(IT_NUM) - iteration);\\n let V = cmul(unityRoot(rootIndex, SIZE, inverse), v2);\\n TEMP[ids.x] = v1 + V;\\n TEMP[ids.y] = v1 - V;\\n\\n workgroupBarrier();\\n }\\n\\n let M = SIZE / WG_SIZE;\\n for (var i = 0u; i < u32(M); i++) {\\n let rowIndex = index + i * WG_SIZE;\\n let idx = getAxisIndex(rowIndex, group, axis);\\n sim[idx] = TEMP[rowIndex] / select(1.0, SIZE, inverse);\\n }\\n\\n storageBarrier();\\n}\\n\\nfn gaussian(dx: vec2f, sigma: float) -> float {\\n return exp(-dot(dx, dx) / (2.0 * sigma * sigma));\\n}\\n\\nfn real_space_update(index: u32)\\n{\\n let id = indexTo2D(index);\\n if(time.frame < 5) \\n {\\n let uv = (vec2f(id.xy) + 0.5) / SIZE;\\n let len = length(vec2f(1.6,1.0)*(uv-vec2f(0.5,0.5)));\\n var color = textureSampleLevel(channel0, bilinear, uv, 0).rgb * 0.5;\\n sim[index] = 1.02*vec2f(1.0 - smoothstep(0.52,0.55,len), 0);\\n } \\n else\\n {\\n //apply potential energy time evolution operator\\n var field = sim[index];\\n for(var i = 0; i < SUBSTEPS; i++)\\n {\\n let dt = TIME_STEP / float(SUBSTEPS);\\n let VA = PotentialAbsorb(field, id.xy);\\n let potentialUpdate = expi(-dt * VA.x / (2*H_PLANK)) * exp(-dt * VA.y / (2*H_PLANK));\\n field = cmul(potentialUpdate, field);\\n }\\n sim[index] = field;\\n }\\n}\\n\\nfn momentum_space_update(index: u32)\\n{\\n let id = indexTo2D(index);\\n let freq2 = vec2i(vec2u(ifftshift(id.x), ifftshift(id.y))) - SIZE/2;\\n let freq = length(vec2f(freq2));\\n //apply kinetic energy time evolution operator\\n let kineticUpdate = expi(-TIME_STEP*H_PLANK*freq*freq/(2.0*MASS));\\n sim[index] = cmul(kineticUpdate, sim[index]) * (1.0 - smoothstep(LOW_PASS_FREQ-1.0,LOW_PASS_FREQ,freq));\\n}\\n\\n#workgroup_count mega_kernel SIZE 1 1\\n@compute @workgroup_size(WG_SIZE)\\nfn mega_kernel(\\n @builtin(local_invocation_index) local_index: u32,\\n @builtin(workgroup_id) workgroup_id: vec3u\\n) { \\n let linearIndex = local_index + workgroup_id.x * SIZE;\\n let currentStep = dispatch.id / DISPATCHES_PER_STEP;\\n var currentDispatch = dispatch.id % DISPATCHES_PER_STEP;\\n\\n if(currentDispatch == 0) { real_space_update(linearIndex); real_space_update(linearIndex+WG_SIZE); return; }\\n currentDispatch-=1;\\n if(currentDispatch < AXIS_COUNT) { fft(local_index, workgroup_id.x, currentDispatch, false); return; }\\n currentDispatch-=AXIS_COUNT;\\n if(currentDispatch == 0) { momentum_space_update(linearIndex); momentum_space_update(linearIndex+WG_SI

@compute.toys