[Question]
  [
This problem will have you analyzing circles drawn on the grid, with the gridlines drawn at integer values of \$x\$ and \$y\$.
Let \$\varepsilon\$ be a very small number (think, \$\varepsilon = 0.0001\$). If we paint a filled-in circle of radius \$\varepsilon\$ centered at the origin, it will require paint to be in \$4\$ boxes; painting a circle of radius \$1 + \varepsilon\$ will send the paintbrush to \$12\$ different boxes (second, red circle); painting a circle of radius \$\sqrt 2 + \varepsilon\$ will send the paintbrush to \$16\$ different boxes (third, green circle).
There's no way to draw a circle based at the origin that lives in exactly \$13\$ or \$14\$ or \$15\$ boxes, so the possible number of boxes that a circle can live in is given by `[4, 12, 16, ...]`.
[![Circles based at the origin containing 4, 12, 16, ... boxes](https://i.stack.imgur.com/9FzcE.png)](https://i.stack.imgur.com/9FzcE.png)
# Rules
This challenge will give you a point \$(x,y)\$ on the plane consisting of two rational numbers \$x\$ and \$y\$. Your job is to output the first sixteen (or more!) terms of the sequence consisting of all numbers \$b\$ such that there is a way to draw a circle centered at \$(x,y)\$ that contains a part of exactly \$b\$ boxes.
Any reasonable input is allowed: you can take the center point as a pair of floats, a complex number, a list containing numerators and denominators in some order, etc.
This is a [code-golf](/questions/tagged/code-golf "show questions tagged 'code-golf'") challenge, so the shortest code wins.
# Examples
In my example, the inputs are written as pairs of rational numbers. The first entry in the table is illustrated above.
```
[input]    | [output]
(x, y)     | sequence
-----------+----------------------------------------------------------------------
(0/1, 0/1) | 4, 12, 16, 24, 32, 36, 44, 52, 60, 68, 76, 80, 88, 104, 112, 120, ...
(5/1, 5/2) | 2,  6,  8, 12, 16, 20, 22, 26, 34, 38, 40, 44, 48,  52,  56,  60, ...
(1/3, 8/1) | 2,  4,  6,  8, 10, 12, 14, 18, 20, 22, 24, 26, 28,  30,  32,  34, ...
(1/2, 3/2) | 1,  5,  9, 13, 21, 25, 29, 37, 45, 49, 61, 69, 77,  81,  89,  97, ...
(1/3,-1/3) | 1,  3,  4,  6,  8,  9, 11, 13, 15, 17, 19, 21, 22,  24,  26,  29, ...
(1/9, 8/7) | 1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12, 13,  14,  15,  16, ...
```
(If it helps, the first example is given by \$4\$ times [A000592](https://oeis.org/A000592).)
  
[Answer]
# [PowerShell](https://github.com/TryItOnline/TioSetup/wiki/Powershell), 317 313 305 298 281 274 231 221 210 207 195 181 136 126 111 109 108 bytes
-an enormous 94 bytes thanks to mazzy!
I must say, I never expected to beat Charcoal in this challenge; and when I posted that 317 byte answer, I certainly never expected to get it all the way down to the low 100s. I'm hopeful I can shave those final 10 to get to sub-100, but this feels so very optimized!
Takes the x and y offsets (center of the circle) as inputs. Generates all points in the range `(-9..9,-9..9)`, then adds points for potential border intersections, which handles offset circles, and handles doubling point values on x,0 or 0,y coordinates when the circle is centered on the origin. Groups all of those points by their squared distance from the origin, then for each group outputs the sum of the group's count and the counts of all groups with distance smaller than that group. Multiplies all of the squared distances by 1e7 then casts to Integer to avoid issues with floating points.
Because it groups the distances, it does output more than the first 16 terms, but as the challenger stated that this is acceptable, the TIO link truncates the output to the first 16 terms to make it more clear which tests succeed.
```
$x,$y=$args|%{,(-9..9+($p=$_%1)|%{($t=$_-$p)*$t*1e7})}
$x|%{$n=$_;$y|%{[Int]($n+$_)}}|group|% c*|%{($s+=$_)}
```
[Try it online!](https://tio.run/##jZNNj9owEIbv/IpRZLoJSdjY@RaKhLRSpV7aquwNIQTUu2UFSRo7KivIb6djJ6Fc2OJD4vnwPO9Y47L4wyvxi@92Z/JS5xu5LfJsALiOZ3JwyHtGVtWrOA2Pjumm43Fqm6TMyHJILfSZROLeJaU1InJEedxYzYAcMEJyDEzIO27nX3K5MEluk6XVNKfXqqjL0xA2I11A2Jnyn/Gc5EI@rQQXkMHU1CKmR/1Ta5uXtVQbFfQcz7pE@KHkG8l/6ogRAGVAI2AB@Az8CIIAQgaRB1ECcQSJB0kC1MM8lcg8o63UOLdgZvhILQe/zLoFxfqQ9GAPGAMWgY8KEgg8pSBIlIgwQh3/59FHH3kJUm/zgpboaSj2kvTcQKFZAr6n@w/uwTHE@R@0RyGEFKgPjAILgaXgxxCEEKQQUYhSiGNIKCQppPG97bnqdxvodx0ilioyDYHGQFMt4dJmeg8t1ZcZfwBjGhciMG6RnqYyDQ40O@pJA2tAnnFSRTY1rauhVUOvjK/1fs2rjDatIKLVOKRH4rSPO9dy3DvbzIqLeidRk/np8hZhqpOtk@A7TAX3ZVsJiXKuis9kpQ61Bsow5sTEV/lWbHPDAcNaGI2lLVBme1DUmw0XIjM6quHy38ZFZJek27Sz@ffZUy1ksf@2fsPwYnp8UJGH7Kpf257MupJ97ckPXfkfohk03dWdPhfVfiXd59V6x89/AQ "PowerShell – Try It Online")
# Ungolfed
At this point, the ungolfed solution is more a representation of the logic behind my solution than it is an ungolfed version of the solution - as the iterations have been made, I haven't been able to keep up with the changes on the ungolfed version. I am leaving it in here, despite the discrepancy, in hopes that it will inspires to answer in their language of choice (or port it to their language)
```
#Take the x-offset, y-offset, and starting radius as a parameter
#radius starts at 0, so that the first element will be the count of points on the origin of the circle
param($xOffset,$yOffset,$radius=0)
#set the x and y offsets to be in the range [0,1)
$xOffset=$xOffset % 1
$yOffset=$yOffset % 1
#do 16 times (for the required 16 elements)
1..16|%{
	#make an array of points
	$points = @()
	#set our limits to the edge of the bounding box of the circle
	#this will normally be within 1 of the actual circle, but offsets of .5
	#will increase the search space due to .NET's banker's rounding
	$upperSearchLimitX = [Int]($xOffset+$radius+1)
	$lowerSearchLimitX = -$upperSearchLimitX 
	$upperSearchLimitY = [Int]($yOffset+$radius+1)
	$lowerSearchLimitY = -$upperSearchLimitY 
	#generate x and y for every combination of values in those ranges
	$lowerSearchLimitX..$upperSearchLimitX | %{
		$x = $_
		$lowerSearchLimitY..$upperSearchLimitY | %{
			$y = $_
			#add the point to our list
			$points+=@{x=$x;y=$y}
			#If x is zero, also add a special point that represents a border intersection to the list
			#this will also make x,0 points count for two when the center is at 0,0
			if (!$x) {
				$points+=@{x=$x+$xOffset;y=$y}
			}
			#If y is zero, also add a special point that represents a border intersection to the list
			#this will also make 0,y points count for two when the center is at 0,0
			if (!$y) {
				$points+=@{x=$x;y=$y+$yOffset}
			}
		}
	}
	#add the center point to the list
	$points += @{x=$xOffset;y=$yOffset}
	#set the minimum point as the closest point to our current circle that's outside of the bounds of the circle
	$minimum =	$points|?{ 
				($_.distance=[Math]::Sqrt(($_.x-$xOffset)*($_.x-$xOffset)+($_.y-$yOffset)*($_.y-$yOffset))) -gt $radius
			}|sort distance -t 1
	#set our new radius to be that point's distance from the x offset
	$radius = $minimum.distance
	#count the number of points that lie inside (but not on the boundary of) this new circle, output implicitly
	@($points|?{$_.distance -lt $radius}).Count
}
```
[Try it online!](https://tio.run/##jZVtb9owEMdfk09xi7I1WQiNnedW0ZAqTZq0rdPoiyGEqpSaNRUBmjgaCPLZme3EAXXQgQTx5R5@/7vYYbn4Q/LiicxmO21azic0XcxjBdhns1smeZLp2up2Oi0I7WpruciTx7QsYttQpDOWC3gPSJGRsVyIu6jXQ/72/UbpaMtFOqcFxNDXDWaWyyXJByTJJ09f0yylv5hn9GVOxy3cbJgm4vEzLvpVvHWkypHSw33p9Tmlh0dLD@GYil7viIYt8I472ooV0u756h/EkbyhzOtoa5kox2bG/c2KDfx6zeZbcUc6Bf2dtjJAZLyOM@UQ9wkHWesTWSLYlENqs9i32j9Akz1BEX0AaBO0LJ2nWZlBLOO3nzYgWLp233tMC5rMJyQefUvo0/jqavCSU517VpYsaHx8ZZvcXluSUvsPbMMA6zeF5pEK0dtikVOQNLAo24qdJoCPtlHZ6lE6fX2v90AoWLO2cmX0bhblnCrVrlI0Sgp6kxSk3tHi9PQ34sI/6XxZUr7gTrvLTo30kNWSTCh5FB7VBYQB@YBdcDA4PrgueBh8G/wQAh9CG8IQkM3ieCC21bpS1T0F071LZHTZLzZOQVl9CCXYBowB@@AwBSG4NlfghlyE5zMd/@ehS4fxQkY9zXNroi2grJdQcl2OxiE4tujfPQeHGc55oz0EHkSAHMAIsAc4AicA1wM3Ah@BH0EQQIggjCAKzm3P4pfTQKfpkGERJyMPUAAoEhLaNqNzaJEYZvAGDAucx4BBjbQFFQuwK9i@JCnsdX3HdmoR87fuftOyV7IwvpfZA8ljVNWCNKGmq0lkzI5Cc@u@J2/WkTkpSnY2YvjQ/odAX8QeVBrQnEXotcGY6khjh9d6ZidN7YJqjNXKEBZws04sysmEFEWsNgjVIi9qq6gJEj2Z8ejH4KYs6CK7fXhm7nF/c8E9F/FBc6Z5PWhKytrXP0XlPaJSqmZO28@LPEuodZc8zMjuLw "PowerShell – Try It Online")
I'm not the best at making images; if anyone can recommend a good tool for things like these, let me know. But here's a rough sketch of all the points that get checked along with the order we iterate to them in for x=3/8, y=1/4
[![First 9 iterations](https://i.stack.imgur.com/XkBKk.png)](https://i.stack.imgur.com/XkBKk.png)
[Answer]
# [Charcoal](https://github.com/somebody1234/Charcoal), 117 bytes
```
ＮθＮηＮζＮε≧×εθ≧×ηζ≧×ηε≔⁰ηＷ‹Ｌυ¹⁶«≔⟦⟧δＦ…·±ηηＦ…·±ηηＦ¬›ΣＸ⟦κλ⟧²ηＦ…·÷⁺θκε±÷⁺θκ±εＦ…·÷⁺ζλε±÷⁺ζλ±εＦ¬№δ⟦νξ⟧⊞δ⟦νξ⟧Ｆ¬№υＬδ⊞υＬδ≦⊕η»Ｉυ
```
[Try it online!](https://tio.run/##jVFLa8JAED7rr5jjLGyhCvbSU7FQhFaC7U08pGaaXUw2uo9YUvrb0wlJrUoqXpb9HvvNzM5axXZdxFldz8w2@HnI38niTtwPj7E6w9UZJsYv8fbBOZ2ahU6Vxzedk5NAEnb/ikpCdUlsYlsFb5lgtFc6I8Bnco4Pk3qFQUgY3QkBX8NBZ16uJCTsHnwUFnBm1llwuqRFbFLCOaWxJx6pSRRwrWVeeHyyxLzF15BjVOz5ttxIyLjaWLTe/ryZ8Y@61AlhxDTuJGxEM52ErlK/oRNJiCuDK27mYnBn6AluxpsWwXhMJCyNhM9Vo0XBqSPm90v/3EFCt4dEHB6ckPzmsOCmf0s5GU9Ju9DvYWQ150xjx2HsrusJjGAC4/qmzH4A "Charcoal – Try It Online") Link is to verbose version of code. Takes input as numerator denominator numerator denominator. Explanation:
```
ＮθＮηＮζＮε
```
Input the co-ordinates.
```
≧×εθ≧×ηζ≧×ηε
```
Increase the scale of the grid so that the gridlines are drawn at intervals given by the product of the denominators.
```
≔⁰η
```
Start enumerating circles with radius \$ \sqrt h + ε \$.
```
Ｗ‹Ｌυ¹⁶«
```
Repeat until 16 box counts have been found.
```
≔⟦⟧δ
```
Start with an empty list of boxes.
```
Ｆ…·±ηηＦ…·±ηη
```
Search points with a distance of up to and including \$ h \$ both horizontally and vertically.
```
Ｆ¬›ΣＸ⟦κλ⟧²η
```
Skip points whose distance from the desired point is greater than \$ h \$.
```
Ｆ…·÷⁺θκε±÷⁺θκ±εＦ…·÷⁺ζλε±÷⁺ζλ±ε
```
Consider all boxes touching this point.
```
Ｆ¬№δ⟦νξ⟧⊞δ⟦νξ⟧
```
Add the box to the list of boxes if it's not already present.
```
Ｆ¬№υＬδ⊞υＬδ
```
Add the count of boxes to the list of box counts if it's not already present.
```
≦⊕η
```
Move on to the next radius.
```
»Ｉυ
```
Print all of the box counts found.
[Answer]
# [MATL](https://github.com/lmendo/MATL), ~~51~~ 50 bytes
```
p-9:9*tJ*!+1ei1GP*1Gp\tFF5Mth3$v-!Z}J*!+u-|X<S&Y'Ys
```
The input consists of an array containing the denominators followed by an array containing the numerators. The output may contain more than 16 terms.
[Try it online!](https://tio.run/##y00syfn/v0DX0spSq8RLS1HbMNU9IFPL0L0gpsTNzdS3JMNYpUxXMaoWJFeqWxNhE6wWqR5Z/P9/tKGCYSxXtIGCQSwA) Or [verify all test cases](https://tio.run/##y00syfmf8L9A19LKUqvES0tR2zDVPSBTy9C9IKbEzc3UtyTDWKVMVzGqFiRXqlsTYROsFqkeWfzf0MxK0yWqIuR/tKGCYSxXtIGCAZA0VDACkqYKpkDSGCxuqGABJI3A4oYKxmBxYzBbFyRtqWAOVQQA) (with output trimmed to 16 terms for convenience)
### Explanation
Let (*x*, *y*) be the input values modulo 1. The code checks the distance from the (*x*, *y*) to each square with lower left corner (*a*, *b*), where *a* and *b* are integers. To do that, 9 points of the square are tested: (*a*+*m*,*b*+*n*) with *m* = 0, *x*, 1 and *n* = 0, *y*, 1, and the minimum of the 9 distances is kept.
A grid of 19×19 squares with *a*, *b* = −9, −8, ..., 9 is generated. The distance from (*x*, *y*) to each square is computed as described in the preceding paragraph. The 361 distances are sorted and run-length encoded, and the cumulative sum of the run lengths is the desired output.
To keep floating-point issues as controlled as possible, instead of using the actual values for *x*, *y*, *a*, *b*, *m*, *n*, the code uses those values multiplied by the product of the two input denominators.
```
p       % Implicit input: array with the denominators. Product
-9:9    % Push [-9, -8, ... 9]
*       % Multiply, element-wise
tJ*     % Duplicate, multiply by j (imaginary unit) element-wise
!+      % Transpose, add with broadcast. This generates the 19×19 grid of
        % complex values
1e      % Linearize into a row vector
i       % Input: array with the numerators
1GP     % Push first input again. Flip it
*       % Multiply, element-wise
1Gp     % Push first input again. Product (*)
\       % Modulo, element-wise (**)
t       % Duplicate
FF      % Push [0 0]
5Mth    % Push (*) again. Duplicate, concatenate with itself horizontally
3$v     % Concatenate three elements vertically. Gives a 3×2 matrix
-       % Subtract (*) with broadcast
!       % Transpose. The result is a 2×3 matrix
Z}      % Split into two rows with 3 elements
J*      % Multiply by j (imaginary unit) element-wise
!+      % Transpose, add with broadcast. This gives a 3×3 matrix with the
        % 9 points to be considered in each square
u       % Unique. This gives a column vector with the (unique) values
-|      % Subtract with broadcast, absolute value. Gives a 9×361 matrix
X<      % Minimum of each column. This gives the 361 distances
S       % Sort
&Y'     % Run-length encoding, producing only the run lengths
Ys      % Cumulative sum. Implicit display
```
[Answer]
# [Perl 5](https://www.perl.org/), 225 bytes
```
sub f{my($x,$y,@a,@r,%b,$n,$i)=(@_,map{[cos,sin],[-cos,sin],[cos,-sin],[-cos,-sin]}map$_/12.566,0..720);while(@r<16){$n+=2e-3;$b{int$x+$n*$$\_[0]+1e4,int$y+$n\*$$_[1]+1e4}++for@a;$b=keys%b;$i=0,push@r,$b if$b>$r[-1]&&$i++>1}@r}
```
[Try it online!](https://tio.run/##TVHbbqswEHzPV/jBbeGwSWxzTTlE/EDf@oZQFCpy6rYBBIlOEOLb012HXJBsxrO7M2NoyvbHPx@7kr2X3eH19a1uS3ZA2CVrP57te5bSKbFmGcvEUgLDLQeWMQ@YVLgCYAqxi9hF7CH2EQcCVwQsRC5CHCGWgobMlBI5y4FEfRL1l8qIYonhAIsexHFYIVaIXTLCmicuRh5i48Z8mgpuonLpouWUlOreg7CYxClM9GDgXUwUibqCNpp0vbso3XFKKskU1wpF0EvhWeFZ4dkNMRliD3GAfIDvMCRzGopWNBXeRVeUNLyJGs8psT@lDi/JjRmll9MNqI@uwSR1ysCI2vFsV7eW@W/2MGP47HuL66oBXp4aO0n5Jr7RJ@A9UVi@kvxffUi@al29wAvsppZLUXcWVYGqFpVTIwlEYst47o4F2w1XYUi3kLbwVACvgGs7sdIN7LfNkH3UHXS6yiGb3yGh@QNr8Ij9fLOUauEHAYjFIlTCjv9/6p/SStu/MrAHXjmJKuduzItBVwd@cnj1h/NNJnJHlh4Q1185abjRcfAjpVscSb7LvnsqYq4TAc2x@8TEvGB6x4s1b7O5zJ@fuXactRzTdjyffwE "Perl 5 – Try It Online")
All brawn and no brain, but the 5 given tests passes in about 11 seconds at link above. It's easy to find new tests that will fail. The function f starts at the given input (x,y) and circles around in ever wider circles counting boxes it falls into. To be (almost) sure it found all boxes it needs two rounds for a new number of boxes to be registered into the return array @r. This solution is bound to have floating point "rounding errors" and "over- or understepping" when the radius is increased. To solve this properly I'd try to sort the distances of the edges of the close by boxes, but I'm too lazy to even think about that now.
```
sub f{
  my($x,$y,@a,@r,%b,$n,$i)=
    ( @_, map {[cos,sin],[-cos,sin],
               [cos,-sin],[-cos,-sin]}
          map $_/12.566, 0..720 );     #12.566 ~ 4*π
  while(@r<16){                        #repeat until 16 elems found
    $n+=2e-3;                          #increase radius
    $b{ int$x+$n*$$\_[0]+1e4, #keys of dict b contains found boxes
 int$y+$n\*$$_[1]+1e4 }++        #+1e4 to avoid rounding up for x<0 or y<0
      for @a;                          #for all angles tested
    $b=keys%b;                         #b = number of boxes so far
    $i=0, push@r,$b                    #register b in solution if
      if !@r || $b > $r[-1] && $i++ > 1#it's > last number and if it's
                                       #found at least twice
  }
  @r                                   #return the 16 counts found
}
```
[Answer]
## Javascript, 199 bytes
My first code golf in a long time!
I originally started with a brute force solution, but couldn't get the step sizes small enough to have it produce the correct answers in a timely manner. So here's a more elegant solution that iterates points in a similar fashion to Zaelin's answer.
```
(x,y)=>{b={0:1};g=w=>b[d=parseInt(w*1e7)]=1+(b[d]|0);for(u=x%1-9;++u<9;){g(u*u);for(v=y%1-9;++v<9;)g(u*u+v*v)||u+1<9||g(v*v)}s=0;return Object.keys(b).sort((a,b)=>a-b).map(v=>(s+=b[v])).slice(0,16)}
```
[Try it online!](https://jsfiddle.net/v1L756kc/1/)
### Ungolfed
```
(x, y) => {
  b = { 0: 1 };             // initialize a map for the lengths
  g = (w) => {              // function g is used to add to the map
    d = parseInt(w * 1e7);  // avoid floating point rounding errors
    b[d] = 1 + (b[d] | 0);  // increment point with length d
  }
  for (u = (x % 1) - 9; ++u < 9; ) {    //loop over "u" coordinates
    g(u * u);                           // add horizontal line
    for (v = (y % 1) - 9; ++v < 9; )    // loop over "v" coords
      g(u * u + v * v);                 // add corner
      (u + 1 < 9) || g(v * v);          // add vertical line (but only once per above for loop)
  }
  s = 0;                        // sum of circles up to this point
  return Object.keys(b)
    .sort((a, b) => a - b)
    .map((v) => (s += b[v]))    // aggregate number of circles visited
    .slice(0, 16);        // first 16 terms only
};
```
It relies on the fact that every new vertical, horizontal and corner node found represents exactly one new grid square found. Note that the x=0 or y=0 cases do not have to be treated specially, it just means that we hit a node and a vertical/horizontal line at the same time, which represents the multiple new grid squares hit.
[Answer]
# [Vyxal](https://github.com/Vyxal/Vyxal/tree/version-2), 15 bytes
```
9:NṡvJΠ$v∆dsĠ@¦
```
[Try it Online!](https://vyxal.pythonanywhere.com/?v=2&c=1#WyIiLCLin6ggMSAzIC8gfCAxTiAzIC8g4p+pOiIsIjk6TuG5oXZKzqAkduKIhmRzxKBAwqYiLCIiLCIiXQ==)
Mostly a port of [Zaelin Goodwin's Powershell answer](https://codegolf.stackexchange.com/a/217459/100664), go look at that if you want to understand why this works. Outputs more than the first 16 terms, although beyond that not all terms may be accurate.
```
   ṡ            # Range from 
9               # 9 to
 :N             # -9
    vJ          # Append x and y to two copies of that
      Π         # Take the cartesian product
       $v∆d     # Distance between each point and the input
           s    # Sort
            Ġ@  # Take lengths of runs of identical items
              ¦ # Cumulative sum
```
]