In this talk, I advertise for a convex variational functional from statistical mechanics, which is particularly suitable for obtaining the free energy of high-dimensional order parameters from simulational sampling. In the numerical minimization of this variational functional, sampling difficulties related to ergodicity break-down are often alleviated. Two applications will be given. The first one is on Monte Carlo renormalization group simulations, where critical slowing down is eliminated in calculating the renormalized Hamiltonian and static critical exponents. The second one is on the nearest-neighbor 3D spin glass problem, where the Edward-Anderson overlap distribution can be obtained at temperatures much lower than the estimated spin glass transition temperature without the aid of parallel tempering. In these two examples, the optimization is fast and robust. Puzzles generated by this optimization process will also be discussed in the talk.