#include <cstdio>
#include <initializer_list>
#define __device__
#define __forceinline__ inline
#include "boundary/planar_velocity.cuh"

int main() {
    const int cx[9]={0,1,-1,0,0,1,-1,1,-1};
    const int cy[9]={0,0,0,1,-1,1,-1,-1,1};
    const float w[9]={4.f/9,1.f/9,1.f/9,1.f/9,1.f/9,1.f/36,1.f/36,1.f/36,1.f/36};
    for (int top : {0,1}) for (float rho : {1.f,1.04f})
    for (float u : {-.02f,0.f,.02f}) for (float v : {-.01f,-.001f,0.f,.001f,.01f})
    for (int neq : {0,1}) {
        float f[9];
        printf("%d %.9g %.9g %.9g %d\n",top,rho,u,v,neq);
        for (int i=0;i<9;++i) {
            float cu=cx[i]*u+cy[i]*v;
            f[i]=w[i]*rho*(1+3*cu+4.5f*cu*cu-1.5f*(u*u+v*v));
            f[i]+=neq*.0001f*(i-4);
            // Unknown inputs must be completely overwritten, never consumed.
            if (cy[i]==(top?-1:1)) f[i]=-7.f-i;
            printf("%.9g ",f[i]);
        }
        puts("");
        if (top) reconstruct_top_velocity(f,u,v);
        else reconstruct_bottom_velocity(f,u,v);
        for (float value : f) printf("%.9g ",value);
        puts("");
    }
}