What is the exchange interaction?
In quantum chemistry, probably the easiest way to understand the "exchange interaction" is within the context of the Hartree-Fock model.
$
\newcommand{\op}{\hat}
\newcommand{\el}{_\mathrm{e}}
\newcommand{\elel}{_\mathrm{ee}}
\newcommand{\elnuc}{_{\mathrm{en}}}
\newcommand{\core}{^{\mathrm{core}}}
\newcommand{\bracket}[3]{\langle{#1}\vert{#2}\vert{#3}\rangle}
$



To reduce the complexity of the electronic Schrödinger equation
$$
    \op{H}\el \psi\el(\vec{q}\el) = E\el \psi\el(\vec{q}\el) \, ,
$$
where
$$
    \op{H}\el
    =
    \op{T}\el + \op{V}\elnuc + \op{V}\elel
    =
    - \sum\limits_{i=1}^{n} \frac{1}{2} \nabla_{i}^{2}
    - \sum\limits_{\alpha=1}^{\nu} \sum\limits_{i=1}^{n} \frac{Z_{\alpha}}{r_{\alpha i}}
    + \sum\limits_{i=1}^{n} \sum\limits_{j &gt; i}^{n} \frac{1}{r_{ij}} \, ,
$$
we would like to separate the electronic coordinates from each other by writing down the many-electron wave function as a product of one-electron ones. Unfortunately, such separation won't work due to the presence of $\op{V}\elel$ term in the electronic Hamiltonian. But if instead of $\op{V}\elel$ potential, which prevents the separation of electronic coordinates, a model potential of the form $\sum\nolimits_{i=1}^{n} v_{\mathrm{MF}}(\vec{r}_{i})$ had entered the electronic Schrödinger equation, it would be reduced to a set of $n$ one-electron Schrödinger equations with the many-electron wave function being just a simple product of their solutions, one-electron wave functions $\psi_{i}(\vec{r}_{i})$.



More precisely, since we have to take spin of electrons and the principle of antisymmetry of the electronic wave function into account, the electronic wave function would be an antisymmetric product of one-electron wave functions $\psi_{i}(\vec{q}_{i})$ termed the Slater determinant,
$$
    \Phi
    =
    \frac{1}{\sqrt{n!}}
    \begin{vmatrix} 
        \psi_{1}(\vec{q}_{1})  &amp;  \psi_{2}(\vec{q}_{1})  &amp;  \cdots  &amp;  \psi_{n}(\vec{q}_{1})  \\
        \psi_{1}(\vec{q}_{2})  &amp;  \psi_{2}(\vec{q}_{2})  &amp;  \cdots  &amp;  \psi_{n}(\vec{q}_{2})  \\
        \vdots                 &amp;  \vdots                 &amp;  \ddots  &amp;  \vdots                 \\
        \psi_{1}(\vec{q}_{n})  &amp;  \psi_{2}(\vec{q}_{n})  &amp;  \cdots  &amp;  \psi_{n}(\vec{q}_{n})
    \end{vmatrix} \, ,
$$
where the one-electron wave functions $\psi_{i}$ of the joint spin-spatial coordinates of electrons $\vec{q}_{i} = \{ \vec{r}_{i}, m_{si} \}$ are referred to as spin orbitals. Physically the model potential of the above mentioned form represent the case when electrons do not instantaneously interact with each other, but rather each and every electron interacts with the average, or mean, electric field created by all other electrons and given by $v_{\mathrm{MF}}(\vec{r}_{i})$.

Minimization of the electronic energy functional
$$
E\el^{\mathrm{HF}} = \bracket{ \Phi }{ \op{H}\el }{ \Phi } = \bracket{ \Phi }{ \op{T}\el }{ \Phi } + \bracket{ \Phi }{ \op{V}\elnuc }{ \Phi } + \bracket{ \Phi }{ \op{V}\elel }{ \Phi }
$$
with respect to variations in spin orbitals subject to constraint that spin orbitals remain orthonormal eventually results in the set of the so-called canonical Hartree-Fock equations which define canonical spin orbitals $\psi_i$ together with the corresponding orbital energies $\varepsilon_i$
$$
    \op{F} \psi_{i}(\vec{q}_{1})
    =
    ε_{i} \psi_{i}(\vec{q}_{1}) \, ,
    \quad
    i = 1, \dotsc, n \, ,
$$
where $\op{F} = \op{H}\core + \sum\nolimits_{j=1}^{n} \big(\op{J}_{j} - \op{K}_{j} \big)$ is the Fock operator with
$$
    \op{H}\core \psi_{i}(\vec{q}_1)
    =
    \Big(
        - \frac{1}{2} \nabla_{1}^2
        - \sum\limits_{\alpha=1}^{\nu} \frac{Z_{\alpha}}{r_{\alpha 1}}
    \Big)
    \psi_{i}(\vec{q}_1) \, , \\
    \op{J}_{j} \psi_{i}(\vec{q}_1) = \bracket{ \psi_{j}(\vec{q}_2) }{ r_{12}^{-1} }{ \psi_{j}(\vec{q}_2) } \psi_{i}(\vec{q}_1) \, , \\
    \op{K}_{j} \psi_{i}(\vec{q}_1) = \bracket{ \psi_{j}(\vec{q}_2) }{ r_{12}^{-1} }{ \psi_{i}(\vec{q}_2) } \psi_{j}(\vec{q}_1) \, .
$$



Here, the first part of the Fock operator, $\op{H}\core$, is known as the one-electron core Hamiltonian and it is simply the Hamiltonian operator for a system containing the same number of nuclei, but only one electron. The second part of the Fock operator, namely, $\sum\nolimits_{j=1}^{n} \big(\op{J}_{j} - \op{K}_{j} \big)$, plays the role of the mean-field potential $v_{\mathrm{MF}}$ that approximates the true potential $\op{V}\elel$ of interactions between the electrons. The first operator $\op{J}_{j}$ rewritten in the following form
$$
    \op{J}_{j}(\vec{q}_{1})
    =
    \int \psi_{j}^*(\vec{q}_{2}) r_{12}^{-1} \psi_{j}(\vec{q}_{2}) \,\mathrm{d} {\vec{q}_{2}}
    =
    \int \frac{ |\psi_{j}(\vec{q}_2)|^{2} }{ r_{12} } \mathrm{d} \vec{q}_2 \, ,
$$
can be clearly interpreted as the the Coulomb potential for electron-one at a particular point $\vec{r}_{1}$ in an electric field created by electron-two distributed over the space with the probability density $|ψ_{j}(\vec{q}_2)|^{2}$, and for this reason it is called the Coulomb operator.
The second operator $\op{K}_{j}$ has no simple physical interpretation, but it can be shown to arise entirely due to anti-symmetry requirement, i.e. if instead of a Slater determinant one uses a simple product of spin orbitals, termed the Hartree product, there will be no $\op{K}_{j}$ terms in the resulting equations (the so-called Hartree equations). It is for that reason that $\op{K}_{j}$ is called the exchange operator.



To quickly understand why exchange terms appear when using a Slater determinant instead of a simple product of spin orbitals look no further than at the Slater rules. For $\op{V}\elel$ which is a two-electron operator we have

$$
    \bracket{ \Phi }{ \op{V}\elel }{ \Phi }
    =
    \sum\limits_{i=1}^{n} \sum\limits_{j&gt;i}^{n}
    \Big(
        \bracket{ \psi_{i}(1) \psi_{j}(2) }{ r_{12}^{-1} }{ \psi_{i}(1) \psi_{j}(2) }
        -
        \bracket{ \psi_{i}(1) \psi_{j}(2) }{ r_{12}^{-1} }{ \psi_{j}(1) \psi_{i}(2) }
    \Big) \, ,
$$
where the second (exchange) part would be absent if $\Phi$ were a simple product of spin orbitals, rather than a Slater determinant.1)



To the question on why do we interpret exchange terms in the Hartree-Fock as a manifestation of some sort of interaction. On the one hand, it must be clear from exchange terms $\op{K}_{j}$ being part of potential energy terms $v_{\mathrm{MF}}$, that they are indeed related to some type of interaction between electrons. However, it must be said that, unlike for the four fundamental interactions, there exists no force based on exchange interaction. Strictly speaking, I don't even think that the term "interaction" has a definite meaning in physics, unless it simply stands for the "fundamental interaction" in which case it is merely a synonym for "fundamental force".

Exchange interaction is not one of the four fundamental interactions, so the meaning of the word "interaction" here is a bit different than that for fundamental interactions. Oxford Dictionary defines interaction as follows,


  1.1 Physics A particular way in which matter, fields, and atomic and subatomic particles affect one another, e.g. through gravitation or electromagnetism.


And exchange interaction indeed is a way by which electrons (of the same spin) "affect one another". This can be readily seen as follows. Consider the simplest many-electron system, a two-electron one, and let us examine the following two important events:


$r_1$ - the event of finding electron-one at a point $\vec{r}_{1}$;
$r_2$ - the event of finding electron-one at a point $\vec{r}_{2}$.


Probability theory reminds us that in the general case the so-called joint probability of two events, say, $\vec{r}_{1}$ and $\vec{r}_{2}$, i.e. the probability of finding electron-one at point $\vec{r}_{1}$ and at the same time electron-two at point $\vec{r}_{2}$, is given by
$$
    \Pr(r_1 \cap r_2)
    = 
    \Pr(r_1\,|\,r_2) \Pr(r_2)
    = 
    \Pr(r_1) \Pr(r_2\,|\,r_1) \, ,
$$
where 


$\Pr(r_1)$ is the probability of finding electron-one at point $\vec{r}_{1}$ irrespective of the position of electron-two;
$\Pr(r_2)$ is the probability of finding electron-two at point $\vec{r}_{2}$ irrespective of the position of electron-one;
$\Pr(r_1\,|\,r_2)$ is the probability of finding electron-one at point $\vec{r}_{1}$, given that electron-two is at $\vec{r}_{2}$;
$\Pr(r_2\,|\,r_1)$ is the probability of finding electron-two at point $\vec{r}_{2}$, given that electron-one is at $\vec{r}_{1}$.


The first two of the above probabilities are referred to as unconditional probabilities, while the last two are referred to as conditional probabilities, and in general $\Pr(A\,|\,B) \neq \Pr(A)$, unless events $A$ and $B$ are independent of each other.

If the above mentioned events $\vec{r}_{1}$ and $\vec{r}_{2}$ were independent, then the conditional probabilities would be equal to their unconditional counterparts 
$$
    \Pr(r_1\,|\,r_2) = \Pr(r_1) \, , \quad
    \Pr(r_2\,|\,r_1) = \Pr(r_2) \, ,
$$
and the joint probability of $\vec{r}_{1}$ and $\vec{r}_{2}$ would simply be equal to the product of unconditional probabilities,
$$
    \Pr(r_1 \cap r_2) = \Pr(r_1) \Pr(r_2) \, .
$$
In reality, however, the events $\vec{r}_{1}$ and $\vec{r}_{2}$ are not independent of each other, because electrons "affect one another",
$$
    \Pr(r_1\,|\,r_2) \neq \Pr(r_1) \, , \quad
    \Pr(r_2\,|\,r_1) \neq \Pr(r_2) \, ,
$$
and consequently the joint probability of $\vec{r}_{1}$ and $\vec{r}_{2}$ is not equal to the product of unconditional probabilities,
$$
    \Pr(r_1 \cap r_2) \neq \Pr(r_1) \Pr(r_2) \, .
$$

Inequalities above hold true for two reasons, i.e. there are two ways electron "affect one another".

First, electrons repel each other by Coulomb forces, and, as a consequence, at small distances between the electrons $\Pr(r_1\,|\,r_2) &lt; \Pr(r_1)$ and $\Pr(r_2\,|\,r_1) &lt; \Pr(r_2)$, while at large distances $\Pr(r_1\,|\,r_2) &gt; \Pr(r_1)$ and $\Pr(r_2\,|\,r_1) &gt; \Pr(r_2)$.
In the extreme case when $\vec{r}_{1} = \vec{r}_{2}$, the Coulomb repulsion between the electrons becomes infinite, and thus, $\Pr(\vec{r}_{1}\,|\,\vec{r}_{1}) = 0$ and consequently $\Pr(\vec{r}_{1} \cap \vec{r}_{1}) = 0$, i.e. the probability of finding two electrons at the same point in space is zero.

Secondly, as a consequence of the Pauli exclusion principle, electrons in the same spin state can not be found at the same location in space, so that for electrons in the same spin state there is an additional contribution into the inequalities between conditional and unconditional probabilities above.
This effect is relatively localized as compared to one due to Coulomb repulsion, but bearing in mind the relationship between probabilities and wave functions and taking into account the continuity of the latter, it is still noticeable when electrons are close to each other and not just at the very same location in space.

Now, for a Hartree product wave function
$$
  \psi_{\mathrm{HP}}(\vec{q}_{1}, \vec{q}_{2}) = \psi_{1}(\vec{q}_{1}) \psi_{2}(\vec{q}_{2}) \, ,
$$
one can show that $\Pr(r_1 \cap r_2) = \Pr(r_1) \Pr(r_2)$ if two electrons are in different spin states as well as when they are in the same spin state2). But for a Slater determinant
$$
    \Phi(\vec{q}_{1}, \vec{q}_{2}) = \frac{1}{\sqrt{2}}
    \Big(
    \psi_{1}(\vec{q}_{1}) \psi_{2}(\vec{q}_{2}) - 
    \psi_{1}(\vec{q}_{2}) \psi_{2}(\vec{q}_{1}) 
    \Big) \, ,
$$
$\Pr(r_1 \cap r_2) = \Pr(r_1) \Pr(r_2)$ will hold only for electrons of unlike spin, while when both electrons are in the same spin state $\Pr(r_1 \cap r_2) \neq \Pr(r_1) \Pr(r_2)$2). And this inequality is a clear indication of exchange interaction.



1) I leave it to OP as an exercise. Derivation for a simple product of spin orbitals is pretty trivial, while for the case of a Slater determinant OP can consult, for instance, Szabo &amp; Ostlund, Modern Quantum Chemistry, Section 2.3.3.

2) I leave it as another exercise.