Flexible intelligent metasurfaces (FIMs) constitute a promising technology that could significantly boost the wireless network capacity. An FIM is essentially a soft array made up of many low-cost radiating elements that can independently emit electromagnetic signals. What's more, each element can flexibly adjust its position, even perpendicularly to the surface, to morph the overall 3D shape. In this paper, we study the potential of FIMs in point-to-point multiple-input multiple-output (MIMO) communications, where two FIMs are used as transceivers. In order to characterize the capacity limits of FIM-aided narrowband MIMO transmissions, we formulate an optimization problem for maximizing the MIMO channel capacity by jointly optimizing the 3D surface shapes of the transmitting and receiving FIMs, as well as the transmit covariance matrix, subject to a specific total transmit power constraint and to the maximum morphing range of the FIM. To solve this problem, we develop an efficient block coordinate descent (BCD) algorithm. The BCD algorithm iteratively updates the 3D surface shapes of the FIMs and the transmit covariance matrix, while keeping the other fixed. Numerical results verify that FIMs can achieve higher MIMO capacity than traditional rigid arrays. In some cases, the MIMO channel capacity can be doubled by employing FIMs.