When pricing a European call option using the Monte Carlo simulation when the interest rate is constant I can pull the discounting factor out of the expected value giving: C=exp(-r*(T-t))*E(max(S ? K, 0)/Ft)And then I can use the Monte Carlo method to solve for the expected value of the future price of the underlying which is the only unknown and plug it in.However when the interest rate is not constant I cannot pull out the discounting factor. Therefore I end up with the exponential of the integral of the interest rate inside the expected value. From looking around online I've noticed that mostly in this case a solution analogous to the Black-Scholes is guessed, plugged into the PDE of the model and then the unknowns are solved for. Is this the only procedure for pricing the call option under stochastic interest rate models? Can I use the Monte Carlo method to somehow solve for the integral of the interest rate?Thanks for any help

The direct answer to your question is that it is easy to simulatesomething like B = exp{-int(0,T) r(t) dt} and C = E[B (S(T) - K)^+] as part of your Monte Carlo.Say you have an sde for r(t) of the formdr = b(r,t) dt + a(r,t) dWLet Y = int(0,T) r(t) dtThen a pseudo-code fragment would beY=0; r = r0; S = S0;For(i=1, i <= nsteps, i++){Y += r dt;Z = draw of random normal variate;r += b(r,t) dt + a(r,t) sqrt(dt) Z;// add previous code stock price simulation here} B = exp(-Y)C = B Max(S-K,0)// repeat over nsims

Alan, thanks for your reply. However, it seems to me that you are taking r(t) and S(t) to be independent so therefore the covariance is zero and E[B (S(T) - K)^+] = E * E[(S(T) - K)^+] . In my model the two are correlated. Am I misconstruing what you did?

I considered that possibility, and it still works in my framework.In the section where I say "add previous code stock price simulation here", you just have to be careful to pick up the correlation, call it rho.For example, if the previous stock price simulation was under GBM, you could use:Z1 = Second draw of independent random normal deviateZ2 = rho Z + sqrt{1- rho^2} Z1S += r S dt + S sigma sqrt(t) Z2p.s. You actually want to insert these lines before advancing r(t), so that r0 is used for the first S advance.