We consider a problem of laser pulse interaction with a nonlinear medium which is accompanied by different nonlinear phenomena. Among them, we highlight the laser pulse self-action, optical bistability realization, formation of laser-induced complicated spatio-temporal structures. For computer modeling of these strongly nonlinear effects, using robust conservative numerical methods is required. Well-known, there are two widely applied approaches for the construction of numerical method: the conservative finite-difference schemes and additive finite-difference schemes (the split-step methods or decomposition methods). The first ones are non-economic, as a rule, while the second type of the methods is economic ones, however they possess well-known disadvantages. In our study, we joint advantages of both approaches by developing an original multi-stage iterative process for the conservative finite-difference scheme realization. Using computer simulation results, we demonstrate the feasibility of the proposed approach for investigating certain nonlinear optical phenomena.