Implement the derivative of a function in stancode in mathematica (cosmology)

Modeling
1

Hi!

I defined some functions in a stan-mathematica interface to estimate parameters of a cosmological inhomogeneous model.

My problem is that i need some derivatives of those functions and to solve a differential equation to obtain the last parameters. For example i defined my most importan function (scale factor) in terms of other previous functions:

 real[] R(real r,         // Comovil Position
		   real t,         // Comovil Time
                 real H0 ,       // Local Hubble Rate
                 real Omegaout,  // Omega out
                 real Omegain,   // Omega in
                 real r0,        // Void size parameter
                 real Deltar     // Change parameter    
) {
	return (R0(r) / contrast(r, Omegaout, Omegain, r0, Deltar)) *
     InverseSerie((2/3) * (contrast(r, Omegaout, Omegain, r0, Deltar) *
     a(r, t, H0, Omegaout, Omegain, r0, Deltar))^(3/2));
  }

And now i need to evaluate this function in some r(z) and t(z) values given by this system of equations that i could solve with mathematica, but i dont know how to implement in the stan code:

DtG[z_] := 
 D[R[rG[z], tG[z]], 
   rG[z]]/((1 + z) D[D[R[rG[z], tG[z]], rG[z]], tG[z]])
DrG[z_] := 
 Sqrt[1 - k[rG[z]]]/((1 + z) D[D[R[rG[z], tG[z]], rG[z]], tG[z]])

trsol = NDSolve[{tG'[z] == -DtG[z], rG'[z] == DrG[z], rG[0] == 0, 
   tG[0] == 1/H0}, {tG, rG}, {z, 0, 1.74}]

It gives me the values to perform the evaluation of the R (scale factor function) in the r(z), t (z) solutions and obtain R as function of z (and this is the function that i need to perform bayesian analisis)

Can i import the function that i have obtained in mathematica via NDSolve? There is a way to implement this in the “StanCode” interface for mathematica? I read some information about autodiff useful, but i dont know how to apply in my case.

I would very grateful if someone can give me an orientation cause im a little lost and i like stan as bayesian MC estimator!!

Thanks!

1 Like
2

Just to close this issue, i found that the best solution is to perform the deruvative in mathematica and then write this in c languaje via Cform[ ] command. It will gives you the function and you can add it to the stan code exporting to a file “func.txt” and then just importing and pasting this text in the stan code like:

“… return” <> file <> ;
}

Otherfunction(…"

Hope it helps to someone

1 Like