Hi! Well, an expert will be @Charles_Driver of cstem package, I think he have deal with models like yours before. And you should check this post it talks of something similar to yours.
Another expert will be Sean Taylor of prophet but I don’t know if he is here, you can see prophet’s Stan code here. But prophet use Generalized additive Models so I think it might be a different model.
For multiple periodicities you can use BRMS and Fourier bases. If I am not wrong @avehtari and @paul.buerkner do it with GP in their article. But check it maybe helps. The used code is here.
Another approach of before (I think) will be using Fourier decomposition base, fourier() function from forecast package does it for you, It can handle multiple periodicities too and you can fit your model in Stan like a simple regression. An idea of it is in Rob’s online book
I hope this works for now! I will try to check your code today. :)