Parcel delivery is a complex logistic service, as it serves many small or medium-sized customers who may send or receive parcels. Modeling such delivery system needs to integrate two different research areas of hub location and vehicle routing. As it totally depends on the network and the linkage of the nodes, in this paper, some door-to-door service providers are taken into account to provide suitable information for modeling parcel deliveries of sparse and wide countries. Since the proposed mixed-integer programming model is NP-hard, a new multi-steps solution method based on a simulated annealing algorithm and local search is presented. The results of the proposed model and the solution method are evaluated based on some small test problems. The performance of the solution method is illustrated by solving a real case with all capital cities of 31 provinces in Iran.