In the subsystem formulation of the density-functional theory a large system is decomposed into several constituting fragments that are treated individually. This approach offers the advantage of focusing the attention and computational cost on a limited portion of the whole system while including all the remaining enviromental effects through an effective embedding potential. Here we refer in particular to the (fully-variational) Frozen Density Embedding (FDE) [168] with the Kohn-Sham Constrained Electron Density (KSCED) equations [169,170].

In the FDE/KSCED method the embedding potential required by an embedded subsystem with
density ρ_{A} to account for the presence of another (frozen) subsystem with density ρ_{B}
is:

| (18.1) |

where v_{ext}^{B}(r) and v_{J}[ρ_{B}](r) are the electrostatic potentials generated by the nuclei and electron
density of the subsystem B, respectively, and

| (18.2) |

| (18.3) |

are the non-additive non-interacting kinetic energy and exchange-correlation energy functionals,
respectively. In the expressions above T_{s}[ρ] is the (unknown) non-interacting kinetic energy
density functional and E_{xc}[ρ] is the exchange-correlation energy functional. Note that, while the
first two terms in Eq. (18.1) refer to classical electrostatics (and could be described
by e.g. external point-charges), the last two terms are related to quantum-mechanical
effects.

Using freeze-and-thaw [171] cycles, the role of the frozen and the embedded subsystem is
iteratively exchanged, till convergence. If expressions (18.2) and (18.3) are computed exactly, then
the density ρ_{A} + ρ_{B} will coincide with the exact density of the total system.

Because the FDE/KSCED was originally developed in the Kohn-Sham framework, using
standard GGA approximations for E_{xc}[ρ], the non-additive exchange-correlation potential
(δE_{xc}^{nadd}∕δρ_{A}(r)) can be computed exactly as a functional of the density, leaving the expression
of the non-additive kinetic energy term as the only approximation (with respect the corresponding
GGA calculation of the total system), because the exact explicit density dependence of T_{s} from
the density is not known. Using GGA approximations for the kinetic energy functional
(T_{s} ≈_{s}^{GGA}) we have:

| (18.4) |

and

| (18.5) |

where ṽ_{T}^{GGA}(r) = δ_{s}^{GGA}∕δρ(r).

The FDE total energy of total system is:

Using the Generalized Kohn-Sham (GKS) theory, also hybrid exchange-correlation functionals can be used in embedding calculations. To obtain a practical computational method, the obtained embedding potential must be approximated by a local expression as shown in Ref. [172]. This corresponds to performing for each subsystem hybrid calculations including the interaction with other subsystems through an embedding potential derived at a semilocal level of theory. When orbital dependent exchange-correlation functionals (e.g. hybrid functional and LHF) are considered within the FDE method, the embedding potential includes a non-additive exchange-correlation term of the form

| (18.7) |

where Φ^{A+B}[ρ_{A} + ρ_{B}] denotes the Slater determinant which yields the total density ρ_{A} + ρ_{B}. Since
such a determinant is not easily available, the non-additive exchange-correlation contribution
cannot be determined directly and the non-additive exchange-correlation term can be
approximated as [173]

| (18.8) |