Geodynamic modeling is often related with challenging computations involving solution of the Stokes and continuity equations under the condition of highly variable viscosity. Based on a new analytical approach we have developed particular analytical solutions for 2-D and 3-D incompressible Stokes flows with both linearly and exponentially variable viscosity. We demonstrate how these particular solutions can be converted into 2-D and 3-D test problems suitable for benchmarking numerical codes aimed at modeling various mantle convection and lithospheric dynamics problems. The Main advantage of this new generalized approach is that a large variety of benchmark solutions can be generated, including relatively complex cases with open model boundaries, non-vertical gravity and variable gradients of the viscosity and density fields, which are not parallel to the Cartesian axes. Examples of respective 2-D and 3-D MatLab codes are provided with this paper.