The reason for your confusion, and the answer to question 4, is that the first equation you give (although you copied it missing crucial primes on the indices or on the variables that serve to distinguish between the x and x' coordinate labels, as Lubos Motl pointed out) is the geodesic equation for a flat space-time in nonlinear coordinates, while the second form of the equation tells you the geodesics in a curved space-time with an arbitrary metric.
These are not a 4 by 4 system of equations, but a system of only 4 equations, since there is only one free (not summed over) index. Please do not be freaked out by the Einstein convention--- it is simple to add the summations explicitly as you did, and if you know how to do this correctly, you are not confused, since there is nothing more to it. The Einstein convention is the only simple way of keeping track of covariant/contravariant quantities, tensor products, and the covariant linear maps allowed between them.
It is good to begin in two dimensions, as you did. The explicit equations are (writing out the sum in full):
$$ \ddot{x}^0 + \Gamma^0_{00} \dot{x}^0 \dot{x}^0 + \Gamma^0_{01} \dot{x}^0 \dot{x}^1 + \Gamma^0_{10} \dot{x}^1\dot{x}^0 + \Gamma^0_{11} \dot{x}^1\dot{x}^1 = F^0 $$
$$ \ddot{x}^1 + \Gamma^1_{00} \dot{x}^0 \dot{x}^0 + \Gamma^1_{01} \dot{x}^0 \dot{x}^1 + \Gamma^1_{10} \dot{x}^1\dot{x}^0 + \Gamma^1_{11} \dot{x}^1\dot{x}^1 = F^1 $$
Where $F^0, F^1$ are an additional non-gravitational force, which I will omit (because it doesn't matter for the gravitational force, you understand that separately--- this extra force is in the equation in the nonhomogenous case).
These are two equations, not a 2 by 2 system of equations (question 5), and in 4 dimensions they are 4 equations, not a 4 by 4 system (although they can be a "4 by 4 system" if you mean that you need to solve 4 linear equations to find the accelerations in a new coordinate, see below). Also, note that the cross terms appear, they are not zero, and this answers question 2, although also note that the two cross terms are equal (because the $\Gamma$ is symmetric in the lower two indices for the usual metric connection), so that you can add the middle two terms in both equations, but you don't want to do this, so that the explicit Einstein contractions in the equation are made obvious.
The answer to question 3 is yes, you can use this convention, with two caveats:
* make sure you know what the derivative means in all cases, in the case you give first, it is differentiating one coordinate with respect to another coordinate (see below).
* Make sure that you don't confuse the quantity ${\partial y^\mu \over \partial x^\nu} = \partial_\nu y^\mu = y^\mu_{,\nu} $ with the inverse quantity ${\partial x^\mu \over \partial y^\nu} = M^\mu_\nu $ where $M^\mu_\nu$ is the inverse matrix to $y^\mu_{,\nu}$.
It is easy to lose sight of the meaning of the equations in the soup of symbols and formal relations, and this is no good. To keep the meaning of everything straight, you need a few examples in the back of your head, and these are useful anyway.
Explicit example
Consider 2d spacetime in the regular x,t rectangular coordinantes, and define the new coordinates
$$t'=t$$
$$ x'= x-t^2$$
These coordinates are useless, the center of the x'-coordinate is accelerating to the right with a constant Galilean acceleration, but they are sufficiently not special to explain equation 1. The inverse map is simple
$$ t = t' $$
$$ x = x' + t'^2 $$
In the original rectangular system of coordinates, the geodesic equation is
$$\ddot{x} = 0$$
$$\ddot{t} = 0$$
You can find the new equations just by substitution (but be sure to understand what this means--- you are substituting the x,t coordinate values as a function of the proper $\tau$ in terms of the x',t' coordinates of the same trajectory).
$${d^2\over d\tau^2} (x' + t'^2) = \ddot{x}' + 2t' \ddot{t}' + 2 \dot{t}'^2 = 0$$
$$\ddot{t'} = 0$$
So you see that to solve for the t' and x' accelerations, you need to invert the coefficients of the second derivative terms, which form the Jacobian matrix of the coordinate transformation, and the coefficients of the first derivative squared terms are the second derivatives of the coordinate transformation. The result is silly n this case:
$$ \ddot{x}' = - 2 \dot{t}'^2 $$
$$ \ddot{t}' = 0 $$
Which solves to
$$ t = a \tau $$
$$ x= b \tau - a^2 \tau^2 $$
which is as expected, it is the transformation of a straight line.
It is easy to work out the general case from this.
General Linear/quadratic transformations
If you look near a point x,t, you can define the general quadratic change of variables
$$ x' = ax + bt + P x^2 + 2Qxt + R t^2$$
$$ t' = cx + dt + S x^2 + 2Txt + U t^2 $$
With constants a,b,c,d P,Q,R,S,T,U. You should study this, because it is the general case! Any coordinate transformation is linear plus quadratic to second order near a point, and the cubic terms are not important for determining the transformation law of the connection coefficients $\Gamma$.
You write this more formally in symbol-soup as
$$ x^{\mu'} = {\partial x^{\mu'} \over \partial x^\mu} x^\mu + {1\over 2} {\partial ^2 x^{\mu'} \over \partial x^{\alpha} \partial x^{\beta} } x^{\alpha} x^{\beta} $$
Where the quantities ${\partial x^{\mu'} \over \partial x^\mu}$ are just the (a,b;c,d) coefficients, and the second derivative quantities are $A,C,B,P,Q,R$. Remember that these are just some numbers at every point. This thing is valid to second order near x=0 (which maps to x'=0 under the coordinate change--- this can be arranged just by additively shifting the x' coordinate).
You can formally write the inverse map too (remember that here and above, primed indices or primed coordinate refer to the new coordinates, while unprimed indices or unprimed coordinates are the old ones).
$$ x^{\mu} = {\partial x^{\mu} \over \partial x^{\nu'} } x^{\nu'} + {1\over 2} {\partial ^2 x^{\mu} \over \partial x^{\nu'} \partial x^{\sigma'}} x^{\nu'} x^{\sigma'} $$
Then you can solve for the new form of the geodesic equation by taking the second $\tau$ derivative along a curve $x^{\mu}(\tau)$, setting it to zero, and inverting the linear coefficients multiplying the second derivatives to solve for the second derivatives, as above (except easier, because you only need to do it near the one point x=0).
$$ {\partial x^\mu \over \partial x^{\mu'} } \ddot{x}^{\mu'} + {\partial^2 {x^\mu}\over \partial x^{\alpha'} \partial x^{\beta'} } \dot{x}^{\alpha'} \dot{x}^{\beta'}$$
If you solve this "system of equations" (I don't like this distinction--- 4 equations are 4 equations, whether the leading linear coefficients are explicitly diagonal or not) for the accelerations using the formal inverse matrix to ${\partial x^\mu \over \partial x^{\mu'}}$, you get the form of the Christoffel symbols in the new coordinates, as your book has it.
The meaning of these formal equations is only fully clarified with appropriate examples, which you should provide for yourself.