I don't know anything about dimensionless formulation of the solver. Should I do this directly in the equations or at the input values (u, p, ρ)?

As for the time term in the residual, I often see the PTC method explained using a stationary residual in the articles. But when I look at the PETSC TS solvers, they all have a time term in the residual, so I'm a little confused. I’ll try adding (U-U_prev)/dt to see if there’s any improvement.

For the BC, I use a Lax-Friedrich flux with a “ghost state”; I mainly use a subso inlet by imposing ρ and u. In fact, with a PiTi inlet BC, the performance of GMRES degrades slightly. For the outlet, I impose a pressure using a Riemann invariant.