How to sample uniformly from 2-dim region with sinusoidal shape I have a pair of continuous random variables $\{X,\,Y\}$ that follow a 2-dim density $f_{XY}$ which is uniform over a region (denoted $\mathcal{R}$) and zero elsewhere. The region $\mathcal{R}$ is the "hill top" of a sinusoidal wave such that
$$\begin{align}
0 &< y < g(x)~ &,&& g(x) &\equiv A\cos x - \lfloor A \rfloor \\
0 &< x < x_0~ &,&& x_0 &\equiv \cos^{-1}\left( \frac{\lfloor A \rfloor}A\right) \quad \text{where}~ g(x_0) = 0
\end{align}$$
where the amplitude $A > 0$ can be any positive value except integers, and $\lfloor A \rfloor$ is the usual flooring, giving the largest integer (including zero) that is less than $A$.
$X$ and $Y$ are correlated simply because $\mathcal{R}$ is not a rectangle.

Here's My Question:

How to sample uniformly from $\mathcal{R}$ described above?

Currently I'm doing acception/rejection, but as usual it is not very efficient: generate $\{x,y\}$uniformly from the rectangle with width $x_0$ and height $g(0)$, and reject those out side of $\mathcal{R}$ where  $y > g(x)$. This rejects about $1/3\ldots$ not too bad, but I wonder if it can be improved. 
It's not clear to me if there's a transformation that can decouple $X$ and $Y$.
Another approached I have considered is to discretize the region (for given desired sample size $n$) into regular lattice points, then it becomes equivalent to a discrete sampling from $\{1,2,\ldots, n\}$. However, this feels like cheap shot, appropriate just for a certain kind of computations.
It's been many years since I last used Gibbs sampling so I'm still studying how to implement it.
 A: Generating a random variable $(X,Y)$ uniformly over that set $\mathcal{R}$ is identical to simulate a random variable $X$ with density the upper bound of $\mathcal{R}$, $g(x)$, and then a Uniform variable $Y\sim\mathcal{U}(0,g(x)$. (In our book we call this equivalence the fundamental theorem of simulation.)
To simulate from $$g(x)\propto(A\cos\{x\}-\lfloor A\rfloor)\mathbb{I}_{(0,x_0(A))}(x)$$one can first attempt to find the associated cdf: for $0<x<x_0(A)$
\begin{align*}G(x)&\propto\int_0^x g(x)\text{d}x=\int_0^x (A\cos\{x\}-\lfloor A\rfloor)\text{d}x\\
&=A\sin\{x\}-\lfloor A\rfloor x\end{align*}
but this is not easily invertible, at least analytically. (Numerically, $G$ is an increasing function, hence can be inverted by standard methods.) Note that
$$G(x_0)\propto A\sin\{x_0\}-\lfloor A\rfloor x_0$$provides the normalisation constant since $G(x_0)=1$.
For instance, here is an R code that inverts $G$:
A=pi
tA=trunc(A)
x0=acos(tA/A)
G=function(x){(A*sin(x)-tA*x)/(A*sin(x0)-tA*x0)}
inG=function(u){
dif=1/2;x=u*x0
while (dif>1e-8){
  if (G(x)>u){ x=(x-dif)*(x>dif)}else{
    x=(x+dif);if (x>x0) x=x0}
  dif=dif/2}
  return(x)}

By drawing first coordinates by cdf inversions,
and second coordinates as $\mathcal{U}(0,g(x))$, one gets uniform points under the curve g:

Another approach is to use the upper bound
$$\cos x \le 1-\frac{x^2}{2}+\frac{x^4}{4!}$$
and figure out a way to simulate from this polynomial density, a secondary accept-reject step being a solution.
The density $g$ being concave, it can also be sandwiched by piece-wise linear functions going through points $(x,g(x))$, suggesting a Gilks and Wild (1992) ARS algorithm. Without the need for a log-transform.

