As you described, we substitute y=0 and x=R into the trajectory equation:
0=H+Rtanθ−R2g2u2sec2θ.(1)
Then, differentiating with respect to θ and setting dRdθ=0:
0=Rmaxsec2θ−R2maxg2u22sec2θtanθ,
which simplifies to
Rmax=u2gcotθ.(2)
Solving (1) and (2) will yield the desired expressions for θ and Rmax.