If the $R_k$ are your only optimization variables, then the constraint
$$
\exp \left[ - (2^{{R}_k } -1) \left(\frac{\tilde{Z} g_{k} p_{\max} + \sigma^2}{g_{pu} p_{pu}} \right) \right] \leq q
$$
is equivalent to
$$
- (2^{{R}_k } -1) \left(\frac{\tilde{Z} g_{k} p_{\max} + \sigma^2}{g_{pu} p_{pu}} \right) \leq \ln q
$$
which is equivalent to
$$
(2^{{R}_k } -1)
\leq
-(\ln q)
\left(\frac{g_{pu} p_{pu}}{\tilde{Z} g_{k} p_{\max} + \sigma^2} \right)
$$
which is equivalent to
$$
2^{{R}_k }
\leq
-(\ln q)
\left(\frac{g_{pu} p_{pu}}{\tilde{Z} g_{k} p_{\max} + \sigma^2} \right)+1
$$
which is equivalent to
$$
{R}_k
\leq
\log_2
\left\{
-(\ln q)
\left(\frac{g_{pu} p_{pu}}{\tilde{Z} g_{k} p_{\max} + \sigma^2} \right)+1\right\}.
$$
This is a linear constraint, and your objective function is already linear. So you ended with a standard linear program, for which there is plenty of software (and if you don't want to use any of that, plenty of simple to implement methods).
The transformation above shows an important principle: one can often transform complicated-looking problems into much simpler problems if one just looks at them closely enough.