FFT design considerations

Super helpful. That’s where I found the formula in enough specificity that I could port it. @stevebronder just helped me figure out what I need for the template guards for types, so I should be able to roll in the custom gradients this week.

The reverse mode is all passing tests, but I’m going to add analytic autodiff. This is another case where I probably should do both var-matrix and matrix-var solutions

1 Like

Howdy folks! Just stumbled upon this thread as we were looking for a discrete FFT in Stan with the “new” complex numbers support. Do we know if and maybe this is still going to make it into a release? Thank you!

1 Like

Do you mean if FFT will be in the 2.30 release? If that is the question, then yes, the next release (RC coming out in a few hours) will include the fft() and inv_fft() as well as the 2D implementation fft2() and inv_fft2().

4 Likes

wOOOOOOOO!!! Thank you @rok_cesnovar @WardBrian @Bob_Carpenter and the rest of the team! Stan never fails to impress!

6 Likes

Does anybody know when FFT will be available in RStan?

Unfortunately not for some time yet. We’ve made a great deal of progress in having 2.26 available in RStan on CRAN, and that shouldn’t be too far off. Once that’s available, future updates should be much easier.

In the interim, if you don’t want to use cmdstanr, you can use the development experimental branch of rstan which provides the latest release:

remotes::install_git("https://github.com/stan-dev/rstan", subdir = "StanHeaders", ref = "experimental") 
remotes::install_git("https://github.com/stan-dev/rstan", subdir = "rstan/rstan", ref = "experimental")

Note that this experimental version is not guaranteed to be stable or without issue though

1 Like

Thank you for your quick reply!
So just to be sure, the experimental version also does not yet contain FFT?

I believe the experimental branch is currently at the 2.30 Release Candidate, so it should include the fft functions

Thank you once again :)

Hi again, I have tried using cmdstanr (version 2.30.0) now, but get the following error:

No matches for: 

  fft(vector)

Function fft not found.

I am probably doing something wrong (first time using CmdStanR and I’m also new to Stan and R) but I can’t see what. I am calling the fft function in my .Stan function that contains the model. Any help would be greatly appreciated!

@BirteT Are you able to share more of your R and Stan code?

Thank you for your quick reply.

The model function (referred to as test_model.stan in the subsequent code) itself is quite long, but this is the relevant section:

vector conv(vector a, vector b) {

    int la = num_elements(a);
    int lb = num_elements(b);

    vector[la+lb-1] out;

    out = inv_fft(dot_product(fft(a), fft(b))); 
    return out;

  }

The relevant code of the R function that runs the sampler is here:

  model_file <- './test_model.stan'
  mod <- cmdstan_model(model_file)
  n_warmup <- 1000
  n_iter <- 2000
  initf <- function() {
    beta_init <- array(c(rnorm(1,0,2*(0.04/1.96))),c(1))
    sigma_init <- runif(1,0,(2*10*(0.04/1.96)))
    
    list(beta = beta_init, sigma = sigma_init)
}

y <- as.numeric(unlist(read.csv(data_path, header=FALSE)))

hrf <- as.numeric(unlist(read.csv('./hrf.csv', header=FALSE)))
u <- as.numeric(unlist(read.csv(file.path(data_dir, '/u.csv'), header=FALSE)))
datal <- list(
  t = (i_ti/N_ti)^5, # TI temperature
  hrf = hrf,
  u = u,
  y = y
)

fit <- mod$sample(
  data = datal,
  chains = 4,
  parallel_chains = 4,
  refresh = 500, # print update every 500 iters
  init = initf,
  warmup = n_warmup,
  iter = n_iter,
  cores = 4,
  control=list(adapt_delta=0.9)
)

Previously, I used rstan instead of cmdstanr. Instead of inv_fft and fft I used a for-loop in conv. That version ran.

I think this code is incorrect. When I try to compile it I get an error:

Ill-typed arguments supplied to function 'inv_fft':
(complex)
Available signatures:
(complex_vector) => complex_vector
  The first argument must be complex_vector but got complex

Additionally, inv_fft will return a complex vector, not a real one

I hate to ask, but are you sure you’re running 2.30?

Right—dot products produce scalars, but FFTs operate over vectors.

Thank you both for your replies.

The output of library('cmdstanr') is:

This is cmdstanr version 0.5.2
- CmdStanR documentation and vignettes: mc-stan.org/cmdstanr
- CmdStan path: /home/birte/.cmdstan/cmdstan-2.30.0
- CmdStan version: 2.30.0

The dot product was something I had forgotten to change when implementing the fft – thank you for pointing this out. I have changed my code as follows:

vector conv(vector a, vector b) {

    int la = num_elements(a);
    int lb = num_elements(b);

    vector[la+lb-1] out;

    out = inv_fft(fft(a) .* fft(b)); 
    return out;

  }

but I still get the original error (Function fft not found.)

Can you try running the following code in a fresh R session and seeing if it compiles?

modcode  = "
vector conv(vector a, vector b) {

    int la = num_elements(a);
    int lb = num_elements(b);

    vector[la+lb-1] out;

    out = inv_fft(fft(a) .* fft(b)); 
    return out;

}
data { real y_mean; } parameters { real y; } model { y ~ normal(y_mean, 0); }
"
library(cmdstanr)
mod <- cmdstan_model(write_stan_file(modcode))

Actually, ignore that. While testing your code I found other errors.

The first:

   -------------------------------------------------
     8:        vector[la+lb-1] out;
     9:    
    10:        out = inv_fft(fft(a) .* fft(b)); 
               ^
    11:        return out;
    12:    
   -------------------------------------------------

Ill-typed arguments supplied to assignment operator =: lhs has type vector and rhs has type complex_vector

Is because you’ve declared out (and the function return) to be type vector, but the inv_fft() function returns a complex_vector.

Once I make those changes, it compiles fine for me:

modcode  = "
functions {
  complex_vector conv(vector a, vector b) {
  
      int la = num_elements(a);
      int lb = num_elements(b);
  
      complex_vector[la+lb-1] out;
  
      out = inv_fft(fft(a) .* fft(b)); 
      return out;
  
  }
}
data { real y_mean; } parameters { real y; } model { y ~ normal(y_mean, 0); }
"
mod <- cmdstan_model(write_stan_file(modcode))

1 Like

Thank you, that was it!

Apologies for the beginner questions.

1 Like

Thanks for posting and absolutely no need to apologize. The typing in Stan trips up a lot of people used to languages like R and Python that try to figure it out automatically.