The modeling of biochemical processes in biofilms is more complex compared to those in suspended biomass due to the existence of substrate gradients. The diffusion and reaction of substrates within the biofilms were simulated in 1D models in the 1970s. The quality of these simulation results was later improved by consideration of mass transfer at the bulk/biofilm interface and detachment of biomass from the surface. Furthermore, modeling of species distribution along the axis perpendicular to the substratum helped to simulate productivity and long-term behavior in multispecies biofilms. Multidimensional models that were able to give a realistic prediction of the surface structure of biofilms were published in the 1990s. The 2D or 3D representation of the distribution of the species in a matrix of extracellular polymeric substances (EPS) helped predict the behavior of multispecies biofilm systems. The influence of shear forces on such 2D or 3D biofilm structures was included by solving the Navier-Stokes equation for the liquid phase above the biofilm. More recently, the interaction between the fluid and biofilm structures was addressed more seriously by no longer considering the biofilm structures as being rigid. The latter approach opened a new door, enabling one to describe biofilms as viscoelastic systems that are not only able to grow and produce but also be deformed or even dislodged if external forces are applied.