We have reformulated the reference interaction site model. self-consistent-held (RISM-SCF) method to calculate the molecular properties of solute in solutions. The key feature is to introduce the Helmholtz free energy given as a sum of the solute electronic energy estimated by ab initio molecular orbital (MO) methods and the hypernetted-chain (HNC) excess chemical potential of solvation for the extended reference interaction site model (XRISM) method. The variational conditions for the multiconfigurational self-consistent-held (MCSCF) wave functions were derived. It was also shown that, in the case of the Hartree-Fock (HF) method, the present expression gives the same Fock matrix defined previously. Moreover the expression for first derivatives of the free energy with respect to the solute nuclear coordinates were obtained. Sample calculations for the chemical equilibrium of cis- and trans-isomers of 1,2-difluoroethylene are presented. With the use of analytical energy gradients for the MCSCF and HF methods, the solute molecular geometries were determined in an aqueous solution.