QM/FQ: Quantum Mechanics/Fluctuating Charges (and Fluctuating Dipoles)
**********************************************************************

The Quantum Mechanics/Fluctuating Charges (QM/FQ) and QM/FQ and
Fluctuating Dipoles (QM/FQFμ) methods are multiscale models designed
to describe the properties of a chemical system perturbed by the
presence of its molecular environment, where the latter is described
using a polarizable classical force field; see for example Refs. [1]
[2] [3] [4] [5]. The QM/FQ model has been applied to the calculation
of spectroscopic properties of molecules in aqueous solution and
different solvents (see [6]), whereas QM/FQFμ has only been applied to
aqueous solutions. The method has also been coupled with the GW many-
body-perturbation-theory method [11] . In QM/FQ each atom of the
environment (e.g. a solvation shell) is endowed with a fluctuating
charge (FQ) which can vary in response to the electrostatic potential
generated by the solute, whereas in QM/FQFμ an additional source of
polarization (an induced dipole moment) is placed on each atom. The FQ
and FQFμ charges can be constrained such that each molecule within the
environment remains electrically neutral. The energy functional of the
FQ and FQFμ systems can be written as:

   F(q,\lambda) = \sum_{\alpha, i} q_{\alpha, i} \chi_{\alpha, i} +
   1/2\sum_{i\alpha j\beta} q_{\alpha, i} T^{qq}_{\alpha i, \beta j}
   q_{\beta, j} + \sum_{\alpha}\lambda_\alpha (\sum_i q_{\alpha, i} -
   Q_\alpha)

   F(q,\mu,\lambda) = \frac{1}{2}\sum_{i\alpha j\beta}
   q_{i\alpha}{\text{T}_{i\alpha,j\beta}^{qq}}q_{j\beta} +\frac{1}{2}
   \sum_{ij}\mu_{i\alpha}^{\dagger}\mathbf{T}_{i\alpha,j\beta}^{\mu\m
   u}{\mu}_{j\beta} + \sum_{ij}q_{i\alpha}\textbf{T}_{i\alpha,j\beta}
   ^{q\mu}{\mu}_{j\beta}^{\dagger} +
   \sum_{i\alpha}q_{i\alpha}\chi_{i\alpha} +
   \sum_{\alpha}\lambda_{\alpha}\left[\sum_{i} q_{\alpha
   i}-Q_{\alpha}\right]

where q_{\alpha, i} and \mu_{\alpha, i} are the i-th FQ charge and the
i-th dipole belonging to the \alpha-th molecule, \chi_{\alpha,i} is
the electronegativity of each atom, \lambda_\alpha is the Lagrange
multiplier associated with the \alpha-th molecule, whose purpose is to
fix the total charge to be Q_\alpha. The matrix T^{qq}_{\alpha i,
\beta j} describes the interaction kernel between the charges: while
the diagonal elements are related to the self-interaction through the
chemical hardnesses eta, the off-diagonal elements may be specified
through different formulations. The T_{i\alpha,j\beta}^{q\mu} and
T_{i\alpha,j\beta}^{\mu\mu} matrices describe the interaction kernel
between charges and dipoles and between the dipoles, respectively. The
electric variables (the FQ charges and the dipoles) can be obtained by
minimizing the energy functional with respect to the charges, the
dipoles and the Lagrange multipliers, which leads to the following set
of linear equations:

   \left(\begin{array}{cc|c} \mathbf{T}^{qq} & \mathbf{1}_{{\lambda}}
   & \mathbf{T}^{q\mu} \\ \mathbf{1}^{\dagger}_{{\lambda}} &
   \mathbf{0} & \mathbf{0}   \\ \hline -\mathbf{T}^{q\mu^{\dagger}} &
   \mathbf{0} & \mathbf{T}^{\mu\mu} \end{array} \right)
   \left({\begin{array}{c} \mathbf{q}\\ {\lambda}\\ \hline {\mu}
   \end{array}}\right) = \left(\begin{array}{c} -{\chi} \\ \mathbf{Q}
   \\ \hline \mathbf{0} \end{array}\right) \qquad \Rightarrow \qquad
   \mathbf{D}\mathbf{Q}_{\lambda} = -\mathbf{C}_Q

where \mathbf{1}_{{\lambda}} is a rectangular matrix which accounts
for the Lagrangians, \mathbf{C}_Q is a vector containing atomic
electronegativities and total charge constraints. \mathbf{Q}_\lambda
is a vector containing charges, dipoles and Lagrange multipliers.

The coupling of the FQ/FQFμ models with a QM Hamiltonian can be done
by introducing the QM/FQ(Fμ) interaction operator as

   H_\text{QM/FQ} = \sum_i^{N} q_i V[\rho](r_i), H_{\text{QM/FQF}\mu}
   = \sum_i^{N} q_i V[\rho](r_i) - \sum_i^{N} \mu_i E[\rho](r_i),

where N is the number of atoms, q_i and \mu_i are the i-th FQ charge
and dipole located at position r_i, respectively. V[\rho](r_i) and
E[\rho](r_i) are the electric potential and field generated by the QM
system at the same point. The introduction of the QM/FQ interaction
leads to a modified set of linear equations for the FQ charges, i.e.

   DQ_\lambda = -C_Q-R[\rho]

where R is an array collecting the electric potential and field
generated by the QM electrons and nuclei at the position of each atom.
The FQ charges (and dipoles) are thus determined self-consistently
along with the ground-state density. Since the charges (and dipoles)
depend on the QM density, explicit terms also appear within response
equations that are solved to simulate spectroscopic and excited-state
properties of the QM system.

Starting from AMS2021.103 a screening is included for the interaction
between MM atoms and the QM density, to avoid unstable results in case
numerical integration points are accidentally close to MM atoms. The
screened 1/r_{ij} has the form:

   \text{ERF: =}  & \frac{erf(ar_{ij})}{r_{ij}} \\ \text{EXP: =}  &
   \frac{(1-exp(-ar_{ij}))^2}{r_{ij}} \\ \text{GAUS: =} &
   \frac{erf(r_{ij}/a)}{r_{ij}}

where a can be changed with the QMSCREENFACTOR keyword.


Gradients and geometry optimization
===================================

Starting from AMS2025 the QM/FQ(Fμ) model can be used when calculating
gradients and geometry optimizations. QM/FQ(Fμ) is a "focused model"
developed to reproduce the properties of the central QM system under
the perturbation of the environment, and not the property of the
entire system (QM + environment). In keeping with this philosophy,
whenever geometries or gradients are calculated, the atoms in the
environment are kept frozen in place and the gradients are only
evaluated for the QM atoms. While in their simplest form these are
purely electrostatic models, when evaluating gradients it is crucial
to include non-electrostatic contributions, particularly repulsion
forces that keep the QM atoms from approaching the MM atoms. This can
be done through the "NonEle" keyword (see below). At present, only
ground-state gradients are implemented.


Coupling with FDE
=================

Starting from AMS2022 the coupling between Frozen Density Embedding
(FDE) and FQ is implemented (see Ref. [7]). When a FDE calculation is
set up, the FQ(Fμ) environment atoms respond to both layers of the FDE
embedding.

In the case of an excited-state calculation performed with the Time-
Dependent DFT method (TDDFT), the atoms belonging to the frozen
density layer do not normally respond dynamically to the external
field. This can create an imbalance because both the non-frozen QM
layer and the fluctuating charge (and dipole) layer both do. A very
cost-effective way to include the response of the frozen-density layer
into the calculation is to use the FQ(Fμ) model to estimate it (see
Ref. [7] for more details) provided the parameters defining the
molecular entities in the frozen layer are available.

This can be accomplished by simply including the FDERESP keyword in
the input (see below).


Input options
=============

The scheme below presents all keywords that can be used in a QM/FQ(Fμ)
calculation. Further down below you can find a detailed explanation of
how to specify the coordinates of the FQ(Fμ) atoms.

   QMFQ
      AtomType
         Alpha float
         Area float
         Charge float
         Chi float
         Dfermi float
         Efermi float
         Eta float
         R0 float
         Sfermi float
         Sigma0 float
         Symbol string
         Tau float
         wFQFMuFile string
      End
      Coords # Non-standard block. See details.
         ...
      End
      FDERESP Yes/No
      Forcefield [FQ | FQFMU | wFQ | wFQ_RQ | wFQFmu | wFQFmu_RQRMU]
      Frozen Yes/No
      Kernel [OHNO | COUL | GAUS]
      LifeFrequency float
      LifeState integer
      MolCharge float
      NonEle [LJ | None]
      QMSCREEN [ERF | EXP | GAUS | NONE]
      QMSCREENFACTOR float
      ReadDynamicMatrix string
   End

"QMFQ"
   Type:
      Block

   Description:
      Block input key for QM/(w)FQ(FMu).

   "AtomType"
      Type:
         Block

      Recurring:
         True

      Description:
         Definition of atomic types in MM environment

      "Alpha"
         Type:
            Float

         Description:
            Polarizability of FQFMU atom

      "Area"
         Type:
            Float

         Description:
            Effective charge-transfer area time between same-type wFQ
            atoms

      "Charge"
         Type:
            Float

         Description:
            MM fixed charge (non-polarizable only)

      "Chi"
         Type:
            Float

         Description:
            Electronegativity of FQ atom

      "Dfermi"
         Type:
            Float

         Description:
            wFQ d parameter for Fermi damping

      "Efermi"
         Type:
            Float

         Description:
            Fermi Energy of wFQ carbon-based structure

      "Eta"
         Type:
            Float

         Description:
            Chemical Hardness of FQ atom

      "R0"
         Type:
            Float

         Description:
            wFQ R0 parameter for Fermi damping

      "Sfermi"
         Type:
            Float

         Description:
            wFQ s parameter for Fermi damping

      "Sigma0"
         Type:
            Float

         Description:
            Static conductivity time of wFQ atom

      "Symbol"
         Type:
            String

         Description:
            Symbol associated with atom type

      "Tau"
         Type:
            Float

         Description:
            Scattering time of wFQ atom

      "wFQFMuFile"
         Type:
            String

         Description:
            The absolute path to the wFQFMu interband polarizability
            file

   "Coords"
      Type:
         Non-standard block

      Description:
         Coordinates and fragment information (FQ only)

   "FDERESP"
      Type:
         Bool

      Default value:
         No

      Description:
         In response calculations (TD), the polarization contribution
         of the FDE part is introduced at the FQ level [See F. Egidi
         et al. J. Chem. Phys. 2021, 154, 164107].

   "Forcefield"
      Type:
         Multiple Choice

      Default value:
         FQ

      Options:
         [FQ, FQFMU, wFQ, wFQ_RQ, wFQFmu, wFQFmu_RQRMU]

      Description:
         Version of the FQ family of polarizable forcefields

   "Frozen"
      Type:
         Bool

      Default value:
         No

      Description:
         Expert option. Do not introduce polarization effect in
         response calculations.

   "Kernel"
      Type:
         Multiple Choice

      Default value:
         OHNO

      Options:
         [OHNO, COUL, GAUS]

      Description:
         Expert option. KERNEL can be used to choose the functional
         form of the charge-charge interaction kernel between MM
         atoms. Recommended is to use the default OHNO. The COUL
         screening is the standard Coulomb interaction 1/r. The OHNO
         choice introduce the Ohno functional (see [K. Ohno, Theoret.
         Chim. Acta 2, 219 (1964)]), which depends on a parameter n
         that is set equal to 2. Finally, the GAUS screening models
         each FQ charge by means of a spherical Gaussian-type
         distribution, and the interaction kernel is obtained
         accordingly. For QM/FQFMU only GAUS SCREEN is implemented.

   "LifeFrequency"
      Type:
         Float

      Default value:
         0.0

      Description:
         Input frequency for calculating TDDFT/wFQ(Fmu) lifetimes

   "LifeState"
      Type:
         Integer

      Default value:
         0

      Description:
         Target excited state for state-specific TDDFT/wFQ(Fmu)
         lifetimes

   "MolCharge"
      Type:
         Float

      Default value:
         0.0

      Description:
         Total charge of each fragment (FQ only)

   "NonEle"
      Type:
         Multiple Choice

      Default value:
         LJ

      Options:
         [LJ, None]

      Description:
         Whether to include non-electrostatic contributions to the
         energy. Default is the Lennard-Jones (LJ) model.

   "QMSCREEN"
      Type:
         Multiple Choice

      Default value:
         GAUS

      Options:
         [ERF, EXP, GAUS, NONE]

      Description:
         Expert option. QMSCREEN can be used to choose the functional
         form of the charge-charge interaction kernel between MM atoms
         and the QM density. The screening types available are ERF
         (error function), EXP (exponential), GAUS (Gaussian), or
         NONE. The default is GAUS.

   "QMSCREENFACTOR"
      Type:
         Float

      Default value:
         0.2

      Description:
         Expert option. Sets the QM/MM interaction kernel screening
         length. Recommended is to use the default value 0.2 with the
         GAUS QM/MM screening function.

   "ReadDynamicMatrix"
      Type:
         String

      Description:
         The absolute path to the decomposed wFQ(FMu) matrix file

The FQ(Fμ) region is specified within the Coords block, using the
following scheme.

   Coords
      elem1 x.xxx y.yyy z.zzz N
      elem2 x.xxx y.yyy z.zzz N
      elem3 x.xxx y.yyy z.zzz N
      ...
   End

You will need to replace 'elem#' with the element you are assigning
parameters to, and after the coordinates you need to include a number
N assigning the element to a specific molecule. The latter number is
important because the total charge of each environment molecule is
constrained to be a specific value (0 by default), which can be
specified with the MolCharge keyword if nonzero (see above scheme).
For instance, the following block specifies two water molecules, where
the oxygen and hydrogen atoms have been given the names OW and HW
respectively:

   Coords
       OW   -0.9957202   0.0160415   1.2422556   1
       HW   -1.4542703  -0.5669741   1.8472817   1
       HW   -0.9377950  -0.4817912   0.4267562   1
       OW   -0.2432343  -1.0198566  -1.1953808   2
       HW    0.4367536  -0.3759433  -0.9973297   2
       HW   -0.5031835  -0.8251492  -2.0957959   2
   End

Parameters must have been assigned for the atom types OW and HW in the
respective AtomType blocks.

In case of water, it is recommended to use the optimized FQ parameters
Chi and Eta for O and H in water proposed in Ref. [9]

   AtomType
     Symbol O
     Chi 0.189194
     Eta 0.523700
   End
   AtomType
     Symbol H
     Chi 0.012767
     Eta 0.537512
   End

Alternatively, in case of water as solvent, one could use the
parameters proposed in Ref. [10]

   AtomType
     Symbol O
     Chi 0.116859
     Eta 0.584852
   End
   AtomType
     Symbol H
     Chi 0.000001
     Eta 0.625010
   End

For other solvents, it is recommended to use the parameters reported
in Ref. [6].

For QM/FQFμ calculations for water as solvent, the parameters reported
in Ref. [4] can be used

   AtomType
     Symbol O
     Chi   0.290840
     Eta   0.562510
     Alpha 2.218790
   End
   AtomType
     Symbol H
     Chi   0.167570
     Eta   0.609320
     Alpha 1.190640
   End

[1] T. Giovannini, F. Egidi, C. Cappelli, *Molecular spectroscopy of
    aqueous solutions: a theoretical perspective*, Chemical Society
    Reviews, 49, 5664 (2020)

[2] C. Cappelli, *Integrated QM/polarizable MM/continuum approaches to
    model chiroptical properties of strongly interacting
    solute–solvent systems*, International Journal of Quantum
    Chemistry, 116, 1532 (2016)

[3] T. Giovannini, F. Egidi, C. Cappelli, *Theory and algorithms for
    chiroptical properties and spectroscopies of aqueous systems*,
    Physical Chemical Chemical Physics, 22, 22864 (2020)

[4] T. Giovannini, A. Puglisi, M. Ambrosetti, C. Cappelli,
    *Polarizable QM/MM approach with fluctuating charges and
    fluctuating dipoles: the QM/FQFμ model*, Journal of Chemical
    Theory and Computation, 15, 2233-2245 (2019)

[5] T.Giovannini, R. R. Riso, M. Ambrosetti, A. Puglisi, C. Cappelli,
    *Electronic transitions for a fully polarizable QM/MM approach
    based on fluctuating charges and fluctuating dipoles: linear and
    corrected linear response regimes*, The Journal of Chemical
    Physics, 151, 174104 (2019)

[6] M.Ambrosetti, S. Skoko, T. Giovannini, C. Cappelli, *Quantum
    Mechanics/Fluctuating Charge Protocol to Compute Solvatochromic
    Shifts*, Journal of Chemical Theory and Computation (2021)

[7] F.Egidi, S. Angelico, P. Lafiosca, T. Giovannini, C. Cappelli, *A
    polarizable three-layer frozen density embedding/molecular
    mechanics approach*, The Journal of Chemical Physics, 154, 164107
    (2021)

[8] K. Ohno, *Some remarks on the Pariser-Parr-Pople method*,
    Theoretica Chimica Acta, 2, 219 (1964)

[9] T. Giovannini, P. Lafiosca, B. Chandramouli, V. Barone, C.
    Cappelli, *Effective yet reliable computation of hyperfine
    coupling constants in solution by a QM/MM approach: Interplay
    between electrostatics and non-electrostatic effects*, Journal of
    Chemical Physics 150 (2019) 124102

[10] S.W. Rick, S.J. Stuart, B.J. Berne, *Dynamical fluctuating charge
     force fields: Application to liquid water*, Journal of Chemical
     Physics 101 (1994) 6141

[11] 7. Nottoli, P. Lafiosca, F.E.Q. Rodríguez, F. Egidi, A. Förster,
        C. Cappelli, *GW Approximation Coupled with Classical
        Fluctuating Charges and Dipoles*, Journal of Chemical Theory
        and Computation 22 (2026) 371–384
