In this paper, we propose a nested simple incomplete LU decomposition (NSILU) method for preconditioning iterative methods for solving largely scale and sparse ill-conditioned hnear systems. NSILU consists of some numerical techniques such as simple modification of Schur complement, compression of ill-condition structure by permutation, nested simple ILU, and inner-outer iteration. We give detailed error analysis of NSILU and estimations of condition number of the preconditioned coefficient matrix, together with numerical comparisons. We also show an analysis of inner accuracy strategies for the inner-outer iteration approach. Our new approach NSILU is very efficient for linear systems from a kind of two-dimensional nonlinear energy equations with three different temperature variables, where most of the calculations centered around solving large number of discretized and illconditioned linear systems in large scale. Many numerical experiments are given and compared in costs of flops, CPU times, and storages to show the efficiency and effectiveness of the NSILU preconditioning method. Numerical examples include middle-scale real matrices of size n = 3180 or n = 6360, a real apphcation of solving about 755418 linear systems of size n = 6360, and a simulation of order n=814080 with structures and properties similar as the real ones.