A general method is suggested for the implementation of polarization in molecular dynamics simulations of small molecules. Induced dipole moments are evaluated on selected polarizability centers and represented by separation of charges. The positive polarization charges reside on the selected atoms. The negative polarization charges are treated as additional particles. The positions of these polarization charges are determined from the electrical fields due to the permanent charges of the system. Thus the induction is treated explicitly, while the higher order contributions, the polarization due to induced dipoles, are taken into account in an average way by modification of potential parameters. The forces can be evaluated for the new charge distribution in the conventional way. As an illustration of this approach initial results are reported for the development of a polarizable water model. The higher order polarization is treated in an average way by slight increase of the permanent charges as compared to the values that would give the gas phase dipole moment. The increase in CPU time is comparable to the addition of one atom per polarizable center.