In recent years, the hybrid quantum mechanics/molecular mechanics (QM/MM) methods have become a rather widely used approach for modeling molecules’ electronic structure in complex environments. Most QM/MM methods partition the molecular system into a small QM region (up to few hundreds of atoms) and a large MM region (thousands of atoms). Full or partial coupling in intermolecular interactions (electrostatics, induction, dispersion, short-range repulsion) is accounted for between the QM and MM parts of the system, where each interaction acts both ways. Non-polarizable or polarizable force fields typically describe the MM part with parameters (point charges, distributed polarizabilities) determined empirically or using quantum chemistry methods. For bonded systems, the MM part of the system’s residues is often subject to bond capping, and the QM and MM parts interface to approximations like atom in-linking. After partitioning the molecular system into QM and MM regions and selecting interaction mode between these systems, the electronic wave function of the QM system can be determined by solving the wave function model specific equation in the presence of external potential generated by the MM region.