In this paper we present a three-dimensional numerical model for the radio emission of Magnetic Chemically Peculiar stars, on the hypothesis that energetic electrons emit by the gyrosynchrotron mechanism. For this class of radio stars, characterized by a mainly dipolar magnetic field whose axis is tilted with respect to the rotational axis, the geometry of the magnetosphere and its deformation due to the stellar rotation are determined. The radio emitting region is determined by the physical conditions of the magnetosphere and of the stellar wind. Free-free absorption by the thermal plasma trapped in the inner magnetosphere is also considered. Several free parameters are involved in the model, such as the size of the emitting region, the energy spectrum and the number density of the emitting electrons, and the characteristics of the plasma in the inner magnetosphere. By solving the equation of radiative transfer, along a path parallel to the line of sight, the radio brightness distribution and the total flux density as a function of stellar rotation are computed. As the model is applied to simulate the observed 5 GHz lightcurves of HD37479 and HD37017, several possible magnetosphere configurations are found. After simulations at other frequencies, in spite of the large number of parameters involved in the modeling, two solutions in the case of HD37479 and only one solution in the case of HD37017 match the observed spectral indices. The results of our simulations agree with the magnetically confined wind-shock model in a rotating magnetosphere. The X-ray emission from the inner magnetosphere is also computed, and found to be consistent with the observations.