For the first time, a meshless R-functions method (RFM) in combination with the Petrov-Galerkin technique (PGRM) is used for finding a numerical-analytical solution to the system of partial differential equations in terms of the temperature, vorticity and stream functions, describing the natural convection in arbitrarily shaped enclosed cavities. Solutions of model problems, obtained with PGRM, are in good consistence with the results of finite-difference and finite-element modeling. PGRM provides solutions, exactly satisfying arbitrary boundary conditions. An accurate solution is represented in the form of generalized Fourier series with respect to the system of a small number of basic functions (polynomials). Appropriate choice of basic functions for the vorticity function permits to obtain the stream function representation without solving the corresponding partial differential equation.