An in silico model for the estimation of volumetric bone mineral density (vBMD) changes at a cortical bone site subjected to mechanobiological bone remodeling is proposed in this manuscript. Mechanisms of cell differentiation, receptor–ligand binding, mechanical signaling, and resorption or deposition of bone matrix were considered, therefore providing a comprehensive description of mechanobiological bone remodeling in the bone microenvironment and enabling the analysis of temporal evolution of disease or therapy scenarios. The proposed model is composed by five modules, namely, bone cells populations, mechanobiology, volume fractions and porosity, mineral density, and structural stiffness. The model is an extension of other models found in the literature because equations for the obtaining of cortical vBMD and the binding of parathyroid hormone (PTH) to parathyroid hormone 1 receptor are included. The proposed model showed a satisfactory agreement with the solutions of other in silico models found in the literature. Simulations of walking and running exercise routines were performed for the evaluation of model capability regarding the control of the numerical error and prediction of vBMD. The computational method used to solve the case study controlled the relative numerical error by less than 1 × 10−7 for approximately 1.7 × 106 time steps. The predicted values correlate with the concept of increasing BMD by vigorous physical activity; however, they contrast with the specific effect of physical activities on cortical vBMD.