Rational Number RNG

JavaScript (ES6), 81 bytes

Returns \$p/q\$ as [p,q].

f=(p=0,q=1)=>Math.random()<(g=(a,b)=>b?g(b,a%b):a)(p,q)/2?f(++p-q?p:!++q,q):[p,q]

Try it online!

Or Try some statistics!

How?

We recursively build the fractions \$p/q\$ in the relevant order, but without any explicit deduplication, until a random number picked in \$[0,1[\$ is greater than or equal to \$gcd(p,q)/2\$.

If \$p\$ and \$q\$ are not coprime, we have \$gcd(p,q)/2\ge1\$, which always triggers the recursive call. This is the expected behavior because the corresponding fraction must be discarded.

If \$p\$ and \$q\$ are coprime, we have \$gcd(p,q)/2=1/2\$, which is the desired probability of triggering the next call.

Commented

f = (                 // f is a recursive function taking:
  p = 0,              //   p = numerator, starting at 0
  q = 1               //   q = denominator, starting at 1
) =>                  //
  Math.random() <     // pick a random number in [0,1[
  (g = (a, b) =>      // g is a recursive function computing the gcd of 2 integers
    b ?               //   if b is not equal to 0:
      g(b, a % b)     //     recursive call with (b, a mod b)
    :                 //   else:
      a               //     return a
  )(p, q) / 2 ?       // if our random number is less than gcd(p, q) / 2:
    f(                //   do a recursive call:
      ++p - q ? p     //     increment p if p < q - 1
              : !++q, //     otherwise: set p to 0 and increment q
      q               //
    )                 //   end of recursive call
  :                   // else:
    [p, q]            //   stop recursion and return [p, q]

Ruby, 74 67 bytes

Probably not the golfiest solution for Ruby, but it was the first thing that sprung to mind. Increment the counter $. (starts at 0) each time rand is at least 0.5 to determine the index, then determine the list of rationals and output the right one.

-7 bytes from histocrat.

d,*a=1r
$.+=1until rand<0.5
a|=(0...d+=1).map{|n|n/d}until p *a[$.]

Try it online!

Julia, 61 49 71 bytes

p>q=rand()<.5 ? (p<q-1 ? (while gcd(p+=1,q)!=1 end;p>q) : 1>q+1) : p//q

In un-golfed form, this reads

function f(p,q)                      # define the function
    if rand() < 0.5                  # flip a coin. 
        if p < (q - 1)               # if the coin gives heads, check if p < (q - 1)
            while gcd(p+1, q) != 1   # if the greatest common divisor or p+1 and q is not equal to 1
                p = p + 1            # we increment p, because otherwise we'd get a fraction we already saw.
            end             
            return f(p, q)           # recursively call f(p, q) (where p has been incremented)
        else
            return f(1, q+1)         # if p is not greater than q - 1, we recursively call f(1, q+1)
        end
    else
        return p//q                  # If the coin flip gives tails, we return the Rational number p//q
    end
end

Now let's check to make sure we did this right:

julia> using UnicodePlots; data = [0>1 for _ in 1:10_000_000]; histogram(data)
                ┌                                        ┐ 
   [0.0 , 0.05) ┤▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇ 5000106   
   [0.05, 0.1 ) ┤ 0                                        
   [0.1 , 0.15) ┤ 1168                                     
   [0.15, 0.2 ) ┤▇ 82600                                   
   [0.2 , 0.25) ┤ 0                                        
   [0.25, 0.3 ) ┤▇▇ 313684                                 
   [0.3 , 0.35) ┤▇▇▇▇▇▇▇▇ 1250427                          
   [0.35, 0.4 ) ┤ 39124                                    
   [0.4 , 0.45) ┤ 310                                      
   [0.45, 0.5 ) ┤ 0                                        
   [0.5 , 0.55) ┤▇▇▇▇▇▇▇▇▇▇▇▇▇▇▇ 2498798                   
   [0.55, 0.6 ) ┤ 168                                      
   [0.6 , 0.65) ┤ 19403                                    
   [0.65, 0.7 ) ┤▇▇▇▇ 625484                               
   [0.7 , 0.75) ┤ 77                                       
   [0.75, 0.8 ) ┤▇ 166131                                  
   [0.8 , 0.85) ┤ 2473                                     
   [0.85, 0.9 ) ┤ 47                                       
                └                                        ┘ 
                                Frequency

julia> for frac in [(0//1), (1//2), (1//3), (2//3), (1//4), (3//4), (1//5), (2//5), (3//5), (4//5), (1//6), (5//6)]
           freq = count(x -> x == frac, data)/length(data)
           println("P($(frac)) ≈ $freq")
       end
P(0//1) ≈ 0.5000106
P(1//2) ≈ 0.2498798
P(1//3) ≈ 0.1250427
P(2//3) ≈ 0.0625484
P(1//4) ≈ 0.0313065
P(3//4) ≈ 0.0156477
P(1//5) ≈ 0.0077772
P(2//5) ≈ 0.0039116
P(3//5) ≈ 0.0019398
P(4//5) ≈ 0.0009654
P(1//6) ≈ 0.0004828
P(5//6) ≈ 0.0002473

Notice that P(5//6) ≈ 0.5 * P(1//6) as required since 2//6, 3//6, and 4//6 are reducible fractions.

Try it online!