QM program calculating forces exerted on point charges?



Hello,
 many QM programs allow to calculate energy and gradients of a molecule under
 the influence of point charges.  For typical QM/MM calculations (coupling term)
 it would be necessary to also calculate the gradients/forces exerted _on_ these
 point charges.  Currently I'm not aware of many programs (maybe I'm missing
 something?) doing that: GAMESS-US for example has been interfaced to CHARMM,
 but obviously the calculation of point charge gradients is a special extension
 for CHARMM (the GAMESS source code has to be prepared for that and requires
 CHARMM to do proper setup).  This makes it somewhat difficult to adapt the
 program for our own purposes.
 Now I wonder if there are any QM programs which output these additional
 gradients by default (or probably by user request).  Since many QM programs
 already allow the calculation of energy/gradients with point charges, would it
 be difficult to enhance the program to compute these additional gradient?
 Thanks for any help,
 Hannes.