This does not seem very trivial and I have not tested this at my end but the following can be tried:
1.) Convert the Neumann BC to Dirichlet type by breaking it into explict type as follows: P_n+1-P_n = c*delta_t*(P_x+1,n-P_x,n)/delta_x. Here P_n+1 is at the current time step and P_n is being used from previous time step. Yes, this will need multiple iterations to converge properly. P_n+1 can be found and applied as Dirichlet BC. The gradient can be found using in-built gradient macros. DEFINE_ADJUST or even DEFINE_PROFILE could be used here.
2.) Freeze the flow field once a single time step has been solved. This can be done by using TUI commands to unselect the flow and turbulence equations. Just have the UDS equation turned on.
3.) Decrease the time step to 1/10th of the existing flow time and then solve only the UDS equation till convergence.
4.) Once it is converged or desired number of UDS iterations have been achieved, unfreeze the flow and turbulence equations again.
5.) Change back the time step to your original flow time step.
6.) Also adjust the flow time accordingly as the UDS iterations were considered to be internal within one flow time. You can find relevant TUI commands to adjust this.
7.) Continue the steps and converge the flow field accordingly.
All of the above steps can be automated using scheme, TUI and Execute commands.
Since I have not tested them, so I cannot guarantee any accuracy or feasibility of the above solution.
You will need to test them using a very simple model and tweak the steps accordingly if needed.
I hope this helps.