#ifndef NAVIERSTOKES_H
#define NAVIERSTOKES_H

#include <complex.h>
#include <fftw3.h>

// to extract elements from array of size S representing a function of momentum, use
//   array[KEXTRACT(kx,ky,size)]
#define KLOOKUP(X,Y,S) (X>=0?X:S+X)*S+(Y>=0?Y:S+Y)


// parameters for the NS equation
typedef struct ns_params {
  // number of modes
  int K;
  // 2*K+1
  int S;
  // 4*K+1
  int N;
  // forcing term
  _Complex double* g;
  // time step
  double h;
  // friction
  double nu;
} ns_params;

// arrays on which the ffts are performed
typedef struct fft_vects {
  fftw_complex* fft1;
  fftw_complex* fft2;
  fftw_complex* invfft;
  fftw_plan fft1_plan;
  fftw_plan fft2_plan;
  fftw_plan invfft_plan;
} fft_vects;

// next time step for Irreversible Navier-Stokes equation
int ins_step(_Complex double* u, ns_params params, fft_vects vects, _Complex double* tmp1, _Complex double* tmp2, _Complex double* tmp3);

// right side of Irreversible Navier-Stokes equation
int ins_rhs(_Complex double* out, _Complex double* u, ns_params params, fft_vects vects);

// compute alpha
_Complex double compute_alpha(_Complex double* u, ns_params params);

#endif