Conservation planning often involves multiple species occupying large areas including habitat sites with varying characteristics. For a given amount of financial resources, designing a spatially coherent nature reserve system that provides the best possible protection to targeted species is an important ecological and economic problem. In this paper, we address this problem using optimization methods. Incorporating spatial criteria in an optimization framework considering spatial habitat needs of multiple species poses serious challenges because of modeling and computational complexities. We present a novel linear integer programming model to address this issue considering spatial contiguity and compactness of the reserved area. The model uses the concept of path in graph theory to ensure contiguity and minimizes the sum of distances between selected sites and a central site in individual reserves to promote compactness. We test the computational efficiency of the model using randomly generated data sets. The results show that the model can be solved quite efficiently in most cases. We also present an empirical application of the model to simultaneous protection of two cohabiting species, Gopher Tortoise and Gopher Frogs, in a military installation in Georgia, USA.