  Cooling of neutron stars with dipole magnetic fields is simulated using a
realistic model of the anisotropic surface temperature distribution produced by
magnetic fields. Suppression of the electron thermal conductivity of outer
stellar layers across the field increases thermal isolation of these layers
near the magnetic equator. Enhancement of the radiative and longitudinal
electron thermal conductivities in quantizing magnetic fields reduces thermal
isolation near the magnetic poles. The equatorial increase of the isolation is
pronounced for $B \ga 10^{10}$ G, while the polar decrease -- for $B \ga
10^{12}$ G. The effects compensate partly each other, and the actual influence
of the magnetic fields on the cooling is weaker than predicted by the
traditional theories where the equatorial effects have been neglected.
