The notation in your equation is ambiguous/confusing. You have the variables "U", "U_z" and "U_r". If these are dependent variables, then a solution is not possible as you have 1 PDE and 3 dependent variables. Let us suppose "U" is the only dependent variable, then one need to have more information about the coefficients "U_r" and "U_z" ; are they constants or known functions of r and t.
If you can be more precise in stating your problem, perhaps someone can help you.
Your answer is still confusing, at least to me. Normally when one says that U_r is the r-component of U, it is implied that U is a vector-valued function, i.e., U(t,r,z)=U_r(t,r,z)e_r+U_z(t,r,z)e_z, where e_r and e_z are orthogonal unit vectors in the r and z coordinate directions. Is this what you mean? Somehow I think not as f(r) looks like a scalar function. Better yet, can you given a reference where this equation comes from, or least say how you determine U_r, and U_z.
OK I now understand what you are trying to do. You can write your equation as a hyperbolic equation of the form U_t+a(t,r,z)U_r+b(t,r,z)U_z=G(r,U), where U(t,r,z) is the dependent variable and a(t,r,z)=u_r, b(t,r,z)=u_z, G(r,U)=f(r)+U^2+K are known functions and U_r=\partial U/\partial r , etc. I have used lower case u to denote the known fluid velocity so not to confuse it with U.
You might try to use the method of characteristics to solve your equation, once you have specified appropriate BCs/IC. Depending on the functional forms for the coefficients a, b and G, you may have to solve the characteristic equations numerically. Most textbooks on advance mathematics (or applied mathematics) have a chapter on the method of characteristics. I have attached some notes that should help get you started.