How do I plot Thomae's function in Mathematica?

Question

I wanted to plot this function

$$f(x) =\begin{cases} 1 & \text{if } x= 0 \\ \tfrac1{q} & \text{if } x = \tfrac{p}{q}\\ 0 & \text{if } x \in \mathbb{R}-\mathbb{Q} \end{cases}$$

so I wrote

FuncThomae[x_] := If[ExactNumberQ[Rationalize[x]], If[x == 0, a = 1, 
                  L = #^-1 & /@ Divisors[Numerator[Rationalize[x]]]], a = 0]

and

ListDomain[xmin_, xmax_] := Table[Outer[List, {x}, FuncThomae[x]],
                            {x,xmin,xmax,0.001}] // Flatten[#, 1] &

My result doesn't take all the real numbers (nor negatives) in its domain, and for $-1$ to $1$, it should have looked like so:

but my function does not cater to negatives, nor does it look like the above plot. It looks like this from 2 to 5:

Nearly similar, but not quite. Can someone help me to perfect the function?

Wanted to duplicate this

Answer 1

I'd suggest producing a list of rational numbers and then plot the function there, like so:

maxq = 100;
fracs = Table[p/q, {q, 2, maxq}, {p, 2, q}] // Flatten // DeleteDuplicates;
pq = {#, 1/Denominator @ #} & /@ fracs;

ListPlot[pq, PlotRange -> {0, 1}]

Answer 2

Another possibility that avoids the generation of fractions not in lowest terms (and thus the use of Union[] or DeleteDuplicates[]) rests on generating a Farey sequence, and then applying the Dirichlet-Thomae function to that:

farey[n_Integer?Positive] := Module[{v = 0, w = 1, p = 1, q = n, t},
  Join[{0, 1/n}, Flatten[Last[
     Reap[While[p < q,
       t = Quotient[n + w, q];
       {{p, v}, {q, w}} = {{t p - v, p}, {t q - w, q}};
       Sow[p/q]]]
     ]]]]

ListPlot[{#, 1/Denominator[#]} & /@ farey[100], Frame -> True, PlotRange -> {0, 1}]

Answer 3

Why not use DiscretePlot directly?

DiscretePlot[1/Denominator[Rationalize[x]], {x, -1, 1, 1/(2*3*4*5*6*7)}, 
       PlotRange -> {0, 1}, Filling -> None, RegionFunction -> Function[{x, y}, Abs[x] != 1]]

The RegionFunction throws out the cases where $\frac{p}{q}=1$.