[upadhyay.1]

Hello Friends!
I was trying to understand the code for rhopimplefoam. I was stuck at Ueqn i.e.

tmp<fvVectorMatrix> tUEqn
(
    fvm::ddt(rho, U) + fvm::div(phi, U)
  + MRF.DDt(rho, U)
  + turbulence->divDevRhoReff(U)
 ==
    fvOptions(rho, U)
);

Thing which I was unable to figure out was that how DDT term in the above equation solves for both U^* and rho^{n+1}(assuming one outer corrector). I checked in programer's guide there this ddt term was solved as ((rho*u)^{n+1} - (rho*u)^{n})/dt.
But in above we don't have values of neither rho^{n+1} nor u^{*}.
Despite from this fact we know that it solves for U^{*}.
So, my question how this is happening and from where is the value of rho^{n+1} is coming into the picture.

I hope someone will reply.
Thank You!
regards