The explanation which I got for mesh current analysis was that it is a way to express branch currents in terms of other variables called "mesh currents" but that was basically it, i.e.:
$$\boldsymbol{i} = \boldsymbol{A}\begin{bmatrix}i_{m1} \\ \vdots \\ i_{mn}\end{bmatrix}$$
The only other thing mentioned was that the mesh currents are not physically measurable. But then what is their meaning and how to actually find the equations which give the branch currents in terms of them?
The explanation you were given was a bit too abstract.
Mesh currents are hypothetical currents circulating around the independent loops of a circuit. They aren’t necessarily real currents in any single branch but are a bookkeeping tool: each branch current is expressed as the algebraic sum of the mesh currents that pass through it. This lets you describe all branch currents in terms of fewer unknowns.
To find the equations for the mesh currents, assign a mesh current to each independent loop, then apply Kirchhoff’s Voltage Law (KVL) around each loop. For branches shared between loops, the branch current is the sum or difference of the mesh currents depending on their directions. Solving the resulting system of linear equations gives the mesh currents, which can then be combined to compute the actual branch currents.
You've selected your answer. But I thought it may be useful to answer your question in a way I've not yet seen here at EESE. It's at a slightly higher level of view than the more "machine code view" that is often given.
(In the following, some of the diagrams were developed using Desmos.)
The only other thing mentioned was that the mesh currents are not physically measurable. But then what is their meaning and how to actually find the equations which give the branch currents in terms of them?
I'll address what I see as three questions:
They are all related.
Truly understanding what a loop is and how they may be uncovered requires a few mathematical foundations. So I'll work to develop that below.
I use mesh all the time when I really mean loop in the academic sense. So I guess the first thing to get right is the terminology.
From a Khan Academy lesson:
Nodal analysis, and its modified form aka MNA, are quite general. But if \$m\$ is taken as the number of edges (2-terminal resistors, for example) in a circuit, then the loops can be readily found by treating each circuit edge as a basis vector. These basis vectors will describe a \$p\$-subspace within the \$m\$-space. Those vectors pass through the origin and that are perpendicular to the \$p\$-subspace identify the loops.
The reason is simple. These are the only cases where all of the nodes between the edges have a current sum of zero. If there were any "leaks" out of, or into, a node other than from the edges involved then it can't be a loop. A loop is only and exactly the case where the current entering the node from one edge is exactly the same as the current leaving the node via another edge, and where all of the edges are ultimately connected.
This may sound obscure. And I apologize. I'll need to set down some of the foundations to clarify.
Let's take the simplest non-trivial case and cast it three ways:
simulate this circuit – Schematic created using CircuitLab
Assume that \$-1\$ means from and that \$+1\$ means to; with \$0\$ meaning not connected, then we can form these incidence matrices to describe the above three cases:
$$\begin{align*} \begin{smallmatrix} \fbox{$A$}&N_1&N_2&N_3\\ \hline R_a:\vphantom{\frac1{R_b}}&\hfill -1&\hfill 1&\hfill 0\\ R_b:\vphantom{\frac1{R_b}}&\hfill 1&\hfill 0&\hfill -1\\ R_c:\vphantom{\frac1{R_b}}&\hfill 0&\hfill -1&\hfill 1 \end{smallmatrix} &&&& \begin{smallmatrix} \fbox{A}&N_1&N_2&N_3\\ \hline R_a:\vphantom{\frac1{R_b}}&\hfill 1&\hfill -1&\hfill 0\\ R_b:\vphantom{\frac1{R_b}}&\hfill -1&\hfill 0&\hfill 1\\ R_c:\vphantom{\frac1{R_b}}&\hfill 0&\hfill 1&\hfill -1 \end{smallmatrix} &&&& \begin{smallmatrix} \fbox{A}&N_1&N_2&N_3\\ \hline R_a:\vphantom{\frac1{R_b}}&\hfill -1&\hfill 1&\hfill 0\\ R_b:\vphantom{\frac1{R_b}}&\hfill -1&\hfill 0&\hfill 1\\ R_c:\vphantom{\frac1{R_b}}&\hfill 0&\hfill -1&\hfill 1 \end{smallmatrix} \\\hline \text{Counter-Clockwise}&&&&\text{Clockwise}&&&&\text{Neither} \end{align*}$$
The row vectors are the basis vectors for the input space (node potentials). The \$R_a\$-basis vector is green. The \$R_b\$-basis vector is red. And the \$R_c\$-basis vector is in purple:
In all three cases there's the exact same null space vector, \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}1\\1\\1\end{smallmatrix}\right]\$. (The fact that there's only this exact null space vector tells you that the circuit is "connected".)
It's obvious that these three basis vectors all lie on a plane. They can only span the space of a plane despite the fact that all three nodes can have independently specified potential values. These basis vectors can be scaled any way you like, but only that planar subspace can be spanned by the row vectors.
This follows because you only need any two of the rows, as the remaining third row can be constructed by scaling (with \$\pm 1\$) and summing the other two.
The rest of the 3D input space can be reached. But only by first taking a vector on that plane and then adding some scaled null space vector, \$m\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}1\\1\\1\end{smallmatrix}\right]\$ with \$m\ne 0\$, to it.
There are an infinite number of parallel planes to the subspace plane and these can all be reached by adding some \$m\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}1\\1\\1\end{smallmatrix}\right]\$ to vectors on that planar space that goes through the origin. But none of these additional planes add anything new because all of them will map, via \$A\$, to the same plane in the output space.
Put another way, you can select any set of three node potentials in the input space. And that point may lay on some plane parallel to the one shown above. But there will always be an equivalent point on the subspace plane shown above, which is the only one of the parallel planes going through the origin. (Possessing the origin is a condition to qualify as a space or subspace. Asserting the unique plane that includes the origin is both a useful and necessary condition.)
(Stop for a moment and just think about the fact that despite initially thinking that each node can have any potential assigned to it, it's also true that you really do not have that freedom. You could just as well subtract the same number from all three nodes and still have the same circuit and necessary fluxes into and out of each node. There's really only two degrees of freedom and these two are represented by the subspace plane going through the origin.)
Likewise, the column vectors are the basis vectors for the output space (potential differences across the three edges). Using the same color coding, here's the output spaces for the three cases:
Again, each of these have a null space vector and a planar subspace. This just means that \$A\$ will map from any of the parallel planes in the node potential space to the indicated planar subspace of the potential difference output space. \$A\$ will not be able to map node potentials to anywhere other than the output plane as shown for each of the three directed graph cases. So all of the other parallel planes in the potential difference output space are unreachable by the map, \$A\$.
Note that only in the case of Neither is the null space vector other than \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}1\\1\\1\end{smallmatrix}\right]\$. In that case it is \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}\hfill 1\\\hfill -1\\\hfill 1\end{smallmatrix}\right]\$.
The output space's null space vector tells you that there is a loop!! In the first two cases, since all three values are \$1\$ it says that the loop matches well with the given choice of current directions. In the third case, Neither, it says that the second current direction, that for \$R_b\$, is wrong and must be reversed to make a loop. (Or else the current directions for both \$R_a\$ and \$R_c\$ must be reversed.)
This is the simplest case and I've shown that the null space vector(s) of the output space tell you the loops. And this is always true for any connected circuit.
I'll try to clarify more about why this occurs, shortly.
In this simplest case, \$A\$ just happens to have the same number of edges as nodes. This is rare. Most of the time \$A\$ will not be a nice, square matrix. It will be rectangular with the number of rows almost always exceeding the number of nodes.
As the cases get more complex, it's very hard to show \$m\$-dimensional or \$n\$-dimensional spaces and their subspaces. So I need to shift over to an algebraic way of expressing things. Just keep in mind that the geometry of all this remains the same.
I'm also going to avoid providing a lot of ways of designating current loops. I did it three ways, here. But that point has been made above and I don't need to work it to death.
Just remember that it doesn't matter how one assigns the current directions. The process will always find the independent loops in a connected circuit, regardless of your choices about current directions when making up the incidence matrix, \$A\$.
So let's look at a more complex case.
Before moving on, the weighted Laplacian \$A^T C A\$ (\$C\$ is the conductance matrix which applies Ohm's Law) is positive semidefinite. By asserting the boundary condition that one of the nodes must be zero and thereby reducing \$A^T C A\$ by removing the associated column and row vectors -- keeping in mind the boundary condition as replacing the equation intended by the removed row -- makes the reduced version of \$A^T C A\$ symmetric-positive-definite. This is important as it allows inversion of the reduced Laplacian.
If we assume that there will be a non-zero current into and out of nodes 3 and 4, then here is a Wheatstone bridge circuit:
Here's it's incidence matrix, \$A\$:
$$\begin{align*} \begin{smallmatrix} \fbox{$A$}&N_1&N_2&N_3&N_4\\ \hline R_a:\vphantom{\frac1{R_b}}&\hfill 1&\hfill 0&\hfill -1&\hfill 0\\ R_b:\vphantom{\frac1{R_b}}&\hfill 0&\hfill 1&\hfill -1&\hfill 0\\ R_c:\vphantom{\frac1{R_b}}&\hfill -1&\hfill 1&\hfill 0&\hfill 0\\ R_d:\vphantom{\frac1{R_b}}&\hfill -1&\hfill 0&\hfill 0&\hfill 1\\ R_e:\vphantom{\frac1{R_b}}&\hfill 0&\hfill -1&\hfill 0&\hfill 1 \end{smallmatrix} \end{align*}$$
There's only one way to sum all the columns together to get \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}0\\0\\0\\0\\0\end{smallmatrix}\right]\$. That's by selecting all of them, once. So the input null space vector is \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}1\\1\\1\\1\end{smallmatrix}\right]\$ and this means it is a connected circuit.
But if we transpose the rows into columns:
$$\begin{align*} \begin{smallmatrix} \fbox{$A^T$}&R_a&R_b&R_c&R_d&R_e\\ \hline N_1:\vphantom{\frac1{R_b}}&\hfill 1&\hfill 0&\hfill -1&\hfill -1&\hfill 0\\ N_2:\vphantom{\frac1{R_b}}&\hfill 0&\hfill 1&\hfill 1&\hfill 0&\hfill -1\\ N_3:\vphantom{\frac1{R_b}}&\hfill -1&\hfill -1&\hfill 0&\hfill 0&\hfill 0\\ N_4:\vphantom{\frac1{R_b}}&\hfill 0&\hfill 0&\hfill 0&\hfill 1&\hfill 1 \end{smallmatrix} \end{align*}$$
Then how many ways can the columns be summed together to get \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}0\\0\\0\\0\end{smallmatrix}\right]\$?
Summing \$R_a\$'s column with the negation of \$R_b\$'s column and then adding \$R_c\$'s column would achieve that. That's one null space vector, \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}\hfill 1\\\hfill -1\\\hfill 1\\\hfill 0\\\hfill 0\end{smallmatrix}\right]\$. Another would be by summing \$R_c\$'s column with the negation of \$R_d\$'s column and then adding \$R_e\$'s column. So that's another null space vector, \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}\hfill 0\\\hfill 0\\\hfill 1\\\hfill -1\\\hfill 1\end{smallmatrix}\right]\$. Yet another would be by summing \$R_a\$'s column with \$R_d\$'s column and then adding in the negation of \$R_b\$'s column and \$R_e\$'s column. Yet another null space vector, \$\left[\vphantom{\begin{smallmatrix}\hfill 0\\-1\\+1\end{smallmatrix}}\begin{smallmatrix}\hfill 1\\\hfill -1\\\hfill 0\\\hfill 1\\\hfill -1\end{smallmatrix}\right]\$.
However, only two of these are independent as the remaining one can be generated by appropriate scaling and adding the other two. Any two is sufficient. But perhaps the first two make more sense.
These are called null space vectors for \$A^T\$ because multiplying any of them by \$A^T\$ will result in the null vector. And if we could visualize it, we'd be able to see that they are each perpendicular to the remaining 3D subspace.
In this case, \$A\$ has \$m=5\$ (rows/edges) and \$n=4\$ (nodes/vertices). \$A\$ has one null space vector, meaning the remaining input subspace is \$\mathbb{R}^3\$. (Only three of its four column vectors are independent.) Likewise, \$A^T\$ has two independent null space vectors, again meaning that the remaining output subspace is also \$\mathbb{R}^3\$. (Only three of its five column vectors are independent.) This is to be expected.
The general solution form is \$A^T C A \vec{x}=\vec{f}\$, where \$\vec{x}\$ is the vector of all the node potentials and \$\vec{f}\$ is all of the node fluxes.
Taken in steps:
In general, \$\vec{f}\ne\vec{0}\$ for any interesting circuit. However, the sum of all the elements in \$\vec{f}\$ must be zero (\$0=\sum f_i\$.)
Normally, we'd expect all the node fluxes to be zero, due to KCL. But for any interesting circuit it must be the case that there is a non-zero flux into one node and the negation of that non-zero flux out from another node for the circuit to have non-zero currents in some of the interior edges. There can be more than two-such, of course. But there must be at least two nodes where the node flux is non-zero for the circuit to have non-trivial interior currents.
\$A^TA\$ is called the Laplacian. And \$A^TCA\$ is called the weighted Laplacian. This calculates the divergence of the gradient. The divergence at each node should be zero, except in at least two cases for any interesting circuit where there will be some non-zero divergence at those nodes.
If one of the values in \$\vec{x}\$ is assigned \$0\$ and removed from \$\vec{x}\$ and the Laplacian is modified to remove the associated column and row then the resulting \$A^TCA\$ is called the reduced Laplacian of the graph.
In practical circuits with at least one voltage source so some of the values in \$\vec{x}\$ will be known. Likewise, KCL will set many of the node fluxes in \$\vec{f}\$ to zero. So a technique called the Schur complement is used to block out the Laplacian and generate two different equations to solve for unknown node voltages and unknown node fluxes. (See here for an example.)
(Further, \$\vec{x}^TA^TCA\vec{x}\$ (aka \$\left[A\vec{x}\right]^TCA\vec{x}\$) is the energy!)
The above structure is the exact same form used in so many different fundamental laws! Probability theory requires that the sum of all possibilities is 1. The individual probabilities are like edge fluxes where the sum of each edge must be 1 at a vertex. Same form, entirely different application space.
In any case, the loops (aka meshes if being sloppy) can all be found in the set of independent left null space vectors of \$A\$ (or the independent right null space vectors of \$A^T\$.) This is because all of the node current fluxes must be zero (KCL) in any loop within the circuit.
Mesh current simply means the current of the individual current loops (i.e. meshes) independent of the other branches (i.e. adjacent or coupled meshes). In other words, it's what the current would be around the loop in isolation of the rest of the circuit.
They aren't physically measurable because the idea is a mathematical abstract of what the mesh current would be in isolation. However, you can emulate and simulate the mesh current if you know what the nodal voltages would be without the other meshes in place. The exception to all this would be the trivial case where you only have one mesh i.e. a series circuit. In this case, the mesh current can be physically measured and is equal to the circuit current.
The whole concept of "mesh analysis" is based on the idea of superposition which is what your matrix product in your question is alluding to. The idea of superposition in mesh analysis is basically the sum of the mesh currents is the net circuit current (for a given branch).
The best way to think about this is to think of mesh currents as, maybe, partial currents? Fictitious currents? Useful variables that can be uniquely defined in specific circumstances?
Where this comes from. In a planar network we can draw a specific number of loops N (where loops are defined in the circuits sense, as closed paths that do not cross). We can associate with these loops a current that is imagined to flow around each loop. Then we can write N KVL equations around each of those loops using these current variables. BY THE PROPERTIES OF PLANAR GRAPHS we will have exactly enough INDEPENDENT simultaneous equations to calculate the N imaginary currents. Once we have those we can find everything we want in the circuit (all branch currents, all node voltages, etc.).
(I am trying to recall what was in an old circuits textbook that made this connection between graphs and circuit analysis methods. That book is long discarded, as once you know how to write the minimal number of independent equations, you can go on to DO things with circuits).
(The reason we should be interested in the general methods (node and mesh analysis) is that they provide a systematic way to write exactly enough independent equations to solve the circuit, and not more.)
So the advice to OP and students is to (1) understand how to identify these mesh current things and how to use them and (2) not worry about it.
But this brings up one of the two difficulties with mesh analysis. The first difficulty is that it requires one to use these partial (?) nonphysical currents in the analysis. The second is that the method is limited to planar circuits, and, many interesting circuits are non-planar.
By the way, on the outer edges of a circuit (or graph) only a single mesh current flows through elements. On these outer edges the mesh currents are total physical currents that can be measured.
