Lattice Boltzmann (LB) algorithms are a mesoscopic approach to the solution of nonlinear macroscopic physics, which can be ideally parallelized on supercomputers to the full number of cores available. This is in contrast to the standard direct solution of the nonlinear physics by computational fluid dynamics codes. In its original formulation, LB algorithm consisted of a split operator method: collision relaxation at a spatial node, followed by appropriate streaming of the post-collision distribution to neighboring lattice sites. To reduce memory constraints, a minimal number of lattice streaming vectors are permitted so that one recovers the given nonlinear physics equations in the Chapman-Enskog limit. As in kinetic theory, the macroscopic transport coefficients are determined by the collisional relaxation time in the LB equation. The algorithm is explicit and typically is second-order accurate. However, the simple LB algorithm is prone to numerical instabilities for high Reynolds number flows. In this paper we review some of the recent attempts at stabilizing the LB algorithm for fluid flows – in particular multiple-relaxation operators, introduction of quasi-equiibria states as well as entropic methods. These methods are then discussed in reference to the LB representation of resistive magnetohydrodynamics. The introduction of a vector distribution function for the magnetic field permits the exact satisfaction of the zero divergence of the magnetic field. As a final example, the LB representation of Landau damping is reviewed.