How to easily create a polynomial function that gives a desired output?

Polynomials in the standard form (as a sum of monomials) are actually the 'wrong' sort of 'encoding' to do interpolation. This post explains how to very efficiently interpolate any polynomial function, without ever computing the polynomial coefficients! In short, we can easily express any given polynomial sequence as a Newton series by computing successive forward differences of the given sequence, and not only can we use this to efficiently compute subsequent terms in the sequence, we can also quickly convert this representation back to standard polynomial form, in the same way the linked post demonstrates for the sum of the first $n$ perfect cubes.

For your example the iterated forward differences are easily computed by hand:

 1,  2,  3,  4,  5,  6,100,659,...
 1,  1,  1,  1,  1, 94,559,...
 0,  0,  0,  0, 93,465,...
 0,  0,  0, 93,372,...
 0,  0, 93,279,...
 0, 93,186,...
93, 93,...

Which also means that you can use this forward difference technique (which works for any polynomial sequence such as $1,4,9,16,25,\cdots$) to justify that you can extend any given sequence polynomially as a reasonable pattern.

That said, if too few terms of a sequence are given (which is often the case), it is actually easier to describe it as simply a list with no pattern, because the intended pattern would have a description that is longer than just listing the few terms given.

As mentioned in a comment, the key is "Lagrange Interpolation". (See, for instance, Wikipedia or MathWorld.) The formal description can be a little opaque, but the idea is actually pretty straightforward: Build a polynomial out of parts, most of which go away with each input, then make the remaining part do what it needs to do.

For example, suppose you want the three (distinct) inputs $a$, $b$, $c$ to give the three outputs $p$, $q$, $r$. Lagrange claims that the polynomial you want has the form

$$h(x-b)(x-c) \;+\; j(x-c)(x-a) \;+\; k(x-a)(x-b) \tag{$\star$}$$

where $h$, $j$, $k$ are some as-yet-undetermined constants.

When $x=a$, two terms of $(\star)$ vanish, leaving $h(a-b)(a-c)$; likewise, $x=b$ yields $j(b-c)(b-a)$, while $x=c$ yields $k(c-a)(c-b)$. So, just choose $h$, $j$, $k$ to give the outputs you desire:

$$\begin{align} h(a-b)(a-c) = p \qquad\to\qquad h = \frac{p}{(a-b)(a-c)} \\[4pt] j(b-c)(b-a) = q \qquad\to\qquad \,j = \frac{q}{(b-c)(b-a)} \\[4pt] k(c-a)(c-b) = r \qquad\to\qquad k = \frac{r}{(c-a)(c-b)} \end{align}$$

Thus, we can re-write $(\star)$ informally as $\require{cancel}$ $$p\;\frac{\cancel{(x-a)}(x-b)(x-c)}{\cancel{(a-a)}(a-b)(a-c)} \;+\;q\;\frac{(x-a)\cancel{(x-b)}(x-c)}{(b-a)\cancel{(b-b)}(b-c)} \;+\;r\;\frac{(x-a)(x-b)\cancel{(x-c)}}{(c-a)(c-b)\cancel{(c-c)}} \tag{$\star\star$}$$

Clearly, this kind of thing can be extended to any number of inputs and outputs.