Since many chemical events occur in solution, it is highly desirable to include solvent effects in the aforementioned accurate quantum chemical calculations. A straightforward approach is to treat solvent molecules explicitly; however, the number of atoms to consider in this manner can be substantial. Even with advancements in the quantum mechanical/molecular mechanical (QM/MM) method, 11 the computational cost of such calculations remains relatively high. One challenge in achieving accurate QM/MM calculations arises from the sampling required to ensure convergence of the free energy, which may be addressed using a reweighting technique. Another class of explicit solvent models includes the effective fragment potential, enabling a more accurate treatment of solvent molecules. On the other hand, solvent effects can also be efficiently considered in an implicit manner. This class includes the generalized Born model, the polarizable continuum model (PCM), the conductor-like screening model (COSMO), and the reference interaction site model (RISM). The first three approaches are known as the dielectric model and have a relatively modest computational demand compared to that of RISM. Among these methods, the PCM approach is widely utilized and has been combined with various computational methodologies. These primarily include single-configuration methods, though some earlier studies have incorporated electron correlation methods. Determining molecular structures using a combination of accurate quantum chemical methods and implicit solvent models, particularly for excited states, remains challenging. One difficulty is developing the analytic derivative with respect to the atomic coordinates. A notable example is the combination of EOM-CCSD or symmetry-adapted cluster-configuration interaction (SAC-CI) with PCM, as demonstrated by Cammi, which has been applied to absorption and emission energies. The integration with multiconfiguration methods is less common, though some earlier studies exist. However, the integration with MRPT is much more limited. Noteworthy developments by Yamazaki (with CASSCF) and Mori employed the RISM method to incorporate non-equilibrium solvent effects and applied it to finding conical intersections (CoIns). Unfortunately, the program code for these methods is not publicly available. Additionally, some studies have successfully incorporated explicit solvents (QM/MM) at the MRPT level, which is useful for including kinetic effects. Technically, the electronic coupling between solute and solvent molecules is typically decoupled in QM/MM calculations, allowing for the use of program code developed for vacuum with minimal modifications. While this approach can be beneficial in many situations, the standard (non-polarizable) QM/MM method employs a fixed charge in the MM region, which is unsuitable when the electronic structure of the solvent varies following an electronic transition. In contrast, implicit solvation models account for a solvent charge dependency on the solute's electronic structure, necessitating additional development for analytic derivatives with PCM and possibly for polarizable QM/MM methods. In this study, we develop analytic first-order derivatives (gradients and non-adiabatic coupling vectors) for the state-averaged (SA) MCSCF and CASPT2 methods in combination with the conductor PCM (C-PCM) within the open-source package OpenMolcas. The CASPT2 energy is obtained using the perturbation to energy (PTE) approach, where the solvent reaction field is generated at the SCF level and remains fixed during the PT2 calculation. While it is theoretically possible to generate the reaction field at the CASPT2 level by iteratively determining the solvent charge and solute density (perturbation to energy and density; PTED), this approach appears impractical for applications and thus is not considered here. This study primarily focuses on the CASSCF-based wavefunction, but the developed method can also be applied to the restricted active space (RAS) SCF and the PT2 extension (RASPT2). The developed method is applied to typical molecules, focusing on solvent effects.