Nonlinear modeling of pulsating stars has two features. Firstly, a self-consistent time-dependent model must include not only the optically thin atmospheric region, but also the sub-photospheric envelope, where the instability regions are located. The physical conditions and spatial scales of these regions are very different, while the accuracy of calculations should be as high as possible, especially for the radiative shock waves. As show numerical experiments, the loss of accuracy in the description of the atmospheric dynamics may seriously affect the calculated motions in the inner region, and vice versa, to say nothing of the synthetic line spectrum which is one of the main goals of the atmosphere modelling. Secondly, to obtain the limit cycle we must run the hydromodel for about 103 to 104 pulsational cycles. It is thus clear that we have to simplify the problem of the radiative hydrodynamics in order to obtain reliable result in a reasonable CPU time. In the present talk we describe one of possible approaches which gives good results at least in a semi-quantitative way.