In various fields of natural science, the chaotic systems of differential equations are considered more than 50 years. The correct prediction of the behaviour of solutions of dynamical model equations is important in understanding of evolution process and reduce uncertainty. However, often used numerical methods are unable to do it on large time segments. In this article, the author considers the modern numerical method and algorithm for constructing solutions of chaotic systems on the example of tumor growth model. Also a modification of Benettin's algorithm presents for calculation of Lyapunov exponents.