The objective of this work is the development of a fully coupled heart model that can serve as a basis for the construction of highly accurate digital twins of the cardiac function. We incorporate in a novel computational framework cardiac electrophysiology, contractile force generation, muscular mechanics, hemodynamics and the circulatory system. The model accounts for electro-mechanical and mechano-electrical feedback, the coupling between muscular deformation and contractile force generation, and fluid-structure interaction (FSI) between the blood and the myocardium. Suitable numerical methods must be employed for the solution of the resulting PDE system. We rely on a geometrically explicit monolithic method for the FSI subproblem, offering a good trade-off among solver robustness, efficiency and accuracy. The coupling of electrophysiology, force generation, and FSI is treated in a segregated-staggered way, leveraging the multiphysics nature of the problem for computational efficiency and flexibility. We implement our solver in a high-performance computing framework. We simulate a realistic human left heart in physiological conditions and compare the numerical results against normal ranges for several biomarkers for ventricular volumes and pressures, flow rates through cardiac valves and the duration of heartbeat phases. The results show that the proposed model reproduces the heart function in healthy conditions. Finally, a whole-heart simulation indicates that previous results can be extended to simulations involving all four cardiac chambers, thus providing an extremely comprehensive representation of the heart. The proposed computational model stands as a milestone towards the development of cardiac digital twins.