Maxwell stress tensor on two steady current carrying wires I am having trouble thinking about the force on a wire via maxwell stress tensor. The problem is as follows; given two  parallel infinite wires separated by a distance of 2a with steady currents in the same direction, find the force per unit length on the two wires. I have calculated the maxwell stress tensor given the field of a single wire, in SI units it's something like $T_{ij}=\frac{I^2\mu_0}{4\pi^2r^2} \begin{pmatrix} \sin^2{\phi}-\frac{1}{2} & -\cos{\phi}\sin{\phi} & 0\\
-\cos{\phi}\sin{\phi} & \cos^2{\phi}-\frac{1}{2} & 0 \\
0 & 0 & \frac{-1}{2}
\end{pmatrix}$
Now if I understand correctly, the trick to making use of the tensor here is to enclose one of the currents by a surface where you let some part of it go off to infinity where the field strength is zero so you are only left with integrating over some nice surface (does it have to be compact?). My real question is where do I start with parameterizing this surface? I thought about making a kind of dome enclose one of the wires and having the non zero curvature part of the dome go out to infinity and therefore im left with integrating over the plane in between and orthogonal to both the wires. Is this the right approach? Is my understanding ok? Do you agree with my calculation thus far?
I would really appreciate some insight and I'm sorry if I'm asking to be spoon fed here in advance.
 A: Good question with a subtle answer.
If you read the argument for the Maxwell stress tensor (I linked Wikipedia but I recommend a book like Griffiths to really understand it), the motivation for it is to have a single object that encodes all the forces in space. The stress tensor works by the equation
$$
\nabla \cdot T = f + \epsilon_{0}\mu_{0}\frac{\partial\vec{S}}{\partial t}
$$
where I am using $T$ to maintain consistency with you (the wiki I linked uses $\sigma$), $f$ is the body force in space (force per unit volume), and $\vec{S}$ is the Poynting vector. The usual method in electrostatics and magnetostatics that you are trying to replicate is to calculate the field at a particle or body from all other charged or current carrying particles or bodies, and then apply that field to get the force or differential force. In essence, the particle or body "does not push on itself." The Maxwell stress tensor, on the other hand, calculates the force from the total field in space, that is, from all the charged and current carrying particles and bodies.
This may seem wrong, like a particle is "pushing on itself," but it is not. You can see this in the argument of its derivation (e.g. from the wiki but also elsewhere):

the force per unit volume is
$$\vec{f} = \rho E + \vec{J} \times \vec{B}$$
Next, $\rho$ and $\vec{J}$ can be replaced by the fields $\vec{E}$ and $\vec{B}$ using Gauss's law and Ampère's circuital law:

What this means is that the charges and currents used to calculate these forces are calculated from the fields, including the fields of the charges and currents upon which the action is happening. If $\vec{F} = q\vec{E}$, the $q$ itself is a function of  $\vec{E}$, so we substitute it and remove $q$ from consideration (and likewise for $\vec{J}$).
If you understood all that, you may be asking "well, can I just include the second wire in my stress tensor then"? The answer is a tentative "yes." That is what you have to do in order to get what you are tying to show, but not directly as you have set this up. The magnetic field of a "wire" (the theoretical wire in these equations) goes to infinity at the wire (at an inverse). This is like the electric field going to infinity as well (at an inverse square). This is not an error in the stress tensor method and is to be expected, as the Maxwell stress tensor is intended to give you force per volume and your infinitely thin wire also suffers from having no volume.
So, to do what you are trying to do, you have to set up two "real" wires (with a non zero radius) and calculate the correct (total) Maxwell stress tensor. Then, you have to take the divergence of that tensor (which will give you a vector). Since the Poynting vector is constant in this set up (it is actually zero), this divergence should equal the force per volume on the wire, but only at the wire which you will have defined based on the assumed radii and locations (centerlines). You should then be able to integrate over the cross sectional area and along the length of one wire and get a result close to the usual "force per length" equation for two current carrying wires (but not exact yet). If you get all that right, you ought to be able to take a limit as the two radii go to zero, and that should then approach the equation you are know and are trying to get.
If I am feeling motivated, I may try to do this and post it back here for you.
A: In a suitably chosen Cartesian coordinate system, the current density (assuming equal currents $I$ in both wires) is given by $$\vec{j}(x,y,z)= I \delta(y)[\delta(x+a)+\delta(x-a)]\vec{e}_z,$$ generating the magnetic field $$\vec{B}(x,y,z)= \frac{\mu_0 I}{2 \pi} \left[\frac{-y \, \vec{e}_x +(x+a) \, \vec{e}_y}{(x+a)^2 +y^2}+\frac{-y \, \vec{e}_x+(x-a) \, \vec{e}_y}{(x-a)^2+y^2} \right].$$
The force $\vec{F}$ acting on the segment $S=\left\{(-a,0,z)|0\le z \le \ell\right\}$ of the left wire is obtained from the Maxwell stress tensor $$T_{ik}= \frac{1}{\mu_0} \left(B_i B_k- \frac{1}{2}\delta_{ik}\,  \vec{B}\!\cdot\! \vec{B} \right)$$ by the surface integral $$F_k= \int\limits_{\partial V} \!d\sigma_i \, T_{i k}, $$ where $\partial V$ is the boundary of an arbitrary volume $V$ containing the segment $S$ but neither any other part of the left wire nor any piece of the right one. As already remarked in a comment by Jan Lalinsky, a convenient choice of the domain $V$ is given by $V=\{(x,y,z)|x^2+y^2\le R^2, \, x\le 0, \, 0\le z \le \ell\}$ in the limit $R\to \infty$. In the surface integral, only the contribution from the slice of the $y$-$z$-plane with $0 \le z \le \ell$ survives. Using $d \vec{\sigma} = dy \, dz \, \vec{e}_x$, we obtain $$F_k  =\frac{1}{\mu_0} \int\limits_0^\ell \! dz \int\limits_{-\infty}^\infty\! dy \,  \left[(B_x(0,y,z) B_k(0,y,z) - \frac{\delta_{x k}}{2} \, \vec{B}(0,y,z) \cdot \vec{B}(0,y,z) \right]. $$ Inserting $B_x(0,y,z)=-\frac{\mu_0 I y}{ \pi (a^2+y^2)}$ and $B_y(0,y,z)=0$, $B_z(x,y,z)=0$, one finds $$F_x= \frac{\mu_0I^2 \ell}{2 \pi^2 a} \int\limits_{-\infty}^\infty \!\! d \eta \, \frac{\eta^2}{(1+\eta^2)^2}= \frac{\mu_0 I^2 \ell }{4  \pi a}, \quad F_y=F_z=0,$$ as to be expected.
