# A Theoretical Perspective on the Thermodynamic Stability of Polymer Blends for Solar Cells: From Experiments to Predictive Modeling

Claudia Caddeo, Jörg Ackermann, Alessandro Mattoni

► To cite this version:

Claudia Caddeo, Jörg Ackermann, Alessandro Mattoni. A Theoretical Perspective on the Thermodynamic Stability of Polymer Blends for Solar Cells: From Experiments to Predictive Modeling. Solar RRL, 2022, pp.2200172. 10.1002/solr.202200172. hal-03763289

HAL Id: hal-03763289

https://hal.science/hal-03763289v1

Submitted on 29 Aug 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# A Theoretical Perspective on the Thermodynamic Stability of Polymer Blends for Solar Cells: From Experiments to Predictive Modeling

Claudia Caddeo,* Jörg Ackermann, and Alessandro Mattoni*

An overview of the theoretical/computational methods that allow the study of the thermodynamic stability of the polymer blends for photovoltaics is provided. After discussing the fundamental concepts of thermodynamic blend stability and solubility, including mixing enthalpy and the Flory-Huggins theory, some experimental approaches to determine the interaction parameter and the stability in organic photovoltaic (OPV) are briefly discussed, and the advances in the modeling of polymer blends based on the use of atomistic simulations and multiscale modeling are reviewed. An outlook on the modeling strategies that can have a strong impact on the design of stable blends and to the commercial OPV technologies is given. In particular, the main directions along which major developments are expected are envisaged: multiscale models with improved accuracy or machine learning methods applied to large ab initio datasets, as well as a judicious combination of the two strategies.

## 1. Introduction

Organic photovoltaic (OPVs) solar cells have attracted great interest in the search for efficient solar cell technologies due to their promising performance and potential for low-cost manufacture.[1,2] Among their characteristics of special relevance are solution processability, flexibility, large scale, and low-cost production. In particular, organic cells can be manufactured over large areas, on lightweight plastic substrates, using high-throughput printing fabrication, potentially resulting in large reductions in production costs and energy payback time.[3]

C. Caddeo, A. Mattoni
Istituto Officina dei Materiali (CNR - IOM)
Cagliari Cittadella Universitaria
I-09042 Monserrato, CA, Italy
E-mail: caddeo@iom.cnr.it; mattoni@iom.cnr.it

J. Ackermann
Aix Marseille University
CINaM UMR CNRS 7325
13288 Marseille, France

<figure><img src="image_2.png" /></figure>

The ORCID identification number(s) for the author(s) of this article
can be found under https://doi.org/10.1002/solr.202200172.

© 2022 The Authors. Solar RRL published by Wiley-VCH GmbH. This is an
open access article under the terms of the Creative Commons Attribution
License, which permits use, distribution and reproduction in any medium,
provided the original work is properly cited.

DOI: 10.1002/solr.202200172

The first examples of organic solar cells (OSCs) featured fullerene derivatives (such as PC60BM and PC70BM) as electron acceptors; polymer-fullerene blends can be considered the prototypical material realizing the bulk heterojunction architecture[4-7] and have been investigated for more than two decades with a rapid increase in power conversion efficiencies (PCE) that surpassed 10% in 2015,[8] since then followed by slower improvements, with a record value of 11.5% for single-junction cells and 13.2% for tandem cells.[9,10] Recently nonfullerene organic semiconductors outperformed the fullerene-based cells in terms of PCEs. In 2019, the discovery of the nonfullerene acceptor (NFA) Y6 led to the realization of a single-junction organic solar cell with PCE > 15%[11] in combination with the large-bandgap donor polymer PM6.[12] Since then, PM6:Y-series blends achieved over 19% PCE.[13] Very efficient charge separation[14] and long carrier diffusion[15] are key aspects of these blends. Several groups have carried out investigations to unravel the reasons behind their exceptional performance. Surface and bulk characterization techniques, as well as theoretical simulations, have been employed to understand how phase separation and other morphological features affect the PCE in these blends.[16-23]

Despite the great progress of organic solar cells, the photovoltaic market is still dominated by crystalline silicon (c-Si) technologies due to the material availability, its excellent electronic charge transport properties, and PCEs >24% for large module areas >1 m². Moreover, the lifetime of c-Si solar cells typically exceeds 25 years[24,25]: thermal and photochemical stability of inorganic semiconductors is in fact generally higher than that of organic counterparts. The operational lifetimes of large-area OPV devices are instead still significantly lower than the market requirements of ≥10 years.[26] When considering the so-called “golden triangle,” used to assess the feasibility of photovoltaic technologies (efficiency, cost, and lifetime),[27,28] it is evident that OPVs fail in meeting the lifetime criteria to be considered competitive with c-Si.

Considerable effort is dedicated to improving OPV efficiency in
academic research, with efficiencies approaching rapidly 20% for
small-area solar blends with the advent of NFAs, as discussed earlier.

On the contrary, the fundamental study of stability has
received much less attention.[29]

Fullerene-based blends have been studied for a long
time[4-7,30,31] and are the systems for which more experimental

<figure><img src="image_3.png" /></figure>
information is available. Remarkable results have been obtained for the PBDTTT-OFT:PC70BM blend, which retained 80% of its initial efficiency under accelerated testing conditions for over 500 h in air at temperature of 85 °C. The storage time at room temperature of the encapsulated devices based on this blend was 2000 h.[32]

Concerning NFA solar cells, several studies have established that they exhibit superior long-term stability than fullerene-based ones[33-36] and further improvement requires a deep understanding of the morphology, thermodynamics, and degradation processes. The UV-induced degradation at the interface between the active layer and the charge transporting layers is one of the main sources of instability. Important results have been obtained for PTB7-th:BT-CIC blends by introducing a buffer between the active and the charge transporting layers, in addition to a UV filter. By performing accelerated aging, it was possible to extrapolate an intrinsic lifetime equivalent to 30 years of outdoor exposure.[37] Photoinduced and thermal degradation modes at the active layer contacts were individuated as responsible for the PCE loss in the absence of the buffer layers.[38]

In general, for both fullerene and nonfullerene solar cells, the lifetime ultimately depends on the morphological, chemical, and photochemical stability of the active layer.[39] At the nanometric scale, a bulk heterojunction (BHJ) is constituted by a good intermix of donor/acceptor domains with similar charge carrier mobilities within the exciton diffusion lengths (nm scale). The BHJ architecture increases the interfacial area between donor and acceptor, overcoming the limits of bilayer junctions in terms of short exciton diffusion length, limited exciton lifetime, and charge separation.[40] Nevertheless, the BHJ structure also presents challenges. The multiple phases and complex interfaces bring along unexpected morphologies and complicated charge dynamics, which are difficult to be observed and controlled. Mixing two or more components and achieving an optimum morphology is challenging[41] and usually involves a trial-and-error approach. Furthermore, bulk heterojunctions are metastable microstructures and are prone to degradation. The factors that contribute to degradation of the active layer include photochemical and photophysical degradation of the active layer materials and interfaces,[42-44] morphological degradation due to high thermal stress under operation, and[45,46] degradation events initiated by H₂O and O₂ in poorly sealed devices.[47-50]

Furthermore, also solvents and additives have important effects on the final morphology and stability. Solvent additives, in fact, are important to manipulate the kinetics of phase separation of the active layer during deposition. A comprehensive review on this issue can be found in the study by McDowell et al.[51] and studies therein.

Roughly speaking, it is possible to distinguish degradation mediated by external factors (chemical contamination, exposure to external excitations such as ultraviolet radiation, mechanical actions) from degradation related to the intrinsic thermodynamic stability of the material at given thermodynamic parameters (e.g., temperature, pressure). Extrinsic degradation mechanisms are discussed in many existing reviews[39,50,52-60] because of their relevance under working conditions. From a theoretical point of view, they involve photochemical reactions that require a specific description of the excitations, interfaces, or chemical impurities.

However, even just considering the intrinsic stability of the material with temperature, we are faced with a topic of great

complexity that represents an important theoretical and applica-
tive scientific challenge and that is the focus of the present article.
While extrinsic factors can be, to some extent, mitigated by
production steps without modifications of the active material
(e.g., encapsulation, contacts), improving the intrinsic stability
requires optimizing molecular interactions and the processed
microstructure of the donor/acceptor blend. The control of
thermodynamic stability of polymer blends is important not only
for the active layer but also for the hole and electron transporting
layers that are often realized with polymeric systems.

Among the factors affecting intrinsic instability, the nanomor-
phological degradation and phase separation of the photoactive
layer are certainly relevant, which give rise to the formation of
nanometer-scaled interpenetrating domains of pure donor, pure
acceptor, and mixed phases. The mixed phases play a fundamen-
tal role in the overall performance of the photovoltaic cell.
A well-mixed phase is in fact needed for efficient phase separa-
tion (and consequently high short-circuit current) Jsc but too
small domains can increase the charge recombination and thus
lower the fill factor (FF).[61-63]

It is necessary to make progress toward the full characteriza-
tion of the morphology and the crystalline order at the nanoscale
for improving charge transport and for blend stability.[64-68]
The determination of the BHJ structure typically requires a com-
bination of complementary analytical techniques, such as X-Ray
diffraction and solid-state nuclear magnetic resonance
(ss-NMR).[69] In particular, ss-NMR is sensitive to short-range
interactions, allowing it to probe diverse material structures,
such as crystallites, lamellar mesophases, and amorphous
regions, and making it particularly suitable for studying
organic materials. Integrated approaches combining density
functional theory (DFT) calculations or molecular dynamics
(MD) simulations to experimental techniques have been also
successfully used.[70,71] In the case of blends based on NFAs
(e.g., PBDB-T:ITIC), however, the direct characterization of the
molecular arrangements at the interface and within the acceptor
phase is challenging (e.g., due to the similar electronic structures
giving overlapping signals). Promising results come from recent
studies based on analytical electron microscopy that enables
material phase identification at the nanometer scale.[64]

The possibility to apply in silico methods for screening
donor/acceptor pairs to identify those that might stabilize the
BHJ blends, analogous with virtual screens to identify hit com-
pounds that bind with a target protein,[72-75] is very appealing.
Molecular dynamics simulations, either with interatomic forces
calculated from first principles, (i.e., ab initio molecular dynam-
ics, AIMD) or derived from classical force fields (i.e., model
potential molecular dynamics, MPMD), can predict, at least in
principle, the full thermodynamics of the blends and their
miscibility.[76] Recently, high-throughput methods based on
machine learning (ML) algorithms have emerged as methods
to screen donor/acceptor molecules which should in principle
provide the best performances in OPV.[77,78] Nevertheless, the
stability is rarely taken into account in these kinds of studies,
which mainly focus on the prediction of electronic properties
(such as, e.g., energy-level alignment) and overall PCE.

An ideal computational-theoretical method for material
screening would use the chemical structure as the only input

<figure><img src="image_4.png" /></figure>
(no need for experimental information) and would be able not only to discriminate between miscible and immiscible materials, but also to predict their volume/weight ratios in the different phases. Computational methods to investigate the miscibility for different classes of materials have been widely used across the literature[79-88] but their potential in the study of blends specific for photovoltaics is still largely unexplored.

Here we provide an overview of the theoretical/computational methods that allow the study of the thermodynamic stability of the polymer blends for photovoltaics. The aim is to describe the state of the art with a perspective on the expected contribution of modeling to the development and design of new highly stable organic polymer mixtures. The article is organized as follows: after discussing the fundamental concepts of thermodynamic blend stability and solubility, we briefly discuss some experimental approaches to the study of stability in OPV, and we review the advances in the modeling of polymer blends based on the use of atomistic simulations and multiscale modeling. We conclude with an outlook on the modeling strategies that can have a strong impact on the design of stable blends and the commercial OPV technologies.

## 2. Fundamentals on the Thermodynamics and Stability of Blends: the Interaction Parameter

The thermodynamics of any chemical or physical process is described by the variation in the free energy between the initial and final states

$$ \Delta G = \Delta H - T\Delta S \quad (1) $$

where ΔH and ΔS are the change in enthalpy and entropy of the system, respectively, and T is the temperature. A reaction occurs spontaneously if ΔG is negative and vice versa. For a donor-acceptor blend the relevant process is the mixing during which the two separate components give rise to the blend. The corresponding free energy of mixing ΔG is the quantity that governs the equilibrium conditions of the blend: the blend is stable if the mixing is thermodynamically favored and the free energy change is negative, while for positive values the blend is unstable and it will tend to separate (or it does not even form).

The thermodynamic stability is in principle dependent only on the free energy difference between the initial (separate components) and final configurations (mixed phase). However, the knowledge of the free energy barrier (i.e., the free energy that must be overcome in order to form or to decompose the blend) is necessary to the study of the kinetics of mixing.

Before discussing modern atomistic approaches to the study of stability, it is necessary to start from the older approaches based on continuum theories of the free energy of mixing ΔGm.

Detailed description and derivation of these theories can be found in many book chapters,[89-91] while we recall here the main concepts for ease of reference for the reader.

Continuum models are expressed in terms of empirical or effective parameters, typically as a function of relative concentrations φ or of the chemical potentials μ = ∂G/∂n and the other intensive thermodynamic parameters T and P. Continuum methods have played a fundamental role in understanding the thermodynamics of solutions (e.g., the theory of regular solutions or the

Flory-Huggins extension to macromolecules and polymers) and provide the formalism and the general concepts necessary to describe stability (stability ranges, spinodal, binodal decomposition, etc.).[92-96] A general expression for ΔGm = ΔHm - TΔSm is given by the regular solution model, a model based on the volume discretization and on mean field approximations. In particular it is based on the following approximations. 1) The available space can be discretized into a lattice of molecular volumes (lattice approximation); in such lattice approximation, the available space is divided into lattice sites of volume Vref, often taken equal to the molecular volume of one of the materials composing the blend. 2) The neighboring sites are independent of each other, so if one site is occupied by a molecule of species A, it is neither more nor less likely that a neighboring site is occupied again by a species A molecule. 3) Molecules interact only with their nearest neighbors in a pairwise additive way. 4) The volume of the mix is independent of the composition.

The entropy of mixing is a measure of the number of ways Ω, in which the components A and B can occupy the sites at the desired relative concentration. For a binary mixture, with nA molecules of species A, nB molecules of species B, and N = nA + nB as the total number of molecules in the mix, we have

$$ \Delta S_m = k_B \ln(\Omega) = k_B \ln\left(\frac{N}{n_A}\right) = k_B \ln\left(\frac{N}{n_B}\right) \\ \approx -k_B \left[n_A \ln\frac{n_A}{N} + n_B \ln\frac{n_B}{N}\right] \quad (2) $$

$$ \frac{\Delta S_m}{V_m} = \frac{\Delta S_m}{N V_{\text{ref}}} = - \frac{k_B}{V_{\text{ref}}} [\phi_A \ln(\phi_A) + \phi_B \ln(\phi_B)] $$

where kB is the Boltzmann constant and φA and φB are the molar fractions of species A and B, respectively, that is, φi = nᵢ/(nA+nB). These are equal to the volume fractions if the volume of the mix is independent on the composition.

The entropy is always positive and favors mixing with a maximum for φA = φB = 0.5.

Assumptions 2) and 3) correspond to the mean field approximation and make it possible to derive an expression for the enthalpy in the reference volume ΔHm/Vm[93] of the form

$$ \frac{\Delta H_m}{V_m} = \phi_A \phi_B \frac{\Delta_{AB}}{V_{\text{ref}}} \quad (3) $$

where the parameter ΔAB = z/2(2εAB - εA - εB) contains the difference of the “molecular energies”, that is, the energy of the mixed A-B interaction (εAB) with respect to homologous A-A (εA) and B-B (εB) interactions, and z is the coordination number. The term φAφB = φA(1 - φA) accounts, within mean field approximation, for the number of mixed interactions per site. For example, at φA = 1 or φB = 1 the probability of mixed interactions is zero. The maximum probability is obtained when φA = φB = 0.5.

In summary, within the mean field assumptions described above the effective model of the free energy per unit volume is

$$ \frac{\Delta G_m}{V_m} = \frac{k_B T}{V_{\text{ref}}} [\phi_A \ln(\phi_A) + \phi_B \ln(\phi_B) + \phi_A \phi_B \chi_{AB}] \quad (4) $$

where here φi = Vref/ni / (Vref(nA+nB)) is the volume fraction of the i-th component (we consider that the volume of the mix is independent

<figure><img src="image_5.png" /></figure>
on the composition), and $\chi_{AB}$ is a dimensionless interaction parameter defined as

$$ \chi_{AB} = \frac{1}{k_B T} \frac{z}{2} (2\epsilon_{AB} - \epsilon_A - \epsilon_B) \quad (5) $$

that characterizes the strength and favorability of mixed (A–B) interactions compared with homologous (A–A and B–B) interactions with respect to the thermal budget of a single degree of freedom at the corresponding temperature $k_B T$.

We remark that the entropic part of the free energy (i.e., $k_B T[\phi_A \cdot \ln(\phi_A) + \phi_B \cdot \ln(\phi_B)]$) increases with temperature and always dominates over the enthalpic term $k_B T[\phi_A \phi_B \chi_{AB}] = \frac{z}{2}(2\epsilon_{AB} - \epsilon_A - \epsilon_B)\phi_A \phi_B$ if the energies $\epsilon_{AB}$, $\epsilon_A$, $\epsilon_B$ are substantially temperature independent. At sufficiently high temperature every blend is soluble within this model.

The free energy curve $\Delta G(\phi)$ for a prototypical polymer solution[97] is reported in Figure 1B as a function of the relative concentration $\phi = \phi_A$ (and $\phi_B = 1 - \phi$) for different values of temperatures. For high temperature ($T_1$) the curve has a single minimum and positive concavity, meaning that the blend will not

<figure><img src="image_6.png" /><figcaption>Figure 1. A) Phase diagram and B) free energy curves of a polymer solution. $\phi_0$ represents the initial composition; $\phi'_2$ and $\phi''_2$ represent the equilibrium phase composition at temperature $T_2$; and $\phi'_3$ and $\phi''_3$ represent the equilibrium phase composition at temperature $T_3$. The black dots on free energy curves represent the inflection points. Reproduced with permission.[97] Copyright 2016, Elsevier.</figcaption></figure>

spontaneously phase separate, as any combination of volumes with different relative concentrations would increase the free energy. For $T_2 < T_1$ we find two minima and a local maximum with negative curvature: in this case the blend is unstable, and it will tend to separate into two volumes, each with the relative concentration $\phi'_2$ and $\phi''_2$ of one of the minima; it is easy to understand that this combination can lower the free energy and decomposition is accordingly spontaneous. These two concentrations are coexisting, and the locus of these concentrations as a function of temperature (and thus of the interaction parameter) is known as the coexistence curve or the binodal (see Figure 1A). The limit of local stability is defined by the points, where the curvature of the free energy changes (inflection points in Figure 1B), and the locus of these points is called the spinodal curve.

The above analysis refers to equilibrium conditions where the increase in temperature always favors mixing. However, out of equilibrium, the effect of temperature can be subtler. For example, metastable blends are commonly synthesized which are kinetically but not thermodynamically stable. In these cases (e.g., polymer-fullerene blends) the increase of temperature can favor phase separation (as typically observed during organic solar cells processing).[86,98] This occurs when demixing requires thermally activated phenomena (e.g., diffusion of fullerenes in the polymer matrix).

In the simplest derivation of free energy, $\Delta_{AB}$ is constant and $\chi_{AB} = \frac{\Delta_{AB}}{k_B T} = \frac{\Delta H_m V_{ref}}{k_B T \phi(1-\phi)}$ is only due to enthalpy and scales with temperature as $\chi_{AB} \approx \frac{1}{T}$. In the most general case, the mixing interaction parameter is a function of $\phi$, $T$, $P$, and it can include a dependence on the local packing of components. This is not yet fully understood and all the deviations from the lattice approximation are put together into an entropic contribution to the interaction parameter. This contribution is an extra term with respect to the combinatorial entropy $\Delta S_m$ described above and that is independent on the specific nature of the interactions, being a sole function of size and relative concentration. $\chi_{AB}$ can thus be described as the sum of two contributions: $\chi_H$ which accounts for enthalpic contributions and $\chi_S$ which accounts for entropy, or other energetic effects, arising from the specific chemical nature of the molecules[99]

$$ \chi_{AB} = \chi_H + \chi_S \quad (6) $$

For example, if a molecule is involved in a hydrogen bonding network with molecules that restrict it into certain constrained positions, this will have a configurational entropic penalty due to the reduced number of configurations available.

The Flory–Huggins theory[94,95] generalizes the expression for the free energy of mixing to macromolecular solutions, such as polymer blends. If only small molecules are involved, any lattice site can be occupied either by the solvent or the solute, but when dealing with high-molecular-weight polymers, the value of $\Delta S_m$ is reduced, because of the constraints induced by the polymer (e.g., the sites occupied by the polymer must be consecutive, self-avoiding, etc.). Flory and Huggins introduced the quantities $r_A$ and $r_B$ representing the volumes of polymer A and polymer B (or polymer and molecule) in units of $V_{ref}$ (e.g., if the latter is chosen equal to the molecular volume of one of the two materials, say, A, then $r_A = 1$). They were able to derive an expression

<figure><img src="image_7.png" /></figure>
for the entropy of polymer blends that is identical to that for the regular solutions

$$-\frac{\Delta S_m}{V_m} = \frac{k_B}{V_{\text{ref}}} [\phi_A \ln(\phi_A) + \phi_B \ln(\phi_B)] \quad (7)$$

provided that the relative concentrations depend on the molecule volumes as

$$\phi_i = \frac{V_{\text{ref}} r_i n_i}{V_{\text{ref}} (r_A n_A + r_B n_B)} \quad (8)$$

In the general case the sizes of the polymers are different (i.e., $r_A \neq r_B$), and the above FH mixing entropy is an asymmetric function of $\phi$. The derivation of Equation (3) for the enthalpy is still valid also for macromolecules, so that Equation (3) and (7) give the formula for the effective FH free energy of polymeric solutions. Further details can be found in classical texts on thermodynamics of the blends.[92,93,100]

It is useful here to relate the interaction parameter to a macroscopic quantity, that is, the enthalpy change upon mixing per unit volume $\frac{\Delta H_m}{V_m}$. The latter quantity can be measured or calculated by atomistic simulations and it can be expressed as the difference between the energy of the mixed and of the unmixed states, divided by the volume of the mix $V_m$.

$$\frac{\Delta H_m}{V_m} = \frac{E_m - E_A - E_B}{V_m} \quad (9)$$

where $E_A$ ($E_B$) is the energy of the same number of particles of the component A (B) present in the mix but obtained in the pure single-component phase. To further simplify we suppose that the volume change upon mixing is negligible (which is true for most blends), so that we can write

$$\frac{\Delta H_m}{V_m} = \frac{E_m}{V_m} - \frac{E_A}{V_A} \frac{V_A}{V_m} - \frac{E_B}{V_B} \frac{V_B}{V_m} = \frac{E_m}{V_m} - \phi_A \frac{E_A}{V_A} - \phi_B \frac{E_B}{V_B} \quad (10)$$

In conclusion $\Delta H_m/V_m$, that is, the mixing enthalpy per unit volume, can thus be related to the cohesive energy densities (i.e., the energy per unit volume needed to vaporize the molecules from the condensed phase) of the i components $CED_i = E_i/V_i$ and of the mix $CED_m = E_m/V_m$ by

$$\frac{\Delta H_m}{V_m} = CED_m - \phi_A CED_A - \phi_B CED_B \quad (11)$$

and, if the entropic contribution to $\chi_{AB}$ is negligible

$$\chi_{AB} = \chi_H + \chi_S \approx \chi_H = \frac{CED_m - \phi_A CED_A - \phi_B CED_B}{k_B T \phi_A \phi_B} \quad (12)$$

It has been evidenced in literature that care should be taken when considering $\chi_{AB} \approx \chi_H$, especially when hydrogen bonds are present in the blend, which may constrain the macromolecules in specific conformations, thus lowering the entropy of the system.[101]

Provided that parameter $\chi_{AB}(T, \phi)$ is known, the free energy of a blend can be calculated giving information on stability, spinodal decomposition, binodal loci, phase behaviors etc.

This approach is discussed for example in the study by Bates et al.[102] and references therein.

In the next sections we discuss how to derive the interaction parameter from experimental data.

### 3. Experimental Determination of the Interaction Parameter

Several studies have tried to determine the interaction parameters for a number of donor/acceptor blends of interest for OPV and establish a correlation between the overall efficiency and the purity of the mixed phase. This has been mainly done to estimate the miscibility of donor polymers with molecular acceptors (fullerene or NFA) or with other polymers (acting, e.g., as sensitizers) in ternary and quaternary devices.[103–105]

A possible way to obtain $\chi$ is through the melting point depression method.[106–108] Using differential scanning calorimetry (DSC) to measure the change in melting temperature of the polymer in the mix ($T_m$) as a function of the acceptor concentration $\phi$, one can estimate the value of the interaction parameter via the relation

$$\frac{1}{T_m} - \frac{1}{T_{m0}} = \frac{R}{\Delta H_f} \frac{V_m}{V_s} (\phi - \chi \phi^2) \quad (13)$$

where $T_m$ is the melting point of the polymer at acceptor volume fraction $\phi$, $T_{m0}$ is the melting point of pure polymer, R is the ideal gas constant, $\Delta H_f$ is the heat of fusion of the polymer, $V_m$ is the monomer molar volume of polymer, and $V_s$ is the acceptor molar volume.

Kozub et al.[109] were among the first groups to study the miscibility of P3HT and PCBM and estimated a value of $\chi$ concluding that miscibility is possible for polymer volume fractions greater than 0.42. Determination of the miscibility is also important when optimizing the choice of the third component in the realization of ternary solar cells. Addition of a sensitizer to the widely studied, highly efficient polymer:NFA-mixed PBDB-T:ITIC has been studied by Yi et al.[110] In their work, PDCBT and P3HT, which have similar chemical structures but different side chains, were chosen as third components. The Flory-Huggins interaction parameters between the sensitizers and both the polymer PBDB-T and the ITIC molecule acceptor were determined via the melting point depression method, and it was found that while both sensitizing polymers are able to form bimolecular crystals with PBDB-T, they possess different miscibility with ITIC, leading to very different morphologies and photovoltaic performances.

Ye et al.[111] developed an approach to determine the dependence on the processing temperature of the FH parameter of polymer-small molecules blends and to establish quantitative relations between $\chi$ and device performance. These authors were among the first to point out the relevance of miscibility in the achievement of good performance in polymer:NFA blends.[112,113] By measuring the polymer concentration $\phi$ in the mixed phase (obtained at different processing temperatures), it is possible to determine the $T - \phi$ phase diagram (i.e., the binodal curve). For each value of temperature, the two values of $\phi$ correspond to minima in the $G_m$, that is, where $\frac{\partial G_m}{\partial \phi} = 0$

<figure><img src="image_8.png" /></figure>
and $\frac{\partial^2 G_m}{\partial\phi^2} > 0$. From these two conditions, it is possible to extract the value of $\chi$ at different temperatures and to obtain $\chi(T)$ by fitting the general formula $\chi = A + \frac{B}{T}$. The authors have then measured the photovoltaic performances of the devices fabricated at different processing temperatures and plotted the fill factor FF versus $\chi$, showing that there is a quantitative relation between the value of the fill factor and the interaction parameter. In particular, they were able to predict the optimum $\chi$ and thus composition leading to the highest device FF for a wide range of binary blends composed of small-molecule acceptors and (amorphous and semicrystalline) polymer donors.

The same group[98] focused on the impact of morphological instabilities on the so-called burn-in (i.e., the efficiency loss in the first few hundred hours of operation of an unstable organic solar cell) in polymer-NFA OPVs, which is typically driven by polymer-NFA demixing in the mixed domains and/or crystallization of the NFA. The crystallization of NFAs is inhibited by quenching the acceptor in an amorphous state which is metastable: the thermodynamically stable state in fact involves NFA crystals. As the metastable state is governed by the binodal composition, if the optimal morphology is quenched far from the binodal, severe burn-in is expected (see Figure 2a). In contrast, a device with a binodal close to the optimal morphology (see Figure 2b) is expected to be relatively stable and thus to exhibit lower burn-in degradation. The authors calculated the FH interaction parameter as in the study by Ye et al,[111] determining the amorphous-amorphous phase diagram (and thus the binodal) at different processing conditions for four polymer-NFA blends (P3HT:EH-IDTBR, FTAZ:EH-IDTBR, FTAZ:ITIC, and PTB7-Th:EH-IDTBR) and delineating thermodynamic drivers and kinetic factors for stability. In particular, they estimated the diffusion coefficients of EH-IDTBR in P3HT, FTAZ, and PTB7-Th, showing that in some blends it is sufficiently high to enable

demixing and crystallization burn-in after no- or low-temperature annealing. They concluded that the optimum morphology must be stabilized by a high degree of vitrification with low diffusion, or it has to be metastable with a mixed composition near the percolation threshold.

## 4. Computational Approaches to the Modeling of Stability of (organic) Blends

Despite the great advances in the experimental techniques to access the nanomorphology of BHJ and to control to some extent the stability of the blends, a truly predictive approach is not currently available and, in principle, only atomistic methods can achieve this goal.

The possibility to predict thermodynamic properties of novel blends and to perform high-throughput screening via in silico methods would be a major step forward in the search for optimum, highly stable blends not only in the field of organic PV but in a wide range of applications, from materials science to pharmaceutics and drug delivery.

To predict the miscibility through computational methods, two strategies are widely used: the first is the determination and comparison of the solubility parameters $\delta$ of only the separate phases of the blend (the so-called Hildebrand and Hansen solubility parameters, HSPs, discussed below) and the other is the simulation of both the mixed phase and the separate pure phases and the determination of the interaction parameter $\chi$.

The solubility parameter approach is based on the idea that materials with similar $\delta$ should be miscible. This approach is derived from an approximation of Equation (5) assuming that the interaction energy $\epsilon_{AB}$ can be approximated as, $\sqrt{\epsilon_A\epsilon_B}$ that is, the geometric average of $\epsilon_A$ and $\epsilon_B$ obtaining

<figure><img src="image_9.png" /><figcaption>**Figure 2.** a) A blend with miscibility that is too low (hypomiscible). Solid lines, including the blue part of the y-axis, correspond to thermodynamic equilibrium, whereas dashed lines are only metastable states. Below $T_g$, which is composition and material dependent, the liquids are frozen into a glass. Equilibrium can still be established by diffusion but is kinetically hindered. The points A, P, and B represent initial average D/A ratio, percolation threshold of the acceptor in the mixed domains (optimal morphology), and the binodal composition at a given temperature, respectively. Point C represents the crystallization of NFAs in a blend. The system generally proceeds from point A to point C during casting and aging. Processing conditions are typically chosen to reach the vicinity of point P, where an optimal tradeoff between charge creation, charge extraction, and charge recombination is achieved if the domain size can be sufficiently controlled. b) A blend of an optimum amorphous-amorphous miscibility. Adapted with permission.[98]
Copyright 2019, Elsevier.</figcaption></figure>

<figure><img src="image_10.png" /></figure>
$$ \chi_{AB} = \frac{1}{k_B T} \frac{z}{2} (2\sqrt{\epsilon_A \epsilon_B} - \epsilon_A - \epsilon_B) \quad (14) $$

By defining the solubility parameter for species i = A, B as $\delta_i = \sqrt{\frac{z \epsilon_i}{2 V_{ref}}} = \sqrt{CED_i}$ from which $\epsilon_i = \frac{z}{2} V_{ref} \delta_i^2$, an approximated expression for $\chi_{AB}$ is finally obtained

$$ \chi_{AB} = \frac{1}{k_B T} V_{\text{ref}} (2\delta_A \delta_B - \delta_A^2 - \delta_B^2) = \frac{V_{\text{ref}}}{k_B T} (\delta_A - \delta_B)^2 \quad (15) $$

that only depends on the difference $\delta_A - \delta_B$. According to Equation (15), the two materials are miscible if their solubility parameters are similar, as in this case $\chi_{AB}$ is small.[114,115] According to Equation (15) the value $\chi_{AB}$ is always positive, and so this approximation applies only if there is an energy cost for the formation of the blends. The theory was developed for nonpolar liquids, and it was extended by Hansen[116] by decomposing $\delta$ into three components, that is, polar, nonpolar, and hydrogen bonding. The solubility parameter can be calculated either via group contribution (GC) methods, where the value of the parameter is given by the sum of the contributions coming from the different parts that make up the molecule, or via molecular dynamics simulations, where the CED is calculated from the full atomic structure. GC methods were first proposed by Small,[117] are relatively easy to use, and have been applied to determine $\delta$ in a few works.[118–120]

The determination of $\delta$ via all-atoms simulations is expected to be more precise, as it takes into account all microscopic interactions between the materials. Nevertheless, the approximations underlying solubility parameter often fail in predicting miscibility when factors such as temperature, concentration, viscosity, ionic interactions are important, such as for example in predicting polymer-drug miscibility for drug delivery applications.[121]

To overcome the above limits, the second strategy is to directly calculate, though at a much higher computational effort, the interaction parameter $\chi$ by the atomistic simulation of both the blends and the separate components.

This has been carried out by different groups in a wide range of applications, such as the prediction of the mixing of impurities with crude oil,[81,82] the study of drug-polymer formulations for pharmaceutical use,[80,122–125] and the study of the solubility of conjugated polymers for organic electronics[83,84] and the thermodynamics of mixing in organic blends for solar cells.[86–88]

The field of drug delivery is the one where the method has found most of its applications: as many newly developed drugs possess poor water solubility, delivery of these drugs within the human body is a major and challenging problem in pharmaceutical technology. Among the approaches to solve this problem there is the incorporation of drugs into nanoparticles, liposomes, emulsions, or polymers. Despite the high number of works on the subject, in each of them only a limited chemical space has been considered, and consequently none of the available studies are able to draw definitive conclusions about the usefulness of these techniques in this field. Turpin et al.[101] in particular conducted a systematic review of works on the in silico determination of drug-polymer $\chi$ for solid dispersions. They have also conducted a comparative study between experiments and different theoretical techniques (namely, determination of $\delta$ via GC or

MD and direct determination of $\chi$ via MD) and have concluded that these theoretical methods have a limited reliability to predict drug-polymer interactions for the purposes of screening or excipient selection. Among the reasons for these deviations there is the possible limited transferability of the available force fields and the incorrect evaluation of the entropic part of $\chi$, which for materials where hydrogen bonds are important can make a large contribution and should not be neglected. In these cases, the calculation of the free energy of mixing and not only the enthalpy should be taken into account. There is extensive literature on the calculation of the free energy difference that relies on thermodynamic integration, which consists of defining a thermodynamic path between initial and final states and integrating over ensemble-averaged enthalpy changes along the path.[126,127] Several methods have been developed such as, for example, umbrella sampling,[128] metadynamics,[129] and free energy perturbation,[130] which have been implemented in most molecular dynamics simulation codes and successfully applied, mainly in the field of biophysics or in other contexts[131–138] but to a less extent to polymer blends for photovoltaics.

The effect of the entropic part of $\chi$ can be minor in systems where hydrogen bonds are less important, such as organic blends for photovoltaics. In these systems, in fact, both the polymers and the small molecules possess large conjugated cores, and large contributions to the intermixing come from by $\pi-\pi$ stacking.[139–142]

Recent works based on MD[86,88] derived the interaction parameters of blends formed by low-bandgap polymers such as PTB7, PTB7-th, and Si-PCPDTBT with fullerene acceptors PC60BM and PC70BM and found good agreement with experiments. Interestingly, it was found that the mixing enthalpy dependence on $\phi$ is not the simple parabolic function $\phi(1-\phi)$ of the mean field approximation and the interaction parameter is not independent on $\phi$ as expected from the ideal Flory-Huggins theory. The dependence of $\chi$ on concentration or on the blend density is typically not taken into account in literature.

A central point for the accurate calculation of the interaction parameter is the quality of the blend model. It is necessary to have realistic microstructures and suitable sampling of the configurational space. In the studies discussed above,[86,88] the atomistic models of polymer blends for photovoltaics were generated by randomly distributing fullerenes within a dilated mesh of polymers chains and letting the system relax by long constant temperature and pressure MD runs under the effect of cohesive forces due to interatomic potentials.[143] For each concentration several initial configurations were considered. Alternative strategies could be to generate models of blends starting from the donor-acceptor molecules in solution and simulating the solvent evaporation. This mimics the experimental process but it is computationally demanding due to the large number of atoms and long simulation time needed. For example, in the study by Gertsen et al.,[144] the solution deposition of NFA thin films via solvent evaporation was simulated, and the resulting structures were compared with experimental results from grazing-incidence X-Ray scattering measurements. Short-range structural properties, such as $\pi$-stacking and alignment effects relative to the substrate, were accurately reproduced, but the time scales affordable were not sufficient to induce significant longer-range order in the thin films. Applications that make use of

<figure><img src="image_11.png" /></figure>
coarse-graining (CG) strategies are needed to reduce the work-
load and extend the scales.[145,146]

The MD results based on the random approach mentioned earlier have shown to provide meaningful morphologies. The value of ΔH_m was found to be dependent on the relative polymer/fullerene concentration with local oscillations (see Figure 3a) and minima which corresponded to the concentrations determined experimentally by scanning transmission electron microscopy combined with spatially resolved spectroscopic imaging (STEM-SI) of low-energy loss (Figure 3c-e). This is a powerful approach able to distinguish between donor and acceptor phases of the polymer blend with nanoscale resolution.[86,147,148]

The analysis of MD blends can also provide information on the
dependence of ΔHm from blend density (see Figure 3b). Indeed,
as polymer cells are processed from solution, fast drying of the
blend can generate nanoscale morphology in nonequilibrium con-
ditions. For fullerene/PTBT/PTB7-th blends it was found that at
low density the polymer chains are not sizably bent, the binding to
fullerene is dominated by the large van der Waals interactions, and
the PTB7-th polymer has a better adhesion to the fullerene than
PTB7, due to the larger conjugated core. As the density of the sys-
tems increases, the polymer chains are constrained by neighbor-
ing molecules and the energy cost of bending increases rapidly for
the more rigid PTB7-th so that the rather flexible PTB7 polymer is
energetically favored and has better miscibility (lower) ΔHm with
fullerenes. From atomistic data, qualitative models to correlate the
stability to the commensurability between polymers and fullerenes
size can be tested and extracted.[88]

The above models and analysis are valid under the condition
that enthalpy dominates over mixing entropy. If this is not the
case, the free energy can in principle be calculated directly by
atomistic methods, making use of suitable thermodynamic inte-
gration along a configurational path connecting separated phases
to mixed phases. However, such an all-atom approach for 3D
blends of macromolecules is very demanding and the direct cal-
culation of free energy is still unexplored for these materials.
The application of such methods to study mixing and solubility
has mainly focused on small molecules or ions such as mixing
water with acetone,[149] small alcohols,[150] and DMSO,[151] to
study the solubility of ions in electrolytes[152] and to model ther-
modynamic mixing properties of solid solutions.[153,154]

To further extend the time and length scales of molecular
dynamics simulations, coarse-grained approaches are
possible in which groups of atoms are replaced by beads with
effective masses and interactions able to reproduce the
average material property but at a much lower computational
cost.[145,146,155–158] However the models developed so far are not
readily transferable, limiting their potential for high-throughput
screening calculations. Recently, Marrink et al.[146,159] proposed a
method based on the transferable Martini CG force field[160,161] to
simulate the morphological evolution during solution processing
of P3HT:PCBM blends. This force field still possesses a degree of
chemical specificity, allowing its use also for the design of com-
pounds which differ only slightly in chemical nature. The force
field was initially developed for lipids, and parameters exist for
water, simple alkanes, organic solvents, surfactants, a wide range
of lipids, cholesterol, proteins,[162] carbohydrates,[163] nucleic

<figure><img src="image_12.png" /><figcaption>**Figure 3.** a) Enthalpy of mixing per unit volume ΔH<sub>m</sub>/V<sub>m</sub> calculated for PTB7:PC60BM and PTB7:PC70BM as a function of the fullerene weight ratio. Error bars on each point are the standard deviation of the ΔH<sub>m</sub>/V<sub>m</sub> calculated on ten systems with different starting configurations. b) ΔH<sub>m</sub>/V<sub>m</sub> at 1:1 polymer: PC70BM weight ratio, calculated for polymers PTB7, PTB7-th, and Si-PCPDTBT as a function of the blend density. c) Conventional annular dark-field STEM micrograph of the as-cast PTB7:PC70BM photoactive layer. d) Plasmon peak map of the same area of observation in (c) showing the nanoscale materials phase (PC70BM rich in red with high plasmon peak energy and PTB7 rich in green with low energy). e) Segmentation of the same area into enriched and mixed domains. a), c–e) Adapted with permission.[86] Copyright 2016, Wiley-VCH. b) Adapted with permission.[88] Copyright 2021, Elsevier.</figcaption></figure>

<figure><img src="image_13.png" /></figure>
acids,[164,165] some carbon nanoparticles,[156] ionic liquids,[166]
and some polymers.[146,167,168]

Among coarse-graining schemes, it is worth citing dissipative
particle dynamics (DPD), which is a method developed to address
mesoscale problems in complex fluids and soft matter in general.
The method allows access to polymer structure and morphology
at a reasonable computational cost, and it has been very successful
in identifying mechanisms in phase separation.[169] DPD simula-
tions have been extensively used to study polymers in solvents but
also polymer melts and blends. Concerning OPVs, the application
of DPD to explore the active layer morphology is limited, includ-
ing P3HT:PCBM[170] and small-molecule OSCs in the presence of
solvents and additives.[171]

Finally, the morphology evolution of the full active layer can be
simulated by continuum approaches, such as phase field (PF)
models.[172-175] These are typically used to model the kinetics of
phase separation of donor:acceptor blends by the Cahn-Hilliard
equation,[176,177] which evaluates the change of local composition
with annealing time. PF models require the knowledge of a
number of effective parameters of the blend, including the
Flory-Huggins parameters and donor/acceptor/solvent interface
energies, which can be obtained from experiments or by
atomistic/CG simulations as described earlier. Applications of
PF methods to OPV include the simulation of the morphology
evolution during solvent evaporation[178] and thermal
annealing.[179] The computed morphologies can be then corre-
lated with optoelectronic and photovoltaic properties.[179,180]

One aspect that is crucial in the calculation of the FH interac-
tion parameter is the accuracy of the force field adopted in
MD/CG simulations. The possibility to adopt more accurate
atomistic calculations at the level of ab initio DFT[181,182] has been

explored within a multiscale strategy. High-quality quantum
chemical calculations have been used for predicting HSPs for
rather simple, low-dimensional molecules.[183] An interesting
approach is represented by hybrid quantum mechanics/molecu-
lar mechanics (QM/MM) simulations: within this scheme, a QM
zone (which is simulated ab initio) is surrounded by a classical
environment (MM zone). This approach is particularly useful
whenever an explicit description of changes in the electronic
structure is necessary as, for example, when modeling charge
transport.[184,185] Obviously, the computational bottleneck of
QM/MM methods is the description of the QM part, and semiem-
pirical methods are often used to reduce the computational
workload.[186–188] A recently proposed alternative approach is rep-
resented by substituting the QM part by machine-learnt models,
but this is challenging in condensed-phase systems due to long-
range interactions and further progress is needed toward this direc-
tion.[189,190] The major obstacle in the application of QM/MM to the
study of polymer blends is that intermixing derive in general from
long-range interactions (dispersive and multipolar) and the linking
sites are distributed within the blend volume. The separation of the
QM and MM zones under these conditions is ill-defined.

In Figure 4 we report a summary of atomistic, coarse-grained, and continuum models of polymer blends that can be used to calculate the free energy of mixing at different levels of approximation. When accurate ab initio methods are employed (downward arrow of Figure 4), their higher computational cost hampers the study of a realistic blend microstructure: the rough- est approximation is to model the blend as a pair formed by one donor and one acceptor molecule or oligomer and use their inter- action in combination with the mean field approximations to esti- mate ΔHm. In the upper direction (upward arrow of Figure 5),

<figure><img src="image_14.png" /></figure>

**Figure 4.** Schematic summary of models (and corresponding methods) that can be used to calculate the thermodynamics of polymer blends for photo-voltlaics. Time and size (logarithmic) scales accessible for different methods are reported on the left. Light blue shadings in the axis represent intermediate scales shared by different models. Ab initio methods (e.g., DFT) that reproduce accurately the molecular interactions that are limited to small systems with oversimplified morphologies (downward arrow). In contrast, blend models with realistic morphologies (upward arrows) can be treated by classical potentials, coarse-grained, and continuum methods but at the price of less accurate description of atomic interactions. The exact free energy requires both accurate interactions and realistic models. Images for DFT: Adapted with permission.[21] Copyright 2020, American Chemical Society. Adapted with permission.[191] Copyright 2018, Wiley-VCH. Image for QM/MM: Reproduced with permission.[184] Copyright 2014, American Chemical Society. Image for MPMD: Reproduced with permission.[88] Copyright 2020, Elsevier Ltd. Image for CG: Adapted with permission.[192] Copyright 2019, Wiley-VCH. Image for DPD: Reproduced with permission under the terms of the Creative Commons CC BY license.[170] Copyright 2015, the Authors. Published by Springer Nature. Image for PF: Reproduced with permission.[180] Copyright 2022, Elsevier.

<figure><img src="image_15.png" /></figure>
model potentials and coarse-grain methods allow to simulate much bigger and realistic portions of the active layer, also allowing, in principle, the determination of the full $\Delta G_m$ but at the cost of less accurate interatomic interactions. A good approximation of the exact free energy requires both accurate atomistic interactions and large-scale realistic models, and this could be obtained by multiscale approaches.

Novel opportunities have been opened by ML techniques. These have been applied mainly to predict PCE and are typically trained on optical properties (such as, e.g., energy levels) but information on solubility and mixing properties is rarely included in the descriptors.

In a recent work, ML coupled with MD simulations has been also employed to design new materials for OPV[71] with improved optoelectronic properties and stability. More than 5000 NFAs were designed (starting from synthesizable building blocks) and screened to be used together with the low-bandgap polymer PTB7-th as a donor. ML models were trained on data collected from the literature and various properties were predicted for the new NFAs, such as energy levels, UV/vis absorption maxima, and PCE. More than 100 NFAs with suitable energy-level alignment, absorption maxima, and PCE > 13% were selected.

<figure><img src="image_16.png" /><figcaption>**Figure 5.** Computational flowchart describing the routine for determining the relative stability capable of describing the microstructure of polymer: fullerene blends from the study by Perea et al.[194] (i) creation of the $\sigma$-profile from the conductor-like screening model (COSMO); (ii) $\sigma$-moments as extracted from COSMO are fed into an ANN to determine HSPs; (iii) HSPs are used to calculate the qualitative FH interaction parameters ($\chi_{1,2}$); (iv) determination of reduced materials properties; (v) spinodal demixing diagrams resulting from polymer blend theory; and (vi) FoM defined as the ratio of the FH intermolecular parameter and spinodal diagram forms the basis of a relative stability metric. Reproduced with permission.[194] Copyright 2017, American Chemical Society.</figcaption></figure>

The selected NFAs were studied using molecular dynamics simulations to predict the Flory-Huggins parameter of the blends. This was obtained by calculating the binding energies and coordination numbers of polymer-polymer, NFA-NFA, and polymer-NFA pairs. Equation (5) was then applied to obtain $\chi$. Finally, 15 NFAs were selected, and the best predicted PCE was found to be over 15%, which is far better than reported results on PTB7-th-based blends.

Another example is the work of Perea et al.[193] who developed a numerical approach to determine the solubility parameters of fullerene mixtures using a mathematical tool based on artificial neural networks (ANNs). They calculated the molecular surface charge density distribution ($\sigma$-profile) by DFT within the framework of a continuum solvation model and then applied ANN to transform it into solubility parameters. The corresponding computational flowchart is reported in Figure 5. The calculated values were validated against experimentally determined HSP, showing excellent agreement. The same authors also conceived a figure of merit (FoM) to describe the mixing thermodynamics (phase evolution) of polymer:fullerene bulk heterojunctions. This new model predicts relative stabilities linking the FH interaction parameter to the spinodal demixing interaction parameter (also referred to as the spinodal line).[194]

The studies reported above are representative of the state of art in the modeling of the thermodynamic stability of polymer blends. To place such studies in a broader context we report in Table 1 a summary of the recent computational works (either atomistic or based on ML analysis of data) on OPV. We report in columns the material investigated, the methodology, and the type of blend model with the corresponding calculated properties.

Most studies focus on the prediction of electronic energy levels or PCE, while stability is taken into account in a few cases. It is clear that the available advanced computational techniques for the modeling of polymer solar cells are still largely underutilized for the study of stability.

As shown, atomistic simulations can be used to determine the parameters necessary for the continuum description of the thermodynamics of mixing and the determination of most stable configurations. However, real systems are not limited to thermodynamically stable mixtures: kinetically stable blends could give rise to de facto stable devices, if demixing happens on a sufficiently long timescale, and should be thus included in the potentially optimum blends for OPV. As such, other attributes of the mix may be important, which cannot be captured in the effective FH or solubility parameters but are directly accessible, in principle, by atomistic models.

## 5. Outlook and Perspectives

Stability remains a fundamental issue in the quest for novel donor:acceptor blends for OPVs. Due to the very large chemical space available for the synthesis of novel donor and acceptor moieties, it is necessary to go beyond the experimental trial-and-error approach. An improvement in this direction can be the extensive use of artificial intelligence and ML techniques to predict the stability of novel blends from available experimental data on existing materials.[77,78] This approach has the advantage of being readily available and relatively fast, but it is not easily generalizable

<figure><img src="image_17.png" /></figure>
<table><caption>**Table 1.** Synoptic table of recent computational works on organic photovoltaics reporting: the investigated material, the method, the type of donor-acceptor model and the calculated properties (thermodynamic, electronic and photovoltaic).</caption><thead><tr><th rowspan="2">Ref</th><th rowspan="2">Material</th><th rowspan="2">Methods</th><th rowspan="2">D/A Model</th><th colspan="3">Calculated Properties</th></tr><tr><th>thermodyn properties</th><th>electronic properties</th><th>PV properties</th></tr></thead><tbody><tr><td>[86]</td><td>polymer/fullerene</td><td>MPMD</td><td>3D blend (thermodyn)</td><td>ΔH<sub>mr</sub>, χ</td><td></td><td></td></tr><tr><td>[87]</td><td>polymer/fullerene</td><td>MPMD</td><td>3D blend (thermodyn)</td><td>ΔH<sub>mr</sub>, χ</td><td></td><td></td></tr><tr><td>[88]</td><td>polymer/fullerene</td><td>MPMD</td><td>3D blend (thermodyn)</td><td>ΔH<sub>mr</sub>, χ</td><td></td><td></td></tr><tr><td>[193]</td><td>fullerenes</td><td>DFT + ML</td><td>A in implicit solvent</td><td>δ</td><td></td><td></td></tr><tr><td>[194]</td><td>polymer/fullerene</td><td>DFT + ML</td><td>D, A molecules in implicit solvent</td><td>χ, δ</td><td></td><td></td></tr><tr><td>[77]</td><td>polymer/fullerene</td><td>ML on existing data</td><td></td><td></td><td></td><td>PCE</td></tr><tr><td>[195]</td><td>polymer/NFA</td><td>DFT + ML</td><td>D/A molecular pairs</td><td></td><td>HOMO, LUMO, abs. spectra, k<sub>CS</sub>, k<sub>CR</sub> (DFT)</td><td>PCE</td></tr><tr><td>[146]</td><td>polymer/fullerene</td><td>CG</td><td>3D blend</td><td>blend morphology</td><td></td><td></td></tr><tr><td>[20]</td><td>polymer/NFA</td><td>EXP + DFT + MPMD</td><td>3D blend</td><td>blend morphology</td><td></td><td></td></tr><tr><td>[78]</td><td>polymer/NFA</td><td>ML on existing data + MPMD</td><td></td><td>χ</td><td>HOMO LUMO absorption maxima (ML)</td><td>PCE</td></tr><tr><td>[196]</td><td>polymer/NFA</td><td>EXP + ML</td><td></td><td></td><td></td><td>PCE</td></tr><tr><td>[197]</td><td>polymer/NFA</td><td>ML on existing data</td><td></td><td></td><td></td><td>PCE</td></tr><tr><td>[198]</td><td>polymer/NFA</td><td>ML on existing data</td><td></td><td></td><td></td><td>PCE</td></tr><tr><td>[199]</td><td>Small molecules (D/A)</td><td>ML on existing data</td><td></td><td></td><td></td><td>PCE</td></tr><tr><td>[200]</td><td>NFA</td><td>ML on existing data + DFT</td><td></td><td>Log(P<sub>o/w</sub>)</td><td>HOMO, LUMO, f<sub>max</sub> (DFT)</td><td></td></tr><tr><td>[201]</td><td>donors</td><td>ML + DFT</td><td></td><td></td><td>HOMO, LUMO (DFT)</td><td>PCE</td></tr><tr><td>[170]</td><td>polymer/fullerene</td><td>MPMD + DPD</td><td>3D blend + solvent + additive</td><td>δ, χ, blend morphology</td><td></td><td></td></tr></tbody></table>

whenever important changes on the molecules involved are committed and it is not easy to extract physicochemical rules for improving the results. The alternative approach is to use atomistic simulations to develop predictive and accurate models of the free energy of mixing and to make use of high-throughput calculations to reduce part of the experimental screening but for the final stages of materials selection.

Provided that fast and reliable models are available, theoretical screening would be possible, saving time and money with respect to simple trial and error.[202] Nowadays, the requested level of accuracy is provided by ab initio methods. In particular DFT calculations are the best choice among first-principles methods as a good compromise between accuracy and computational cost.[203–205] However, ab initio methods are still too computationally expensive for the study of stability through the use of realistic (i.e., large-scale) atomistic models. The perspective of ab initio-based predictive models of free energy of blends requires necessarily a multiparadigm/multiphysics/multiscale approach. We envisage at least two directions that are still ongoing and that can provide important support to experiments: 1) predictions based on ML techniques trained on ab initio data calculated on representative small-scale models and 2) development of more accurate multiscale methods (e.g., based on a hierarchy of coarse-graining steps with minimum human effort) for the generation of atomistic models of blends and the calculation of corresponding free energy.

The first strategy requires intensive use of ab initio methods to perform calculations on simplified models that still retain a

minimum amount of information on the donor-acceptor and mixing interactions (e.g., at least the most relevant chemical groups of the donor-acceptor mixtures and binding). By suitable standardization of the system generation and calculated properties (it is necessary to automatically explore the compositional and configurational space of donor:acceptors), it is possible to perform high-throughput calculations organized into clusters of embarrassingly parallel calculations on independent materials.[206–208] ML techniques with suitable training of the algorithms can be used to predict stability and to drive additional data generation. The efficiency and the role of strategies based on artificial intelligence are expected to grow rapidly in next years, also considering likely improvements of ML algorithms, which is a very active rising research field.[195,209–214] Furthermore, several freely available databases have been made easily accessible through web interfaces, such as AFLOW, Materials Project, NOMAD, Materials Cloud.[215–218] Nevertheless, when dealing with stability, the complexity of the blends and entropic effects plays an important role and it is difficult to include them in the atomistic models affordable by only ab initio methods. Furthermore, ML predictions do not necessarily reveal the key physicochemical factors or provide rules able to drive the synthesis of improved stable blends. We believe that ML will play a major role in the development of OPVs, but for more accurate predictions not only electronic properties should be included within the descriptors, but also morphology-related properties (e.g., intermolecular packing and mixing, extent of phase separation, crystallinity, size and purity of domains, χ parameter) are needed.[219]

<figure><img src="image_18.png" /></figure>
The second path for the study and in silico design of polymer blends is the use of multiscale models, which can in principle simulate realistic donor:acceptor blends. As discussed earlier, complex blends can be simulated by molecular dynamics and the corresponding mixing enthalpy calculated. Simulating atom-istic models of interfaces, both crystalline and amorphous, under experimental conditions (temperature, pressure, stress, presence of solvents...) can give important information on the blend morphology and density, as well as enthalpy of mixing.[86,88] The relevant open challenges within these methods are related to the limited accuracy and transferability of the interatomic force fields (which are required to be possibly at ab initio accuracy levels) and to the need to extend the length and time scales to allow including more realistic microstructures and thermodynamics with better sampling of phase space. Concerning accuracy, model potentials trained on ab initio through ML (which promise to reach the accuracy of quantum mechanical computations at a substantially reduced computational cost) constitute a promising approach to solve large-scale problems in materials science.[220-224] Further progress is necessary, as currently available models are still computationally more demanding than typical classical MD simulations.[225]

The possibility to simulate the synthesis/mixing processes to gain access to entropic effects is also intriguing. In this direction, coarse-grained molecular dynamics simulations have been used to simulate organic materials for OPVs.[145,156-158] However, the atom-to-bead mapping and the determination of the topology is in general not automatic, and the coarsening process typically has to be carried out manually by expert researchers. Programs which facilitate automatic mapping for CG simulations exist for biomolecular simulations[226-228] while in most cases they are missing for other systems including blends for photovoltaics, with few exceptions.[229-232] Interestingly, ML

methods have also been applied to parameterize bottom-up coarse-grained force fields.[233-235]

Finally, developing methods that allow a smooth transition from discrete (i.e., all atom or CG) to continuum (e.g., PF[172]) models with minimum information losses would be highly desirable to better link atomistics to experimental scales. As a final remark, we point out that important advances are needed to allow the simulation of more complex structures which are already being synthesized, such as, for example, ternary blends, or to include the effects of charge transport layers into the stability of organic solar cells.[236,237] To our knowledge, a full simulation of a ternary blend has not been carried out yet, though theoretical insight has been obtained from atomistic simulations of binary blends.[88]

We summarize in **Figure 6** the schematic roadmap with the two modeling strategies discussed earlier that are expected to play an important role in predicting materials with optimal stability. Starting from the chemical structures *x* of the blend constituents, the stability properties *P*(*x*) (e.g., the interaction parameter *χ* or the free energies) are calculated by either one or the other strategy: multiscale models with improved accuracy (left path) or ML methods applied to large ab initio datasets (right arrow), as well as a judicious combination of the two strategies (horizontal arrow). Both are rooted necessarily on predictive ab initio calculations or experiments. It seems very likely that in the near future the two strategies will reach a level of reliability and computational efficiency to allow routinely calculations of the structure–property relation of polymer blends, that is, *P*(*x*) (blue box).

Provided that this goal is achieved, the design of the blends $x_{opt}$ that optimizes the stability $P(x)$ can be obtained by a suitable optimization algorithm in which the chemical structure is iteratively varied ($dx$) and the new $P(x)$ calculated until a stationary point is obtained ($dP(x)/dx = 0$).

<figure><img src="image_19.png" /><figcaption>**Figure 6.** Roadmap of the possible paths to achieve the chemical structure–property relation $P(x)$ (arrows in blue rectangle) and predict the best stability (optimization/design right arrow) in complex blends for OPVs. Multiscale models (MD/CG) can be derived from experimental data and DFT (also with the use of ML as in ML force fields). ML can also be directly used on experimental/DFT results to calculate $\Delta G_m$. MD image is adapted with permission.[88] Copyright 2021, Elsevier. CG image is adapted with permission.[146] CC-BY-ND-NC license. DFT image is adapted with permission.[240] Copyright 2019, American Chemical Society. ML image is adapted with permission.[77] Copyright 2019, American Chemical Society.</figcaption></figure>

<figure><img src="image_20.png" /></figure>
The practical implementation of such a conceptual optimization/design scheme will depend on the availability of effective procedures for mapping the chemical molecular structures (e.g., based on the SMILES representation[238,239]) and by efficient optimization algorithms (nonlinear optimizers, evolutionary algorithms[199]).

The methodological progress of recent years, the constant
increase in high-performance computing resources (GPU,
exascale computing), and the exploitation of artificial intelligence
suggest that computational methods will have a strong impact in
the near future for the design of stable OPV blends fulfilling
commercial requirements.

Acknowledgements

The authors acknowledge CNR for funding under Bilateral agreement CNR-RFBR Project CUP B55F21000620005 “Fullerene containers for BeNCT”; MUR for funding under Project PON04a2_00490 “M2M Netergit” CUP B22112000280001; and CINECA for awarding access to high-performance computing resources under the ISCRA initiative.

Open Access Funding provided by Consiglio Nazionale delle Ricerche within the CRUI-CARE Agreement.

Conflict of Interest

The authors declare no conflict of interest.

Keywords

density functional theory, Flory-Huggins, machine learning, mixing
enthalpy, molecular dynamics, organic photovoltaics

Received: April 26, 2022

Revised: May 28, 2022

Published online:

<figure><img src="image_21.png" /></figure>

<footnote>[1] S. Park, T. Kim, S. Yoon, C. W. Koh, H. Y. Woo, H. J. Son, *Adv. Mater.* **2020**, *32*, 2002217.</footnote>

<footnote>[2] C. J. Brabec, A. Distler, X. Du, H. Egelhaaf, J. Hauch, T. Heumueller, N. Li, *Adv. Mater.* **2020**, *10*, 2001864.</footnote>

<footnote>[3] M. Moser, A. Wadsworth, N. Gasparini, I. McCulloch, *Adv. Mater.* **2021**, *11*, 202100056.</footnote>

<footnote>[4] J. Peet, A. J. Heeger, G. C. Bazan, *Acc. Chem. Res.* **2009**, *42*, 1700.</footnote>

<footnote>[5] N. S. Sariciftci, L. Smilowitz, A. J. Heeger, F. Wudl, *Science* **1992**, *258*, 1474.</footnote>

<footnote>[6] N. S. Sariciftci, D. Braun, C. Zhang, V. I. Srdanov, A. J. Heeger, G. Stucky, F. Wudl, *Appl. Phys. Lett.* **1993**, *62*, 585.</footnote>

<footnote>[7] L. Smilowitz, N. S. Sariciftci, R. Wu, C. Gettinger, A. J. Heeger, F. Wudl, *Phys. Rev. B* **1993**, *47*, 13835.</footnote>

<footnote>[8] J. Huang, C.-Z. Li, C.-C. Chueh, S.-Q. Liu, J.-S. Yu, A. K.-Y. Jen, *Adv. Mater.* **2015**, *5*, 1500406.</footnote>

<footnote>[9] National Renewable Energy Laboratory (NREL), *Best Research-Cell Efficiency Chart*, https://www.nrel.gov/pv/cell-efficiency.html (accessed: March 2022).</footnote>

<footnote>[10] A. Polman, M. Knight, E. C. Garnett, B. Ehrler, W. C. Sinke, *Science* **2016**, *352*, aad4424.</footnote>

<footnote>[11] J. Yuan, Y. Zhang, L. Zhou, G. Zhang, H. L. Yip, T. K. Lau, X. Lu, C. Zhu, H. Peng, P. A. Johnson, M. Leclerc, Y. Cao, J. Ulanski, Y. Li, Y. Zou, *Joule* **2019**, *3*, 1140.</footnote>

<footnote>[12] M. Zhang, X. Guo, W. Ma, H. Ade, J. Hou, *Adv. Mater.* **2015**, *27*, 4655.</footnote>

<footnote>[13] Y. Cui, Y. Xu, H. Yao, P. Bi, L. Hong, J. Zhang, Y. Zu, T. Zhang, J. Qin, J. Ren, Z. Chen, C. He, X. Hao, Z. Wei, J. Hou, *Adv. Mater.* **2021**, *33*, 202102420.</footnote>

<footnote>[14] L. Perdigón-Toro, H. Zhang, A. Markina, J. Yuan, S. M. Hosseini, C. M. Wolff, G. Zuo, M. Stolterfoht, Y. Zou, F. Gao, D. Andrienko, S. Shoae, D. Neher, *Adv. Mater.* **2020**, *32*, 1906763.</footnote>

<footnote>[15] N. Tokmoldin, S. M. Hosseini, M. Raoufi, L. Q. Phuong, O. J. Sandberg, H. Guan, Y. Zou, D. Neher, S. Shoae, *J. Mater. Chem. A* **2020**, *8*, 7854.</footnote>

<footnote>[16] A. Karki, J. Vollbrecht, A. L. Dixon, N. Schopp, M. Schrock, G. N. M. Reddy, T. Nguyen, *Adv. Mater.* **2019**, *31*, 1903868.</footnote>

<footnote>[17] G. Zhang, X.-K. Chen, J. Xiao, P. C. Y. Chow, M. Ren, G. Kupgan, X. Jiao, C. C. S. Chan, X. Du, R. Xia, Z. Chen, J. Yuan, Y. Zhang, S. Zhang, Y. Liu, Y. Zou, H. Yan, K. S. Wong, V. Coropceanu, N. Li, C. J. Brabec, J.-L. Bredas, H.-L. Yip, Y. Cao, *Nat. Commun.* **2020**, *11*, 3943.</footnote>

<footnote>[18] Y. Cui, H. Yao, J. Zhang, K. Xian, T. Zhang, L. Hong, Y. Wang, Y. Xu, K. Ma, C. An, C. He, Z. Wei, F. Gao, J. Hou, *Adv. Mater.* **2020**, *32*, 1908205.</footnote>

<footnote>[19] A. Karki, J. Vollbrecht, A. J. Gillett, S. S. Xiao, Y. Yang, Z. Peng, N. Schopp, A. L. Dixon, S. Yoon, M. Schrock, H. Ade, G. N. M. Reddy, R. H. Friend, T.-Q. Nguyen, *Energy Environ. Sci.* **2020**, *13*, 3679.</footnote>

<footnote>[20] B. R. Luginbuhl, P. Raval, T. Pawlak, Z. Du, T. Wang, G. Kupgan, N. Schopp, S. Chae, S. Yoon, A. Yi, H. Jung Kim, V. Coropceanu, J. Brédas, T. Nguyen, G. N. M. Reddy, *Adv. Mater.* **2022**, *34*, 2105943.</footnote>

<footnote>[21] W. Zhu, A. P. Spencer, S. Mukherjee, J. M. Alzola, V. K. Sangwan, S. H. Amsterdam, S. M. Swick, L. O. Jones, M. C. Heiber, A. A. Herzing, G. Li, C. L. Stern, D. M. DeLongchamp, K. L. Kohlstedt, M. C. Hersam, G. C. Schatz, M. R. Wäsielewski, L. X. Chen, A. Facchetti, T. J. Marks, *J. Am. Chem. Soc.* **2020**, *142*, 14532.</footnote>

<footnote>[22] V. Coropceanu, X.-K. Chen, T. Wang, Z. Zheng, J.-L. Brédas, *Nat. Rev. Mater.* **2019**, *4*, 689.</footnote>

<footnote>[23] C. Li, J. Zhou, J. Song, J. Xu, H. Zhang, X. Zhang, J. Guo, L. Zhu, D. Wei, G. Han, J. Min, Y. Zhang, Z. Xie, Y. Yi, H. Yan, F. Gao, F. Liu, Y. Sun, *Nat. Energy* **2021**, *6*, 605.</footnote>

<footnote>[24] M. Green, E. Dunlop, J. Hohl-Ebinger, M. Yoshita, N. Kopidakis, X. Hao, *Prog. Photovoltaics: Res. Appl.* **2021**, *29*, 3.</footnote>

<footnote>[25] G. Hou, H. Sun, Z. Jiang, Z. Pan, Y. Wang, X. Zhang, Y. Zhao, Q. Yao, *Appl. Energy* **2016**, *164*, 882.</footnote>

<footnote>[26] G. Bernardo, T. Lopes, D. G. Lidzey, A. Mendes, *Adv. Mater.* **2021**, *11*, 202100342.</footnote>

<footnote>[27] Z. Liu, S. E. Sofia, H. S. Laine, M. Woodhouse, S. Wieghold, I. M. Peters, T. Buonassisi, *Energy Environ. Sci.* **2020**, *13*, 12.</footnote>

<footnote>[28] L. Meng, J. You, Y. Yang, *Nat. Commun.* **2018**, *9*, 5265.</footnote>

<footnote>[29] Q. Burlingame, M. Ball, Y.-L. Loo, *Nat. Energy* **2020**, *5*, 947.</footnote>

<footnote>[30] N. S. Sariciftci, A. J. Hebger, *Int. J. Mod. Phys. B* **1994**, *08*, 237.</footnote>

<footnote>[31] S. Morita, A. A. Zakhidov, K. Yoshino, *Solid State Commun.* **1992**, *82*, 249.</footnote>

<footnote>[32] X. Xu, K. Fukuda, A. Karki, S. Park, H. Kimura, H. Jinno, N. Watanabe, S. Yamamoto, S. Shimomura, D. Kitazawa, T. Yokota, S. Umezu, T.-Q. Nguyen, T. Someya, *Proc. Natl. Acad. Sci.* **2018**, *115*, 4589.</footnote>

<footnote>[33] A. Karki, A. J. Gillett, R. H. Friend, T. Nguyen, *Adv. Mater.* **2021**, *11*, 2003441.</footnote>

<footnote>[34] C. Yan, S. Barlow, Z. Wang, H. Yan, A. K.-Y. Jen, S. R. Marder, X. Zhan, *Nat. Rev. Mater.* **2018**, *3*, 18003.</footnote>

<footnote>[35] S. Holliday, R. S. Ashraf, A. Wadsworth, D. Baran, S. A. Yousaf, C. B. Nielsen, C.-H. Tan, S. D. Dimitrov, Z. Shang, N. Gasparini, M. Alamoudi, F. Laquai, C. J. Brabec, A. Salleo, J. R. Durrant, I. McCulloch, *Nat. Commun.* **2016**, *7*, 11585.</footnote>

<footnote>[36] S. Li, L. Zhan, F. Liu, J. Ren, M. Shi, C.-Z. Li, T. P. Russell, H. Chen, *Adv. Mater.* **2018**, *30*, 1705208.</footnote>
[37] Y. Li, X. Huang, K. Ding, H. K. M. Sheriff, L. Ye, H. Liu, C.-Z. Li, H. Ade, S. R. Forrest, Nat. Commun. **2021**, 12, 5419.

[38] K. Ding, Y. Li, S. R. Forrest, ACS Appl. Mater. Interfaces **2022**, 14, 5692.

[39] L. Duan, A. Uddin, Adv. Sci. **2020**, 7, 1903259.

[40] F. Zhao, H. Zhang, R. Zhang, J. Yuan, D. He, Y. Zou, F. Gao, *Adv. Mater.* **2020**, *10*, 202002746.

[41] A. Wadsworth, Z. Hamid, J. Kosco, N. Gasparini, I. McCulloch, *Adv. Mater.* **2020**, *32*, 2001763.

[42] Y. Wang, J. Han, L. Cai, N. Li, Z. Li, F. Zhu, *J. Mater. Chem. A* **2020**, *8*, 21255.

[43] Y. Wang, M. J. Jafari, N. Wang, D. Qian, F. Zhang, T. Ederth, E. Moons, J. Wang, O. Inganäs, W. Huang, F. Gao, *J. Mater. Chem. A* **2018**, *6*, 11884.

[44] H. K. H. Lee, A. M. Telford, J. A. Röhr, M. F. Wyatt, B. Rice, J. Wu, A. de Castro Maciel, S. M. Tuladhar, E. Speller, J. McGettrick, J. R. Searle, S. Pont, T. Watson, T. Kirchartz, J. R. Durrant, W. C. Tsoi, J. Nelson, Z. Li, *Energy Environ. Sci.* **2018**, *11*, 417.

[45] Y.-C. Huang, W.-S. Liu, C.-S. Tsao, L. Wang, ACS Appl. Mater. Interfaces **2019**, *11*, 40310.

[46] C. J. Schaffer, C. M. Palumbiny, M. A. Niedermeier, C. Jendrzejewski, G. Santoro, S. V. Roth, P. Müller-Buschbaum, *Adv. Mater.* **2013**, *25*, 6760.

[47] P. Boldrighini, A. Fauveau, S. Thérias, J. L. Gardette, M. Hidalgo, S. Cros, *Rev. Sci. Instrum.* **2019**, *90*, 014710.

[48] X. Wang, C. Xinxin Zhao, G. Xu, Z.-K. Chen, F. Zhu, *Sol. Energy Mater. Sol. Cells* **2012**, *104*, 1.

[49] V. M. Drakonakis, A. Savva, M. Kokonou, S. A. Choulis, *Sol. Energy Mater. Sol. Cells* **2014**, *130*, 544.

[50] P. Cheng, X. Zhan, *Chem. Soc. Rev.* **2016**, *45*, 2544.

[51] C. McDowell, M. Abdelsamie, M. F. Toney, G. C. Bazan, *Adv. Mater.* **2018**, *30*, 1707114.

[52] N. Grossiord, J. M. Kroon, R. Andriessen, P. W. M. Blom, *Org. Electron.* **2012**, *13*, 432.

[53] M. Jørgensen, K. Norrman, F. C. Krebs, *Sol. Energy Mater. Sol. Cells* **2008**, *92*, 686.

[54] S. Rafique, S. M. Abdullah, K. Sulaiman, M. Iwamoto, *Renewable Sustainable Energy Rev.* **2018**, *84*, 43.

[55] H. Cao, W. He, Y. Mao, X. Lin, K. Ishikawa, J. H. Dickerson, W. P. Hess, *J. Power Sources* **2014**, *264*, 168.

[56] I. Fraga Domínguez, A. Distler, L. Lüer, *Adv. Mater.* **2017**, *7*, 1601320.

[57] E. M. Speller, A. J. Clarke, J. Luke, H. K. H. Lee, J. R. Durrant, N. Li, T. Wang, H. C. Wong, J.-S. Kim, W. C. Tsoi, Z. Li, *J. Mater. Chem. A* **2019**, *7*, 23361.

[58] E. K. Lee, M. Y. Lee, C. H. Park, H. R. Lee, J. H. Oh, *Adv. Mater.* **2017**, *29*, 1703638.

[59] M. Jørgensen, K. Norrman, S. A. Gevorgyan, T. Tromholt, B. Andreasen, F. C. Krebs, *Adv. Mater.* **2012**, *24*, 580.

[60] M. Hösel, R. R. Søndergaard, M. Jørgensen, F. C. Krebs, *Adv. Mater.* **2014**, *4*, 1301625.

[61] S. Roland, M. Schubert, B. A. Collins, J. Kurpiers, Z. Chen, A. Facchetti, H. Ade, D. Neher, *J. Phys. Chem. Lett.* **2014**, *5*, 2815.

[62] S. Mukherjee, C. M. Proctor, G. C. Bazan, T.-Q. Nguyen, H. Ade, *Adv. Mater.* **2015**, *5*, 1500877.

[63] L. Ye, X. Jiao, S. Zhang, H. Yao, Y. Qin, H. Ade, J. Hou, *Adv. Mater.* **2017**, *7*, 1601138.

[64] W. Köntges, P. Perkhun, J. Kammerer, R. Alkarsifi, U. Würfel, O. Margeat, C. Videlot-Ackermann, J. J. Simon, R. R. Schröder, R. R. Schröder, M. Pfannmöller, *Energy Environ. Sci.* **2020**, *13*, 1259.

[65] Q. Liang, J. Han, C. Song, X. Yu, D. M. Smilgies, K. Zhao, J. Liu, Y. Han, *J. Mater. Chem. A* **2018**, *6*, 15610.

[66] X. Wang, Y. Yang, Z. He, H. Wu, Y. Cao, *J. Mater. Chem. C* **2019**, *7*, 14861.

[67] H. Jiang, X. Li, H. Wang, G. Huang, W. Chen, R. Zhang, R. Yang, ACS Appl. Mater. Interfaces **2020**, *12*, 26286.

[68] G. Gryn'ova, K.-H. Lin, C. Corminboeuf, *J. Am. Chem. Soc.* **2018**, *140*, 16370.

[69] M. Seifrid, G. N. M. Reddy, B. F. Chmelka, G. C. Bazan, *Nat. Rev. Mater.* **2020**, *5*, 910.

[70] P. Zhan, W. Zhang, I. E. Jacobs, D. M. Nisson, R. Xie, A. R. Weissen, R. H. Colby, A. J. Moulé, S. T. Milner, J. K. Maranas, E. D. Gomez, *J. Polym. Sci., Part B: Polym. Phys.* **2018**, *56*, 1193.

[71] N. C. Miller, E. Cho, M. J. N. Junk, R. Gysel, C. Risko, D. Kim, S. Sweetnam, C. E. Miller, L. J. Richter, R. J. Kline, M. Heeney, I. McCulloch, A. Amassian, D. Acevedo-Feliz, C. Knox, M. R. Hansen, D. Dudenko, B. F. Chmelka, M. F. Toney, J. Brédas, M. D. McGehee, *Adv. Mater.* **2012**, *24*, 6071.

[72] A. Kumar, K. Y. J. Zhang, *Methods* **2015**, *71*, 26.

[73] A. Lavecchia, C. Giovanni, *Curr. Med. Chem.* **2013**, *20*, 2839.

[74] S. Kar, K. Roy, *Expert Opin. Drug Discovery* **2013**, *8*, 245.

[75] E. Lionta, G. Spyrou, D. Vassilatis, Z. Cournia, *Curr. Top. Med. Chem.* **2014**, *14*, 1923.

[76] K. Do, M. K. Ravva, T. Wang, J.-L. Brédas, *Chem. Mater.* **2017**, *29*, 346.

[77] S. Nagasawa, E. Al-Naamani, A. Saeki, *J. Phys. Chem. Lett.* **2018**, *9*, 2639.

[78] A. Mahmood, A. Irfan, J.-L. Wang, *J. Mater. Chem. A* **2022**, *10*, 4170.

[79] L. Huynh, C. Neale, R. Pomès, C. Allen, *Nanomed. Nanotechnol. Biol. Med.* **2012**, *8*, 20.

[80] T.-X. Xiang, B. D. Anderson, *J. Pharm. Sci.* **2013**, *102*, 876.

[81] D. C. Santos, S. D. Filipakis, M. P. Rolemburg, E. R. A. Lima, M. L. L. Paredes, *Fuel* **2017**, *199*, 606.

[82] S. Alimohammadi, J. Sayyad Amin, E. Nikooee, *Neural Comput. Appl.* **2017**, *28*, 679.

[83] C. Caddeo, A. Mattoni, *Macromolecules* **2013**, *46*, 8003.

[84] C. Caddeo, D. Fazzi, M. Caironi, A. Mattoni, *J. Phys. Chem. B* **2014**, *118*, 12556.

[85] S. Bellani, M. Porro, C. Caddeo, M. I. Saba, P. B. Miranda, A. Mattoni, G. Lanzani, M. R. Antognazza, *J. Mater. Chem. B* **2015**, *3*, 6429.

[86] S. Ben Dkhil, M. Pfannmöller, M. I. Saba, M. Gaceur, H. Heidari, C. Videlot-Ackermann, O. Margeat, A. Guerrero, J. Bisquert, G. Garcia-Belmonte, A. Mattoni, S. Bals, J. Ackermann, *Adv. Mater.* **2017**, *7*, 1601486.

[87] S. Ben Dkhil, P. Perkhun, C. Luo, D. Müller, R. Alkarsifi, E. Barulina, Y. A. Avalos Quiroz, O. Margeat, S. T. Dubas, T. Koganezawa, D. Kuzuhara, N. Yoshimoto, C. Caddeo, A. Mattoni, B. Zimmermann, U. Würfel, M. Pfannmöller, S. Bals, J. Ackermann, C. Videlot-Ackermann, ACS Appl. Mater. Interfaces **2020**, *12*, 28404.

[88] C. Caddeo, A. Filippetti, A. Bosin, C. Videlot-Ackermann, J. Ackermann, A. Mattoni, *Nano Energy* **2021**, *82*, 105708.

[89] J. D. Perea Ospina, S. Langner, T. Ameri, C. J. Brabec, in *Encyclopedia of Physical Organic Chemistry* (Eds.: Z. Wang, U. Wille, E. Juaristi), John Wiley & Sons, Inc., Hoboken NJ **2017**, p. 697.

[90] Y. S. Lipatov, A. E. Nesterov, *Thermodynamics of Polymer Blends*, CRC Press, Boca Raton, FL **2020**.

[91] B. A. Wolf, in *Polymer Thermodynamics* (Eds.: B.A. Wolf, S. Enders), Springer, Berlin, Heidelberg, **2011**, pp. 1–66.

[92] R. H. Colby, M. Rubinstein, *Polymer Physics*, Oxford University Press, Oxford **2003**.

[93] R. A. L. Jones, *Soft Condensed Matter*, Oxford University Press, New York, **2002**.

[94] P. J. Flory, *J. Chem. Phys.* **1942**, *10*, 51.

[95] M. L. Huggins, *J. Chem. Phys.* **1941**, *9*, 440.

[96] H. Hildebrand, *Proc. Natl. Acad. Sci.* **1927**, *13*, 267.

[97] J. Zhu, X. Lu, R. Balieu, N. Kringos, *Mater. Des.* **2016**, *107*, 322.

<figure><img src="image_22.png" /></figure>
[98] M. Ghasemi, H. Hu, Z. Peng, J. J. Rech, I. Angunawela, J. H. Carpenter, S. J. Stuard, A. Wadsworth, I. McCulloch, W. You, H. Ade, *Joule* **2019**, *3*, 1328.

[99] J. H. Hildebrand, R. L. Scott, *The Solubility of Nonelectrolytes*, 3rd Ed., Reinhold, New York, **1950**.

[100] H. C. van Ness, *Classical Thermodynamics Of Nonelectrolyte Solutions : With Applications To Phase Equilibria*, McGraw-Hill, New York, NY **1982**.

[101] E. R. Turpin, V. Taresco, W. A. Al-Hachami, J. Booth, K. Treacher, S. Tomasi, C. Alexander, J. Burley, C. A. Laughton, M. C. Garnett, *Mol. Pharmaceutics* **2018**, *15*, 4654.

[102] F. S. Bates, *Science* **1991**, *251*, 898.

[103] P. P. Khlyabich, A. E. Rudenko, B. C. Thompson, Y.-L. Loo, *Adv. Funct. Mater.* **2015**, *25*, 5557.

[104] M. Ghasemi, L. Ye, Q. Zhang, L. Yan, J. H. Kim, O. Awartani, W. You, A. Gadisa, H. Ade, *Adv. Mater.* **2017**, *29*, 1604603.

[105] J. A. Emerson, D. T. W. Toolan, J. R. Howse, E. M. Furst, T. H. Epps, *Macromolecules* **2013**, *46*, 6533.

[106] T. K. Kwei, T. Nishi, R. F. Roberts, *Macromolecules* **1974**, *7*, 667.

[107] T. Nishi, T. T. Wang, *Macromolecules* **1975**, *8*, 909.

[108] D. Lin, Y. Huang, *Int. J. Pharm.* **2010**, *399*, 109.

[109] D. R. Kozub, K. Vakhshouri, L. M. Orme, C. Wang, A. Hexemer, E. D. Gomez, *Macromolecules* **2011**, *44*, 5722.

[110] N. Yi, Q. Ai, W. Zhou, L. Huang, L. Zhang, Z. Xing, X. Li, J. Zeng, Y. Chen, *Chem. Mater.* **2019**, *31*, 10211.

[111] L. Ye, H. Hu, M. Ghasemi, T. Wang, B. A. Collins, J. H. Kim, K. Jiang, J. H. Carpenter, H. Li, Z. Li, T. McAfee, J. Zhao, X. Chen, J. L. Y. Lai, T. Ma, J. L. Bredas, H. Yan, H. Ade, *Nat. Mater.* **2018**, *17*, 253.

[112] L. Ye, W. Zhao, S. Li, S. Mukherjee, J. H. Carpenter, O. Awartani, X. Jiao, J. Hou, H. Ade, *Adv. Mater.* **2017**, *7*, 1602000.

[113] L. Ye, B. A. Collins, X. Jiao, J. Zhao, H. Yan, H. Ade, *Adv. Mater.* **2018**, *8*, 1703058.

[114] H. G. Harris, J. M. Prausnitz, *Ind. Eng. Chem. Fundam.* **1969**, *8*, 180.

[115] J. H. Hildebrand, R. L. Scott, *Regular Solutions*, Prentice-Hall, Englewood Cliffs, N.J., **1962**.

[116] C. M. Hansen, *Hansen Solubility Parameters*, CRC Press, Boca Raton, **2007**.

[117] P. A. Small, *J. Appl. Chem.* **1953**, *3*, 71.

[118] K. Hoy, *J. Paint Technol.* **1970**, *42*, 76.

[119] P. Hoftyzer, D. van Krevelen, *Properties of Polymers*, 1st Ed., Elsevier, Amsterdam, **1976**.

[120] R. F. Fedors, *Polym. Eng. Sci.* **1974**, *14*, 147.

[121] K. DeBoyace, P. L. D. Wildfong, *J. Pharmaceutical Sci.* **2018**, *107*, 57.

[122] S. Patel, A. Lavasanifar, P. Choi, *Biomacromolecules* **2008**, *9*, 3014.

[123] S. K. Patel, A. Lavasanifar, P. Choi, *Biomacromolecules* **2009**, *10*, 2584.

[124] A. O. Kasimova, G. M. Pavan, A. Danani, K. Mondon, A. Cristiani, L. Scapozza, R. Gurny, M. Möller, *J. Phys. Chem. B* **2012**, *116*, 4338.

[125] T. X. Xiang, B. D. Anderson, *J. Pharmaceutical Sci.* **2017**, *106*, 803.

[126] T. P. Straatsma, H. J. C. Berendsen, *J. Chem. Phys.* **1988**, *89*, 5876.

[127] T. P. Straatsma, H. J. C. Berendsen, J. P. M. Postma, *J. Chem. Phys.* **1986**, *85*, 6720.

[128] G. M. Torrie, J. P. Valleau, *J. Comput. Phys.* **1977**, *23*, 187.

[129] A. Barducci, G. Bussi, M. Parrinello, *Phys. Rev. Lett.* **2008**, *100*, 020603.

[130] R. W. Zwanzig, *J. Chem. Phys.* **1954**, *22*, 1420.

[131] F. Baftizadeh, P. Cossio, F. Pietrucci, A. Laio, *Curr. Phys. Chem.* **2012**, *2*, 79.

[132] S. Angioletti-Uberti, M. Ceriotti, P. D. Lee, M. W. Finnis, *Phys. Rev. B* **2010**, *81*, 125416.

[133] Q. Wei, W. Zhao, Y. Yang, B. Cui, Z. Xu, X. Yang, *ChemPhysChem* **2018**, *19*, 690.

[134] T. Mandal, P. H. Koenig, R. G. Larson, *Phys. Rev. Lett.* **2018**, *121*, 038001.

[135] A. de Izarra, C. Choi, Y. H. Jang, Y. Lansac, *J. Phys. Chem. B* **2021**, *125*, 1916.

[136] T. Ludwig, A. R. Singh, J. K. Nørskov, *J. Phys. Chem. C* **2020**, *124*, 26368.

[137] R. Réocreux, É. Girel, P. Clabaut, A. Tuel, M. Besson, A. Chaumonnot, A. Cabiac, P. Sautet, C. Michel, *Nat. Commun.* **2019**, *10*, 3139.

[138] H. J. Heelweg, R. A. de Souza, *Phys. Rev. Mater.* **2021**, *5*, 013804.

[139] L. Zhong, L. Gao, H. Bin, Q. Hu, Z.-G. Zhang, F. Liu, T. P. Russell, Z. Zhang, Y. Li, *Adv. Mater.* **2017**, *7*, 1602215.

[140] Y. Ma, D. Cai, S. Wan, P. Yin, P. Wang, W. Lin, Q. Zheng, *Natl. Sci. Rev.* **2020**, *7*, 1886.

[141] R. Yu, H. Yao, Z. Chen, J. Xin, L. Hong, Y. Xu, Y. Zu, W. Ma, J. Hou, *Adv. Mater.* **2019**, *31*, 1900477.

[142] C. B. Nielsen, S. Holliday, H.-Y. Chen, S. J. Cryer, I. McCulloch, *Acc. Chem. Res.* **2015**, *48*, 2803.

[143] J. Wang, R. M. Wolf, J. W. Caldwell, P. A. Kollman, D. A. Case, *J. Comput. Chem.* **2004**, *25*, 1157.

[144] A. S. Gertsen, M. K. Sørensen, J. W. Andreasen, *Phys. Rev. Mater.* **2020**, *4*, 075405.

[145] C.-K. Lee, C.-W. Pao, *J. Phys. Chem. C* **2014**, *118*, 11224.

[146] R. Alessandri, J. J. Uusitalo, A. H. de Vries, R. W. A. Havenith, S. J. Marrink, *J. Am. Chem. Soc.* **2017**, *139*, 3697.

[147] A. Guerrero, M. Pfannmöller, A. Kovalenko, T. S. Ripolles, H. Heidari, S. Bals, L.-D. Kaufmann, J. Bisquert, G. Garcia-Belmonte, *Org. Electron.* **2015**, *16*, 227.

[148] S. Ben Dkhil, M. Pfannmöller, R. R. Schröder, R. Alkarsifi, M. Gaceur, W. Köntges, H. Heidari, S. Bals, O. Margeat, J. Ackermann, C. Videlot-Ackermann, *ACS Appl. Mater. Interfaces* **2018**, *10*, 3874.

[149] A. Pinke, P. Jedlovszky, *J. Phys. Chem. B* **2012**, *116*, 5977.

[150] A. Idrissi, P. Jedlovszky, *J. Mol. Liq.* **2021**, *338*, 116777.

[151] A. Idrissi, B. Marekha, M. Barj, P. Jedlovszky, *J. Phys. Chem. B* **2014**, *118*, 8724.

[152] S. Gupta, T. M. Lim, S. H. Mushrif, *Electrochim. Acta* **2018**, *270*, 471.

[153] X. Liu, V. L. Vinograd, X. Lu, E. V. Leonenko, N. N. Eremin, R. Wang, B. Winkler, *Am. Mineral.* **2016**, *101*, 1197.

[154] D. Y. Jung, V. L. Vinograd, O. B. Fabrichnaya, A. R. Oganov, M. W. Schmidt, B. Winkler, *Earth Planet. Sci. Lett.* **2010**, *295*, 477.

[155] S. J. Marrink, D. P. Tieleman, *Chem. Soc. Rev.* **2013**, *42*, 6801.

[156] L. Monticelli, *J. Chem. Theory Comput.* **2012**, *8*, 1370.

[157] M. Modarresi, J. F. Franco-Gonzalez, I. Zozoulenko, *Phys. Chem. Chem. Phys.* **2018**, *20*, 17188.

[158] N. Rolland, M. Modarresi, J. F. Franco-Gonzalez, I. Zozoulenko, *Comput. Mater. Sci.* **2020**, *179*, https://doi.org/10.1016/j.commatsci.2020.109678.

[159] R. Alessandri, S. Sami, J. Barnoud, A. H. Vries, S. J. Marrink, R. W. A. Havenith, *Adv. Funct. Mater.* **2020**, *30*, 2004799.

[160] S. J. Marrink, A. H. de Vries, A. E. Mark, *J. Phys. Chem. B* **2004**, *108*, 750.

[161] S. J. Marrink, H. J. Risselada, S. Yefimov, D. P. Tieleman, A. H. de Vries, *J. Phys. Chem. B* **2007**, *111*, 7812.

[162] L. Monticelli, S. K. Kandasamy, X. Periole, R. G. Larson, D. P. Tieleman, S.-J. Marrink, *J. Chem. Theory Comput.* **2008**, *4*, 819.

[163] C. A. López, A. J. Rzepiela, A. H. de Vries, L. Dijkhuizen, P. H. Hünenberger, S. J. Marrink, *J. Chem. Theory Comput.* **2009**, *5*, 3195.

[164] J. J. Uusitalo, H. I. Ingólfsson, P. Akhshi, D. P. Tieleman, S. J. Marrink, *J. Chem. Theory Comput.* **2015**, *11*, 3932.

<figure><img src="image_23.png" /></figure>
[165] J. J. Uusitalo, H. I. Ingólfsson, S. J. Marrink, I. Faustino, *Biophys. J.* **2017**, *113*, 246.

[166] L. I. Vazquez-Salazar, M. Selle, A. H. de Vries, S. J. Marrink, P. C. T. Souza, *Green Chem.* **2020**, *22*, 7376.

[167] H. Lee, R. G. Larson, *J. Phys. Chem. B* **2008**, *112*, 7778.

[168] G. Rossi, L. Monticelli, S. R. Puisto, I. Vattulainen, T. Ala-Nissila, *Soft Matter* **2011**, *7*, 698.

[169] P. Español, P. B. Warren, *J. Chem. Phys.* **2017**, *146*, 150901.

[170] C. Du, Y. Ji, J. Xue, T. Hou, J. Tang, S.-T. Lee, Y. Li, *Sci. Rep.* **2015**, *5*, 16854.

[171] X. Xu, Y. Ji, C. Du, T. Hou, Y. Li, *RSC Adv.* **2015**, *5*, 70939.

[172] K. Bergermann, *GAMM Arch. Stud.* **2019**, *1*, 18.

[173] R. Rabani, H. Sadafi, H. Machrafi, M. Abbasi, B. Haut, P. Dauby, *Colloids Surf.*, **A** *2021*, *612*, 126001.

[174] R. C. Ball, R. L. H. Essery, *J. Phys. Condens. Matter* **1990**, *2*, 10303.

[175] I. Steinbach, *Annu. Rev. Mater. Res.* **2013**, *43*, 89.

[176] A. Choudhury, in *Handbook of Solid State Diffusion*, Vol. 1, Elsevier, **2017**, pp. 363–389.

[177] J. W. Cahn, *J. Chem. Phys.* **1965**, *42*, 93.

[178] O. Wodo, B. Ganapathysubramanian, *Comput. Mater. Sci.* **2012**, *55*, 113.

[179] B. Ray, P. R. Nair, M. A. Alam, *Sol. Energy Mater. Sol. Cells* **2011**, *95*, 3287.

[180] F. Kaka, M. Keshav, P. C. Ramamurthy, *Sol. Energy* **2022**, *231*, 447.

[181] P. Hohenberg, W. Kohn, *Phys. Rev.* **1964**, *136*, B864.

[182] W. Kohn, L. J. Sham, *Phys. Rev.* **1965**, *140*, A1133.

[183] G. Járvás, C. Quellet, A. Dallos, *Fluid Phase Equilib.* **2011**, *309*, 8.

[184] P. Friederich, F. Symalla, V. Meded, T. Neumann, W. Wenzel, *J. Chem. Theory Comput.* **2014**, *10*, 3720.

[185] S. Roosta, F. Ghalami, M. Elstner, W. Xie, *J. Chem. Theory Comput.* **2022**, *18*, 1264.

[186] E. Brunk, U. Rothlisberger, *Chem. Rev.* **2015**, *115*, 6217.

[187] A. S. Christensen, T. Kubař, Q. Cui, M. Elstner, *Chem. Rev.* **2016**, *116*, 5301.

[188] J. J. P. Stewart, *J. Mol. Model.* **2013**, *19*, 1.

[189] F. Häse, S. Valleau, E. Pyzer-Knapp, A. Aspuru-Guzik, *Chem. Sci.* **2016**, *7*, 5139.

[190] L. Böselt, M. Thürlemann, S. Riniker, *J. Chem. Theory Comput.* **2021**, *17*, 2641.

[191] G. Han, Y. Guo, X. Ma, Y. Yi, *Sol. RRL* **2018**, *2*, 1800190.

[192] Y. Li, V. Agrawal, J. Oswald, *J. Polym. Sci., Part B: Polym. Phys.* **2019**, *57*, 331.

[193] J. D. Perea, S. Langner, M. Salvador, J. Kontos, G. Jarvas, F. Winkler, F. Machui, A. Görling, A. Dallos, T. Ameri, C. J. Brabec, *J. Phys. Chem. B* **2016**, *120*, 4431.

[194] J. D. Perea, S. Langner, M. Salvador, B. Sanchez-Lengeling, N. Li, C. Zhang, G. Jarvas, J. Kontos, A. Dallos, A. Aspuru-Guzik, C. J. Brabec, *J. Phys. Chem. C* **2017**, *121*, 18153.

[195] M.-Y. Sui, Z.-R. Yang, Y. Geng, G.-Y. Sun, L. Hu, Z.-M. Su, *Sol. RRL* **2019**, *3*, 1900258.

[196] Y.-C. Lin, Y.-J. Lu, C.-S. Tsao, A. Saeki, J.-X. Li, C.-H. Chen, H.-C. Wang, H.-C. Chen, D. Meng, K.-H. Wu, Y. Yang, K.-H. Wei, *J. Mater. Chem. A* **2019**, *7*, 3072.

[197] Y. Wu, J. Guo, R. Sun, J. Min, *npj Comput. Mater.* **2020**, *6*, 120.

[198] K. Kranthiraja, A. Saeki, *Adv. Funct. Mater.* **2021**, *31*, 2011168.

[199] Z. W. Zhao, M. del Cueto, Y. Geng, A. Troisi, *Chem. Mater.* **2020**, *32*, 7777.

[200] Z. W. Zhao, Ö. H. Omar, D. Padula, Y. Geng, A. Troisi, *J. Phys. Chem. Lett.* **2021**, *12*, 5009.

[201] W. Sun, Y. Zheng, Q. Zhang, K. Yang, H. Chen, Y. Cho, J. Fu, O. Odunmbaku, A. A. Shah, Z. Xiao, S. Lu, S. Chen, M. Li,

B. Qin, C. Yang, T. Frauenheim, K. Sun, *J. Phys. Chem. Lett.* **2021**, *12*, 8847.

[202] Y. Cui, P. Zhu, X. Liao, Y. Chen, *J. Mater. Chem. C* **2020**, *8*, 15920.

[203] J. Klimeš, A. Michaelides, *J. Chem. Phys.* **2012**, *137*, 120901.

[204] M. Stöhr, T. van Voorhis, A. Tkatchenko, *Chem. Soc. Rev.* **2019**, *48*, 4118.

[205] C. Li, X. Zheng, N. Q. Su, W. Yang, *Natl. Sci. Rev.* **2018**, *5*, 203.

[206] G. Pilania, C. Wang, X. Jiang, S. Rajasekaran, R. Ramprasad, *Sci. Rep.* **2013**, *3*, 2810.

[207] K. T. Butler, D. W. Davies, H. Cartwright, O. Isayev, A. Walsh, *Nature* **2018**, *559*, 547.

[208] G. R. Schleder, B. Focassio, A. Fazzio, *Appl. Phys. Rev.* **2021**, *8*, 031409.

[209] A. Mahmood, J.-L. Wang, *Energy Environ. Sci.* **2021**, *14*, 90.

[210] R. Olivares-Amaya, C. Amador-Bedolla, J. Hachmann, S. Atahan-Evrenk, R. S. Sánchez-Carrera, L. Vogt, A. Aspuru-Guzik, *Energy Environ. Sci.* **2011**, *4*, 4849.

[211] L. Pattanaik, C. W. Coley, *Chem* **2020**, *6*, 1204.

[212] N. Marzari, A. Ferretti, C. Wolverton, *Nat. Mater.* **2021**, *20*, 736.

[213] C. Duan, F. Liu, A. Nandy, H. J. Kulik, *J. Phys. Chem. Lett.* **2021**, *12*, 4628.

[214] Y. Miyake, A. Saeki, *J. Phys. Chem. Lett.* **2021**, *12*, 12391.

[215] A. Jain, S. P. Ong, G. Hautier, W. Chen, W. D. Richards, S. Dacek, S. Cholia, D. Gunter, D. Skinner, G. Ceder, K. A. Persson, *APL Mater.* **2013**, *1*, 011002.

[216] S. Curtarolo, W. Setyawan, S. Wang, J. Xue, K. Yang, R. H. Taylor, L. J. Nelson, G. L. W. Hart, S. Sanvito, M. Buongiorno-Nardelli, N. Mingo, O. Levy, *Comput. Mater. Sci.* **2012**, *58*, 227.

[217] S. Kirklin, J. E. Saal, B. Meredig, A. Thompson, J. W. Doak, M. Aykol, S. Rühl, C. Wolverton, *npj Comput. Mater.* **2015**, *1*, 15010.

[218] L. Talirz, S. Kumbhar, E. Passaro, A. V. Yakutovich, V. Granata, F. Gargiulo, M. Borelli, M. Uhrin, S. P. Huber, S. Zoupanos, C. S. Adorf, C. W. Andersen, O. Schütt, C. A. Pignedoli, D. Passerone, J. VandeVondele, T. C. Schulthess, B. Smit, G. Pizzi, N. Marzari, *Sci. Data* **2020**, *7*, 299.

[219] T. Wang, G. Kupgan, J. L. Brédas, *Trends Chem.* **2020**, *2*, 535.

[220] J. Behler, M. Parrinello, *Phys. Rev. Lett.* **2007**, *98*, 146401.

[221] J. Behler, *Angew. Chem., Int. Ed.* **2017**, *56*, 12828.

[222] A. P. Bartók, M. C. Payne, R. Kondor, G. Csányi, *Phys. Rev. Lett.* **2010**, *104*, 136403.

[223] Z. Li, J. R. Kermode, A. de Vita, *Phys. Rev. Lett.* **2015**, *114*, 096405.

[224] R. Jinnouchi, F. Karsai, G. Kresse, *Phys. Rev. B* **2019**, *100*, 014105.

[225] P. Friederich, F. Häse, J. Proppe, A. Aspuru-Guzik, *Nat. Mater.* **2021**, *20*, 750.

[226] Y. Qi, H. I. Ingólfsson, X. Cheng, J. Lee, S. J. Marrink, W. Im, *J. Chem. Theory Comput.* **2015**, *11*, 4486.

[227] G. T. Johnson, L. Autin, M. Al-Alusi, D. S. Goodsell, M. F. Sanner, A. J. Olson, *Nat. Methods* **2015**, *12*, 85.

[228] T. A. Wassenaar, H. I. Ingólfsson, R. A. Böckmann, D. P. Tieleman, S. J. Marrink, *J. Chem. Theory Comput.* **2015**, *11*, 2144.

[229] T. Bereau, K. Kremer, *J. Chem. Theory Comput.* **2015**, *11*, 2783.

[230] M. A. Webb, J.-Y. Delannoy, J. J. de Pablo, *J. Chem. Theory Comput.* **2019**, *15*, 1199.

[231] F. Grünewald, R. Alessandri, P. C. Kroon, L. Monticelli, P. C. T. Souza, S. J. Marrink, *Nat. Commun.* **2022**, *13*, 68.

[232] T. D. Potter, E. L. Barrett, M. A. Miller, *J. Chem. Theory Comput.* **2021**, *17*, 5777.

[233] A. E. P. Durumeric, G. A. Voth, *J. Chem. Phys.* **2019**, *151*, 124110.

[234] W. Wang, R. Gómez-Bombarelli, *npj Comput. Mater.* **2019**, *5*, 125.

<figure><img src="image_24.png" /></figure>
[235] J. Wang, S. Olsson, C. Wehmeyer, A. Pérez, N. E. Charron, G. de Fabritiis, F. Noé, C. Clementi, ACS Cent. Sci. **2019**, *5*, 755.

[236] R. Bhowmik, R. J. Berry, M. F. Durstock, B. J. Leever, ACS Appl. Mater. Interfaces **2017**, *9*, 19269.

[237] R. Bhowmik, R. J. Berry, V. Varshney, M. F. Durstock, B. J. Leever, J. Phys. Chem. C **2015**, *119*, 27909.

[238] N. M. O'boyle, A. Dalke, DeepSMILES: An Adaptation of SMILES for Use in Machine-Learning of Chemical Structures, ChemRxiv, Cambridge Open Engage, Cambridge **2018**.

[239] G. A. Pinheiro, J. Mucelini, M. D. Soares, R. C. Prati, J. L. F. da Silva, M. G. Quiles, *J. Phys. Chem. A* **2020**, *124*, 9854.

[240] L. Benatto, M. Koehler, *J. Phys. Chem. C* **2019**, *123*, 6395.

<figure><img src="image_25.png" /></figure>

**Claudia Caddeo** is staff researcher at the Italian Research Council, Istituto Officina dei Materiali. She received her M.S. in electronics engineering in 2009 and Ph.D. in physics at the University of Cagliari in 2013. She has been working for Centre National d'Études Spatiales (French Space Agency) in Toulouse, University of Cagliari, and Catholic University of Brescia. Her research interests focus on multiscale modeling of hybrid and polymer-based nanomaterials for applications in third-generation photovoltaics, electronics, and sensors.

<figure><img src="image_26.png" /></figure>

**Jörg Ackermann** is a research director at the CNRS in the Interdisciplinary Center of Nanoscience of Marseille (CINaM) in France. He received his physics diploma from the Saarland University, Germany, in 1998. After working as an engineer at the Bavarian Center for Applied Energy Research, he obtained his Ph.D. in physical chemistry in 2002 from the University of Aix-Marseille 2 in France. After two years as a postdoc fellow, he started his academic career at the CNRS in 2004. His research focuses on the development of inorganic, organic, and hybrid semiconducting nanomaterials, corresponding ink formulations, and printing of optoelectronic devices such as solar cells.

<figure><img src="image_27.png" /></figure>

**Alessandro Mattoni** is a research director at the Italian Research Council and director of the unit of Cagliari of the Istituto Officina dei Materiali, where he coordinates a theory group on the multiscale modeling of nanomaterials. He received a master's degree in physics at the University of Perugia and a Ph.D. in solid-state physics at the University of Padova. His research focuses on polymers and hybrid perovskites for photovoltaics and nanoelectronics. He developed the first interatomic force field for classical molecular dynamics of hybrid perovskites.

<figure><img src="image_28.png" /></figure>
<figure><img src="image_29.png" /></figure>