Why are quantum mechanical calculations so computationally expensive? What exactly is the process that makes it so time taking?
What I mean is, one can perform classical dynamics simulations on a large system like proteins with no difficulty unlike ab initio calculations. Why? What type of calculations are to be performed in ab initio that makes it expensive? Matrix diagonalization?
It basically comes down to scaling with the size of the system, which will correspond to the number of electrons (that is, occupied valence orbitals) and the size of the basis set (which must have a sufficient quality, that is, size to yield meaningful results).
Edit: I should note that not all quantum chemists consider DFT to be ab-initio. Given the prevalence of such calculations, I have included it in my answer.
Let's consider density-functional theory (DFT) and Hartree-Fock (HF) theory first. These are a ground-state methods that describe a system with a single "wave function" called a Kohn-Sham or Slater determinant, respectively. This is valid for many, but not all systems - especially when there are several d-block metals closely interacting, this can be a bad approximation. For DFT, the cheapest useful specific methods (called functionals) in the most common implementations (i.e. programs) will scale as $\mathcal{O}(N^3)$ (where $N$ is some measure of the system size) due to the calculation of the electron-electron interaction.$^1$ At some point, a matrix diagonalization is performed, which has the same scaling. For HF, the scaling is worse, because of the calculation of a term called exchange (that arises from the specific form of the wave function ansatz), which scales as $\mathcal{O}(N^4)$.
At this point, I should mention that efforts to lower the prefactor of the computational cost function of these workhorses of quantum chemistry are ongoing and some have reached maturity. For "small" systems, the cost is then often dominated by terms other than the highest scaling ones. However, you will eventually hit the scaling wall. There are also efforts to lower the scaling behavior by introducing more approximations, smarter discarding of near-zero terms etc. Again, you may find out the hard way, that in the system you wish to calculate, the approximations do not hold or are not effective.
For excited states, one can assume a $+1$ scaling behavior and a worse prefactor. When considering the ground state with more accuracy, one then turns to post-DFT methods such as double-hybrid density functionals or post-HF methods such as Configuration Interaction, Møller–Plesset perturbation theory, Coupled-Cluster theory (e.g. CISD, MP$n$, CCSD(T)). In the standard implementations, they scale at least as $\mathcal{O}(N^5)$, but can go up to 6 and 7. Again, smart selection of the involved terms can lead to significant speed-ups, leading to an empirical linear or near-linear scaling at the expense of increased prefactors and less certainty.
The comments to the question mention classical dynamics of proteins. The scaling wall even for a single-point energy is hit well before approaching 100 $\ce{CNO}$ atoms (in my personal experience, which is a bit dated. Make it 200, then). The necessary derivatives with respect to the nuclear coordinates also tend to add to the prefactor or the scaling in a significant fashion.
$^1$ This scaling applies if a particular additional approximation called density fitting is made, which is usually a very good one. Note that better DFT methods will have worse scaling.
@TAR86 already described whats makes DFT and HF are computationally expensive. Here I would just like to add, that classical molecular dynamics simulations no explicit treatment of quantum mechanics is considered. This means the fundamental theory (Newton's equations of motion) are much simpler, and therefore easier to calculate. But of course this sacrifices accuracy. As a side note, those two approaches are actually simulating different things: DFT and HF are only about the electrons, while molecular dynamics is about where the atoms are (and move).
Also I want to try to answer the question from a more abstract perspective:
The problem about the equations we need to solve (in both, classical molecular dynamics and quantum mechanics) is that they aredifferential equations. In general, differential equations have no analytic solutions. There are some special cases which have (e.g. particle-in-a-box, harmonic-oscillator, Hydrogen atom), but when it comes to multiple electrons and nuclei, we have a $N$-particle system and their pairwise interactions (be that the Coulomb force in chemistry/physics or gravity in astrophysics), it results in a type of differential equations which have no general analytical solution.
The only remaining option to solve this anyway is using numerical techniques. This essentially means splitting the whole problem in many tiny small steps. The more steps, the higher the accuracy and the higher the computational cost. A simple example would the integrate a 1D function, by splitting its area into many small rectangle or triangle (the above mentioned "steps"), for which we can use an analytic formula, and then add everything up. In HF and DFT the "steps" would be the basis set of orbitals, and in molecular dynamics small steps in time are considered, which yield small steps in the movement of the atoms.
Computers happen to be an excellent tool for numerical approaches, because for each step the same algorithm has to be repeated. So we can easily scale up the number of steps. However if the scaling for the method is higher then linear (which is usually the case), then the we are limited in the size of the problem. Faster or more computers can shift that limit a bit up, but fundamentally the limit remains. This means we cannot just wait for better computers to tackle larger problems, we really need methods with a better scaling behavior.
