This answer uses a very similar approach to that of egreg, but it may be easier to follow for readers who are comfortable with single-variable differentiation but who aren't yet familiar with partial derivatives.
But first I'll show how we can use finite differences to verify that $(\frac12, \frac12, 2)$ gives a stationary value for $A=\frac{125}{16}$.
Recall that for small $\Delta x$,
$$f(x+\Delta x) \approx f(x) + f'(x)\Delta x$$
At a stationary point (a minimum, maximum, or point of inflexion), $f'(x)=0$.
We are given
$$A = (x^2+1)(y^2+1)(z^2+1)$$
with the constraint
$$x+y+z=3$$
So differentiating,
$$\Delta x +\Delta y +\Delta z = 0$$
For $f(x)=x^2+1, f'(x)=2x$, so
$$A+\Delta A = (x^2+1 + 2x\Delta x)(y^2+1 + 2y\Delta y)(z^2+1 + 2z\Delta z)$$
We want to show that $\Delta A=0$.
Note that when doing arithmetic with these deltas that terms like $\Delta x\Delta y$ and $(\Delta x)^2$ vanish.
Plugging in $(\frac12, \frac12, 2)$,
$$\begin{align}
A+\Delta A = & (\frac54 +\Delta x)(\frac54 +\Delta y)(5 +4\Delta z)\\
= & (5 +4\Delta x)(5 +4\Delta y)(5 + 4\Delta z)/16\\
= & (125 + 100(\Delta x+ \Delta y +\Delta z))/16\\
= & 125/16 = A\\
\end{align}$$
Thus $\Delta A$ is zero, and $A$ is indeed stationary at $(\frac12, \frac12, 2)$.
If we plug in $(1, 1, 1)$, we get
$$8+\Delta A = 8 + 2(\Delta x+ \Delta y +\Delta z)$$
i.e., $\Delta A = 0$, so it's also a stationary point. In fact, it's a local maximum.
Now I'll show how to find these points. Let's temporarily treat one of the variables, eg $z$, as if it were a constant value. Geometrically, we want to know how the A function behaves on the $z=k$ plane as $x$ and $y$ vary. So we look at the function
$$B=(x^2+1)(y^2+1)$$
with $A=B(k^2+1)$.
If $A$ is stationary at $(x, y, z=k)$ then $B$ must also be stationary at $(x, y)$.
At $z=k$, we have
$$y=(3-k)-x$$
so $y$ and hence $B$ are simple functions of $x$, and we can use plain single-variable calculus on them. Firstly,
$$\frac{dy}{dx}=-1$$
Using the product and chain rules,
$$\frac{dB}{dx} = 2x(y^2+1) + 2y(x^2+1)\frac{dy}{dx}$$
$$\frac{dB}{dx} = 2x(y^2+1) - 2y(x^2+1)$$
At a stationary point, $$\frac{dB}{dx}=0$$
so
$$x(y^2+1) = y(x^2+1)$$
and for $x,y \ne 0$,
$$(y+1/y) = (x+1/x)$$
This is a quadratic, with solutions $x=y$ and $x=1/y$. We can see that they're the only solutions by factoring:
$$x(y^2+1) = y(x^2+1)$$
$$xy^2+x - yx^2-y = 0$$
$$xy(y-x) + x-y = 0$$
$$(xy-1)(y-x) = 0$$
So we've established that $B$ is stationary when $x=y$ or $xy=1$. That's true for any given $z=k$, regardless of whether A is also stationary at that point.
By symmetry, we can apply the same argument to a fixed $x$ or $y$ plane, and so we get a stationary point on $A$ when a point $(x,y,z)$ is stationary on all 3 planes simultaneously.
So now we just need to combine each of those solutions:
$x=y$ or $xy=1$,
$x=z$ or $xz=1$,
$y=z$ or $yz=1$
with $x+y+z=3$ to find the stationary points of A. We get either
$$x=y=z=1$$
or
$xy=1$ and $y=z$, and its symmetrical variations.
$$x+y+z=3$$
$$xy+2y^2=3y$$
$$2y^2-3y+1=0$$
$$(y-1)(2y-1)=0$$
So either
$$y=z=1,x=1$$
or
$$y=z=\frac12, x=2$$
and the symmetrical variations thereof.
Thus we end up with these four stationary points:
$(1,1,1),$
$(2,\frac12,\frac12),$
$(\frac12,2,\frac12),$
$(\frac12,\frac12,2)$
It would be nice to look at the graph of this function. Of course, that's a little difficult since we aren't four-dimensional. ;) But we can use our constraint $$x+y+z=3$$ to turn it into a viewable 3D function. We eliminate $z$, and use the 3rd axis for $A$.
Here's a graph of
$$A = (x^2+1)(y^2+1)((3-x-y)^2+1)$$
in the region $-0.5 \le x, y \le 3.5$. The function $A$ grows fairly quickly outside that region, and I've sliced off the boring bits above $A=10$ to focus on the interesting region containing the stationary points, which I've marked with coloured dots, red for the maximum, and blue for the minima.

Click on the image for an interactive version that you can zoom, rotate, etc.