The following answer presumes some familiarity with the limitations of floating-point arithmetic (rounding, truncation and representation errors), which I will lump together as “rounding error”. It is a trimmed down version of a longer blog post [1] with more detail and some astute comments from subject matter experts.
$M$ in the constraints
The constraint $a'x\le b+My$ illustrates the situation where $M$ shows up as a coefficient in the constraint matrix. This is where serious misadventures can occur. If $M$ appears in the basis matrix, rounding errors can make the basis look singular or make a column that should be ineligible for entry to the basis look eligible. The rounding errors can also cause severe loss of precision in the computation of the basic feasible solution. In technical terms, $M$ (or $1/M$) appearing in the basis matrix can make it ill-conditioned.
Even if ill-conditioned basis matrices do not occur, large values of $M$ can cause branch-and-bound solvers to make slow progress solving the mixed-integer programming model (MIP). Consider the constraint $x_{ij}-My_{j}\le0$, where $x_{ij}$ is the flow of some commodity $i$ across and arc $j$, and $y_{j}$ is a binary variable indicating the presence or absence of that arc. There will be an associated penalty cost (say, $c_{j}y_{j}$) in the objective function for providing arc $j$. Now suppose that, at some point, the solver is considering a node in the search tree where $x_{1j}=2$, $x_{2j}=5$ and $x_{ij}=0$ for other values of $i$. Logically, we know this requires $y_{j}=1$ and incurs cost $c_{j}$; but in the LP-relaxation of the node problem, $y_{j}=5/M$ is sufficient, incurring a cost of just $5c_{j}/M$. For large values of $M$, this substantially underestimates the true cost, leading to loose node bounds. Loose bounds at nodes, in turn, make it hard for the solver to prune nodes based on objective value, and so more nodes need to be examined, slowing the solution process.
The previous example also reveals another problem. Solvers have a tolerance for how close a variable value needs to be in order to consider it an integer. If $5/M$ is small enough, a solver may consider $y_{j}=5/M$ to be an integer (zero) and accept the solution as integer-feasible, when clearly it is not.
Choosing a value for $M$
There are several considerations when choosing values for $M$.
• You want the smallest value of M that is large enough to avoid cutting off a legitimate solution. If $M$ is chosen too small in $a'x\le b+My$, a valid choice of $x$ may violate the constraint even when $y=1$.
• On paper, “big M” models usually use a single symbol $M$ in every constraint. This is usually laziness by the authors. In practice, constraint-specific values of $M$ should be chosen.
• Calibrating $M$ is generally problem-specific. For constraints like $a'x\le b+My$, one possibility is to relax the integrality constraints and maximize $a'x-b$ subject to the other constraints (relaxing others by setting the binary variables to whichever value makes the constraints looser). Done properly, the maximum objective value should provide a safe choice for $M$ (though not necessarily the tightest safe choice).
• There may be context-specific ways to choose a safe value for $M$. An example of one, in the context of discriminant analysis, is in one of my old papers [2].
Alternatives to $M$ models
For situations where constraints need to be turned on or off using binary variables, a possible alternative is “combinatorial Benders decomposition” [3]. This is not automatically superior, but can be effective in some cases.
Some solvers (including but not limited to CPLEX and SCIP) provide “indicator constraints” or “if-then” constraints. The preceding example, posed as an indicator constraint, would essentially be $y_{j}=0\implies x_{ij}=0$. It is up to the solver designers how to process such constraints, but frequently they will simply turn into “big M” constraints with the solver selecting a value for $M$. Even if that does not happen, indicator constraints may weaken LP relaxations.
References
[1] Rubin, Paul. “Perils of 'Big M'.” OR in an OB World, 11 July 2011, https://orinanobworld.blogspot.com/2011/07/perils-of-big-m.html.
[2] Rubin, P. A. (1990). Heuristic Solution Procedures for a Mixed-Integer Programming Discriminant Model. Managerial and Decision Economics 11, 255-266.
[3] Codato, G. & Fischetti, M. (2006). Combinatorial Benders' Cuts for Mixed-Integer Linear Programming. Operations Research 54, 756-766.