## 1 Introduction  
Many phenomena occurring in nature for their investigation can be described via mathematical models based on time-dependent nonlinear diffusion equations [1]. Examples include genetics [2, 3], image processing [4], quantum mechanics [5], and laser-material interactions [6]. Although during the last decades a big effort has been undertaken to find efficient numerical schemes for solution of the corresponding mathematical problem, some of the applications are still a challenging task. Specifically, the efforts in model implementation as well as their demands on the computational power during processing can substantially hinder the theoretical interpretation of the investigated problem. In this work, we consider an application of the nonlinear parabolic diffusion equation to describe the response of solids to an ultrashort laser pulse irradiation. Apart from the insights into the material structure, this topic is important for the description of laser machining [7, 8, 9] and nanostructuring experiments [10, 11, 12] with applications in biotechnology [13] and IT [14, 15]. For metals,the problem may be mathematically formulated in the form of frequently used Two-Temperature Model (TTM) [16], whereas for semiconductors, a similar TTM-like approach has been proposed [17]. The latter is based on the system of partial differential equations, reflecting the conservation laws in the atomic subsystem of a solid and its electronic subsystem. Although it is relatively simple to apply an explicit finite-difference numerical scheme to solve such systems in metals [18] or semiconductors [19], the corresponding stability criteria demand the integration time steps to be small, causing high computational costs as a result. The main restriction on the time step often comes from the nonlinear diffusion equation describing the carrier heat conduction process [20]. One of the possibilities to increase the time step of diffusion equation is to use implicit or semi-implicit integration schemes. For instance, the Crank-Nicolson semi-implicit scheme [21; 22] provides unconditionally stable solution when applied to linear diffusion equations. However, this approach is not directly applicable when nonlinear terms play an important role.  
In this work, we present a semi-implicit finite-difference method for the solution of a system of differential equations--one of which is a diffusion equation with nonlinear terms--and apply it to model short laser pulse interaction with silicon. The presented approach is based on Crank-Nicolson method with predictor-corrector algorithm and provides high stability and precision. It has been already successfully applied for the investigation of ultrashort laser interaction with metals [23] and semiconductors [24; 25; 26]. Section 2 describes the theoretical model and presents the system of equations where a nonlinear diffusion equation results in strong restriction on the time step for the explicit integration algorithm. In Section 3, we give a detailed description of semi-implicit numerical solution scheme, which was modified with predictor-corrector algorithm to account for the nonlinearity in the diffusion equation. Further, in Section 4, the calculation results for a particular set of parameters are presented and the energy conservation versus the applied iteration time step is investigated. Section 5 mentions the existing works, in which this approach has been successfully utilized, and suggests possible improvements for the application of the presented approach in two-dimensional (2D) and three-dimensional (3D) cases. Finally, in Section 6, we give a summary of our results.