## 1 Introduction
In this model, the rate of change of cell density can be written as a sum of the invasion of cells into nearby tissue (heterogeneous diffusion coefficient times the cell density) and the proliferation of cells (in terms of logistic growth). The resulting equations are similar to the reaction-diffusion equation with radiotherapy and chemotherapy effects [18]. We note that we neglect the effects of radiotherapy and chemotherapy in this paper because, for the time being, we aim to introduce an accurate and effective solution methodology and test its convergence and efficiency on a simplified model. However, we will consider a more realistic model with all additional effects included in a follow-up study [17]. We also note that a _one-dimensional_ simplified model was numerically studied by the co-author of this work in [19]. The major drawback of [19] is that it relies on the construction, storage, and inversion of the Jacobian matrix of the nonlinear system. This approach can be computationally very inefficient and very difficult to implement in a situation where discontinuous or solution-dependent diffusion coefficient exists together with nonlinear source terms, which is the case herein. In this paper, we utilize a Jacobian-free Newton-Krylov (JFNK) method with new effective preconditioning strategies to solve the growth model. The JFNK methods are synergic combinations of Newton type-methods for the super linearly convergent solution of nonlinear equations and Krylov subspace methods for solving the linear equations [20]. The implementation of a JFNK method requires the solution of a linear system that comes from the Newton correction step. Generally, the linear system is solved based on a particular Krylov subspace method and the generalized minimal residual (GMRES) iterative method [20; 21]. The main advantage of a Jacobian-free Newton-Krylov method is that it does not require the formation, storage, and inversion of the Jacobian matrix. Perhaps, this itself is the main reason why the JFNK methods became very popular in the scientific community. They have been successfully used to solve boundary or initial valued problems coming from diverse areas; compressible or incompressible fluid dynamics problems, heat transfer engineering and thermal hydraulic problems, neutronics applications, solid mechanics problems, astrophysical applications, and optimization problems from industry to name a few [20; 22; 23; 24; 25; 26; 27; 28]. As part of the JFNK method, we need to include a preconditioning step to accelerate the convergence of the Krylov block. We note that the preconditioning step became almost standard for any type of iterative scheme [21; 29; 30]. Here, we introduce very effective preconditioning techniques that share similarities to those presented in [31; 32; 33; 34; 35] in such a way that they all rely on the underlying physics. In the literature, these kinds of methodsare referred to as physics-based or physical preconditioning. The main benefit of this preconditioning strategy is that one can target a particular part of the model that creates stiffness or causes slow convergence and only apply the procedure to such a part rather than to the entire system. The numerical results in Section 4 indicate that the new preconditioners perform very well, decreasing the average number of GMRES and the Newton iterations to very low numbers. The organization of this paper is as follows: In Section 2, we define the governing equations. In Section 3, we describe the Jacobian-free Newton-Krylov (JFNK) method together with the new physics-based preconditioning techniques. In Section 4, we present our numerical results to verify the numerical convergence of the proposed schemes and to show the performance and capabilities of the new preconditioning strategies. In Section 5, we summarize our concluding remarks and hint at possible future works.