One of the critical issues in today's society is improvement of life expectancy. Because of increase in the population of old people and quality of life (welfare), preference of patients to be cured at home and also the constraint of the hospital capacity for patients, requests for Home Health Care (HHC) services are rising. To improve the operation of an HHC provider, with limited resources, we can apply appropriate optimization methods to reduce costs and improve the satisfaction of patients and caregivers. Sometimes an HHC organization with multiple independent depots across a city, might be unable to give on time services. But integrating these depot may decrease the total travel costs and increase patients satisfaction level. In this study, a novel linear mathematical bi-objective model for HHC routing and scheduling problem is presented. The first Objective functio minimizes the total traveling time. The second one minimizes the dissatisfaction level of patients and caregivers. We assume each depot has serval caregivers with various skills. In order to increase caregivers satisfaction level, to balance in total working times of caregivers and their overtimes the second objective function is introduced. Some patient needs more than one service in a day that some of these services are time dependent. Thus each service is considered as a node in this study and each service has a soft time window where the deviation of these time windows lead to increase in patients dissatisfaction level. Since routing and scheduling problems are known NP-Hard problems, a novel multi-objective variable neighborhood search is employed in this study to solve the model. The model results of this algorithm for small-size instances are compared with the epsilon-constraint method solution. For large-scale instances, the results are compared and verified using the solutions of NSGA-II.