# Production of $X_{cs\bar{c}\bar{s}}$ in heavy ion collisions

Yuanyuan Hu  Guangdong Provincial Key Laboratory of Nuclear Science, Institute
of Quantum Matter, South China Normal University, Guangzhou 510006, China
Guangdong-Hong Kong Joint Laboratory of Quantum Matter, Southern Nuclear
Science Computing Center, South China Normal University, Guangzhou 510006,
China  Hui Zhang<EMAIL_ADDRESS>Guangdong Provincial Key
Laboratory of Nuclear Science, Institute of Quantum Matter, South China Normal
University, Guangzhou 510006, China  Guangdong-Hong Kong Joint Laboratory of
Quantum Matter, Southern Nuclear Science Computing Center, South China Normal
University, Guangzhou 510006, China

###### Abstract

The yields of $X_{cs\bar{c}\bar{s}}$ with its two possible configurations,
i.e., the hadronic molecular state and tetraquark state, for Pb-Pb collisions
at $\sqrt{s_{NN}}=5.02~{}\rm{TeV}$ is studied. A volume effect is found from
the centrality distribution of $X_{cs\bar{c}\bar{s}}$, which could help to
distinguish the inner structure of $X_{cs\bar{c}\bar{s}}$. We also show the
rapidity and the transverse momentum distributions of $X_{cs\bar{c}\bar{s}}$
production as well as its elliptic flow coefficient as a function of the
transverse momentum.

Exotic

##  I Introduction

Quarks and gluons are the fundamental degrees of freedom of quantum
chromodynamics (QCD). Because of the nonperturbative feature of QCD, we can
only observe confined colorless hadrons. A normal hadron has two modes: a
meson is made up of one quark and one antiquark, and a baryon is made up of
three (anti)quarks. Multiquark hadrons made up of more than three quarks were
proposed at the beginning of the construction of the quark model by Gell-Mann
and Zweig Gell-Mann:1964ewy ; Zweig:1964ruk ; Zweig:1964jf . However, the
existence of tetraquarks and pentaquarks was not proven until the observation
of $XYZ$ states ParticleDataGroup:2020ssz , hidden-charm $P_{c}$ states
LHCb:2015yax ; LHCb:2019kea ; LHCb:2020jpq ; LHCb:2021chn , doubly-charm
$T_{cc}^{+}$ LHCb:2021vvq ; LHCb:2021auc  and fully-charm tetraquark states
LHCb:2020bwg , etc. Chen:2016qju ; Esposito:2016noz ; Guo:2017jvc ;
Liu:2019zoy ; Brambilla:2019esw ; Chen:2022asf ; Meng:2022ozq

Five $J/\psi\phi$ structures $X(4140)$, $X(4274)$, $X(4500)$, $X(4685)$ and
$X(4700)$ in the $B^{+}\to J/\psi\phi K^{+}$ decay process were observed by
the LHCb Collaboration LHCb:2016nsl ; LHCb:2016axx ; LHCb:2021uow , CDF
Collaboration CDF:2009jgo ; CDF:2011pep , CMS Collaboration CMS:2013jru , D0
Collaboration D0:2013jvp  and BARAR Collaboration BaBar:2014wwp . $X(4140)$
and $X(4274)$ are considered as the $cs\bar{c}\bar{s}$ tetraquark ground
states, whereas $X(4500)$ and $X(4700)$ are considered as the
$cs\bar{c}\bar{s}$ tetraquark excited states, in various theoretical methods
Chen:2010ze ; Stancu:2009ka ; Chen:2016oma ; Ebert:2008kb ; Lu:2016cwr ;
Wu:2016gas ; Liu:2021xje ; Deng:2017xlb ; Wang:2016gxp . In Ref. Wang:2021ghk
, $X(4685)$ was also considered as the axial vector 2S radial excited
$cs\bar{c}\bar{s}$ tetraquark state. In Refs. Reinders:1984sr ; Shifman:1978bx
; Turkan:2021ome ; Yang:2022zxe , the mass spectra of the S-wave and D-wave
$cs\bar{c}\bar{s}$ tetraquarks in different excitation structures are
calculated using the QCD sum rules method. There have been many theoretical
studies on the inner structure of these $X$’s, such as the molecular states
Liu:2009ei ; Mahajan:2009pj ; Wang:2009ue ; Wang:2009ry ; Albuquerque:2009ak ;
Branz:2009yt ; Ding:2009vd ; Zhang:2009st ; Shen:2010ky ; Liu:2010hf ;
Finazzo:2011he ; He:2011ed ; Hidalgo-Duque:2012rqv ; Wang:2014gwa ; Ma:2014ofa
; Ma:2014zva ; Karliner:2016ith ; MartinezTorres:2016cqv ; Mutuk:2022ckn ,
compact or diquark-antidiquark states Stancu:2009ka ; Drenska:2009cd ;
Ebert:2010zz ; Chen:2010ze ; Vijande:2014cfa ; Patel:2014vua ;
Anisovich:2015caa ; Wang:2015pea ; Zhou:2015frp ; Padmanath:2015era ;
Lebed:2016yvr , cusp effects vanBeveren:2009dc ; Luo:2022xjx ;
Nakamura:2021bvs , dynamically generated resonances Molina:2009ct ;
Branz:2010rj , conventional charmonium Liu:2009iw , and hybrid charmonium
states Mahajan:2009pj ; Wang:2009ue . However, overall, the inner structure of
$X(4140)$, $X(4274)$, $X(4500)$, $X(4685)$ and $X(4700)$ remains an open
question.

In the molecular picture, a $X_{cs\bar{c}\bar{s}}$ is formed by a strange-
charmed meson $D_{s}^{+}$ ($D_{s}^{-}$) and a $D_{s}^{*-}$ ($D_{s}^{*+}$),
while a $X(3872)$ is formed by a charmed meson $D_{0}$ ($D_{0}^{*}$, $D^{+}$,
$D^{-}$) and a $\bar{D_{0}}$ ($\bar{D_{0}^{*}}$, $D^{*-}$, $D^{*+}$). In the
tetraquark picture, a $X_{cs\bar{c}\bar{s}}$ is formed by a spin triplet
diquark $[cs]_{1}$ (spin singlet diquark $[cs]_{0}$) and a spin singlet
antidiquark $[\bar{c}\bar{s}]_{0}$ (spin triplet antidiquark
$[\bar{c}\bar{s}]_{1}$), while a $X(3872)$ is formed by a diquark $[cq]_{1}$
($[cq]_{0}$) and a $[\bar{c}\bar{q}]_{0}$ ($[\bar{c}\bar{q}]_{1}$), $q$ for
$u$ and $d$ quarks. Although light quarks $u$ and $d$ in $X(3872)$ are
replaced with $s$ quarks in $X_{cs\bar{c}\bar{s}}$, their inner structures may
or may not be the same. This motivates the present study, in which we examine
whether the approach we proposed in Ref. Zhang:2020dwn  can also be applied to
the $X_{cs\bar{c}\bar{s}}$ case and thus find a way to distinguish the two
internal structures with heavy ion measurements. In this work, we try to
distinguish the two aforementioned possible inner structures of
$X_{cs\bar{c}\bar{s}}$, i.e., a loose hadronic molecule or a compact
tetraquark, by studying its production in heavy ion collisions. The remainder
of this paper is organized as follows. In section II, we introduce the
generation mechanism of $X_{cs\bar{c}\bar{s}}$ into the AMPT model
corresponding to its two possible inner structures following the production of
$X(3872)$ described in Ref. Zhang:2020dwn . In section III, we examine the
production of $X_{cs\bar{c}\bar{s}}$ as a function of centrality, transverse
momentum, and rapidity. A volume effect is found, which can be a probe of the
inner structure of $X_{cs\bar{c}\bar{s}}$. A summary and outlook are presented
in section IV.

##  II Framework

In this study, we generate a total of one million minimum bias events for Pb-
Pb collisions at $\sqrt{s_{NN}}=5.02~{}\rm{TeV}$ by using the framework
developed in Ref. Zhang:2020dwn . We introduce the production mechanism to
produce $X_{cs\bar{c}\bar{s}}$ for its two possible configurations, i.e., the
hadronic molecular configurations and the tetraquark configurations into the
default version (v1.26t9b) of AMPT transport model Lin:2004en . Given that
$X_{cs\bar{c}\bar{s}}$ contains (anti-)charm quarks and (anti-)strange quarks,
we need to generate a reasonable number of individual charm and strange quarks
in the partonic phase. On top of the default version of AMPT, we modify the
factor of $K$ (Zhang:2019utb ) to enhance the initial $c$ and $\bar{c}$
spectra because of a lack of some channels related to initial heavy quarks.
The AMPT calculation gives a reasonable (order-of-magnitude) description of
the experimental data ALICE:2021rxa  for the total yield of $D^{+}+D^{+*}$ in
the low $p_{T}$ region (see upper panel of Fig. 1). For the strange quarks, an
upper limit on the relative production of strange to non-strange quarks in
AMPT is set to 0.6 because of the strangeness enhancement effect (see
Lin:2014tya ), and our calculations also give a reasonable (order-of-
magnitude) description of the experimental data ALICE:2021kfc  for the yield
of $D_{s}^{+}$ meson (see lower panel of Fig. 1). The main purpose of this
work is to distinguish two inner structures of $X_{cs\bar{c}\bar{s}}$ through
their significantly different production rates. The difference of $D$ and
$D_{s}^{+}$ mesons production between our calculation and experimental data
should not influence the relative yield between two inner structures and thus
cannot change the qualitative results.

Figure 1: (color online) Upper panel: total production of $D^{+}+D^{+*}$ from
the ALICE Collaboration ALICE:2021rxa ; Lower panel: the production of
$D_{s}^{+}$ from the ALICE Collaboration ALICE:2021kfc . The bands reflect the
uncertainty due to constituent composition as discussed around Eq. (1) that
are obtained from varying the composition fraction by $\pm 10\%$.

We use the same production mechanism developed in Ref. Zhang:2020dwn  for the
hadronic molecule and tetraquark configurations of the $X_{cs\bar{c}\bar{s}}$.
For the molecular picture, the charmed-strange mesons are collected after the
hadronization process. Then, $D_{s}^{+}$ ($D_{s}^{-}$) and $D_{s}^{-*}$
($D_{s}^{+*}$) are coalesced (similar to the hadronization process mentioned
in Lin:2004en ) to form the “molecule” $X_{cs\bar{c}\bar{s}}$ according to the
following conditions: the relative distance within the region
$[5\rm{fm},~{}7\rm{fm}]$ and invariant mass within the region
$[2M_{D_{s}^{+}},~{}2M_{D_{s}^{+*}}]$. For the tetraquark picture, the “tetra”
$X_{cs\bar{c}\bar{s}}$ is formed via two steps. (i) First, diquarks ($cs$) and
diquarks ($\bar{c}\bar{s}$) are formed by matching a (anti-)charm quark with
the nearest (in both position space and momentum space) (anti-)strange quark
in the parton. (ii) Then, these (anti)diquarks are coalesced to form the
$X_{cs\bar{c}\bar{s}}$ according to the following conditions: the relative
distance $<1\rm{fm}$ and invariant mass within the region
$[2M_{[cs]_{1}},~{}2M_{[cs]_{0}}]$ (the spin triplet and singlet diquark
masses are defined in Refs. Ebert:2008kb ; Lu:2016cwr ). Owing to a lack of
spin information in the AMPT model for the formation of the charmed-strange
mesons and (anti)diquarks, the relative yield ratios are estimated using the
thermal model:

$R(\frac{A}{B})\equiv\frac{\text{Yield}(A)}{\text{Yield}(B)}=e^{-(m_{A}-m_{B})/T_{\text{freezeout}}},$
(1)

where $m_{A}$ and $m_{B}$ represent the masses of hadrons A and B,
respectively. Here, $T_{\text{freezeout}}\simeq 160~{}\mathrm{MeV}$ is the
freeze-out temperature. For the hadronic picture, A and B are the $D_{s}^{+}$
and $D_{s}^{+*}$ mesons, respectively. For the tetraquark picture, A and B are
the spin triplet and singlet diquark, respectively. This estimate indicates a
composition of $30\%$($70\%$) for $D_{s}^{+}$($D_{s}^{+*}$) and a composition
of $35\%$($65\%$) for spin triplet(singlet) diquarks. We also vary the
composition between $20\%$($80\%$) and $40\%$($60\%$) to show the uncertainty
bands.

##  III Results and Discussions

Within this simulation framework, we use the Monte Carlo method to generate a
total of one million minimum bias events for Pb-Pb collisions at
$\sqrt{s_{NN}}=5.02\ \rm{TeV}$. The inclusive yield of $X_{cs\bar{c}\bar{s}}$
is found to be approximately 42000 in the molecular picture and approximately
200 in the tetraquark picture. As a benchmark for comparison, we also estimate
the yield of $X(3872)$ within the same framework (see the production mechanism
in Ref. Zhang:2020dwn , the yield should be multiplied a factor $\tfrac{1}{4}$
owing to wavefunction normalization for both the molecular and tetraquark
pictures). The inclusive yield of $X(3872)$ is found to be approximately
171000 in the molecular picture and approximately 600 in the tetraquark
picture. The yield of $X_{cs\bar{c}\bar{s}}$ is approximately $\tfrac{1}{4}$
of that of $X(3872)$. Compared with the experimental data of $X(3872)$
measured by the CMS collaboration for Pb-Pb collisions at $\sqrt{s_{NN}}=5.02\
\rm{TeV}$ CMS:2021znk , our finding suggests that an observable signal of
$X_{cs\bar{c}\bar{s}}$ could be measured in heavy ion collisions at the LHC
energy.

One can also find the production in the molecular picture significantly
exceeds that in the tetraquark picture, by a factor of 200 — a 2-order-of-
magnitude difference. This result may be understood as follows: the
$c-\bar{c}$ and $s-\bar{s}$ quarks must be pair produced in the initial
conditions of heavy ion collisions and then expand and cool with the bulk
flow; the molecule $X_{cs\bar{c}\bar{s}}$ needs a large volume to be formed,
while the tetraquark $X_{cs\bar{c}\bar{s}}$ needs a compact volume to be
formed; thus, the probability of the formation of hadron molecules is far
higher than that for the tetraquark state.

Figure 2: (color online) Centrality dependence of the $X_{cs\bar{c}\bar{s}}$
in Pb-Pb collisions at $\sqrt{s_{NN}}=5.02~{}\rm{TeV}$ for hadronic molecular
configuration (red solid boxes) and tetraquark configuration (blue shaded
boxes). The bands reflect both statistical uncertainty from our simulations
and the uncertainty due to constituent composition as discussed around Eq.
(1), and are obtained from varying the composition fraction by $\pm 10\%$.

We plot the $X_{cs\bar{c}\bar{s}}$ production as a function of centrality in
Pb-Pb collisions at $\sqrt{s_{NN}}=5.02~{}\rm{TeV}$ for the hadronic molecular
state and tetraquark state in Fig. 2. One can find the yield of the
$X_{cs\bar{c}\bar{s}}$ in the molecular picture is 2 orders of magnitude
larger than that in the tetraquark picture. From the central collision region
to the peripheral collision region, the production first increases then
decreases for both the molecular state and the tetraquark state, and the slope
of the decrease is far larger in the molecular state than in the tetraquark
state. This results from a competing effect between the volume of the bulk
system and the size of $X_{cs\bar{c}\bar{s}}$. For central collisions, the
number of (anti-)charm and (anti-)strange quarks is large, the bulk volume is
large, and its evolution time is long; thus, the (anti-)charm and
(anti-)strange quarks separate sufficiently, which benefits the production of
a large-size molecular state. For the peripheral collisions, both the number
of (anti-)charm and (anti-)strange quarks and the size of the fireball are
small; as such, the evolution time of the fireball is short, which benefits
the production of small-sized tetraquark states. This size effect could help
to explore the internal structure of $X_{cs\bar{c}\bar{s}}$ through different
collision systems, e.g., Pb-Pb, Au-Au, Xe-Xe, Cu-Cu, O-O, and $d-A$/$p-A$.

Figure 3: (color online) Rapidity $y$ and transverse momentum $p_{T}$
distribution of the $X_{cs\bar{c}\bar{s}}$ yield in Pb-Pb collisions at
$\sqrt{s_{NN}}=5.02~{}\mathrm{TeV}$ for hadronic molecular configuration (red
solid boxes) and tetraquark configuration (blue shaded boxes). The bands are
determined as described in Fig. 2.

In Fig. 3, we present the rapidity and the transverse momentum distributions
of $X_{cs\bar{c}\bar{s}}$. One can find that the distribution for both the
hadronic molecular state and the tetraquark state is similar to that of the
usual hadrons CMS:2011aqh ; ALICE:2013jfw . We also show the elliptic flow
coefficient $v_{2}$ of $X_{cs\bar{c}\bar{s}}$ as a function of the transverse
momentum $p_{T}$ in Fig. 4. The elliptic flow is sensitive to the geometry of
the initial fireball and the generation mechanism of $X_{cs\bar{c}\bar{s}}$.

Figure 4: (color online) Elliptic flow coefficient $v_{2}$ versus transverse
momentum $p_{T}$ for produced $X_{cs\bar{c}\bar{s}}$ in minimum bias Pb-Pb
collisions at $\sqrt{s_{NN}}=5.02~{}\mathrm{TeV}$, predicted from our
computation for the hadronic molecule picture. The bands are similarly
determined as described in Fig. 2.

##  IV Summary and Outlook

In this work, we studied the yields of $X_{cs\bar{c}\bar{s}}$ for Pb-Pb
collisions at $\sqrt{s_{NN}}=5.02~{}\rm{TeV}$ by introducing the production
mechanism of two possible configurations, i.e., the hadronic molecular state
and tetraquark state, into the AMPT model. We found that the production in the
molecular picture exceeds in the tetraquark picture by two orders of
magnitude. The centrality distribution of the yields of $X_{cs\bar{c}\bar{s}}$
shows a strongly decreasing trend for the hadronic molecular state and a mild
change for the tetraquark state. This system size dependence could be a good
probe for the inner structure of $X_{cs\bar{c}\bar{s}}$. We also showed the
rapidity and the transverse momentum distributions of $X_{cs\bar{c}\bar{s}}$
production, as well as its elliptic flow coefficient, as a function of the
transverse momentum, which can be tested in the future experimental
measurements. In Ref. ALICE:2021kfc , a strangeness enhancement effect in
heavy ion collisions was found by ALICE Collaboration, which could be evidence
for quark-gluon plasma. We expect a similar effect to be found in the ratio of
$X_{cs\bar{c}\bar{s}}$ to $X(3872)$, which will be studied in our future work.

###### Acknowledgements.

The authors would like to thank Dr. J. Liao, E. Wang, Q. Wang, and H. Xing for
the helpful discussion. This work is partly supported by Guangdong Major
Project of Basic and Applied Basic Research No. 2020B0301030008, the National
Natural Science Foundation of China with Grant No. 12105107, Science and
Technology Program of Guangzhou No. 2019050001.

## References

  * (1) M. Gell-Mann, Phys. Lett. 8, 214-215 (1964) doi:10.1016/S0031-9163(64)92001-3
  * (2) G. Zweig, CERN-TH-401.
  * (3) G. Zweig, CERN-TH-412.
  * (4) P. A. Zyla et al. [Particle Data Group], PTEP 2020, no.8, 083C01 (2020) doi:10.1093/ptep/ptaa104
  * (5) R. Aaij et al. [LHCb], Phys. Rev. Lett. 115, 072001 (2015) doi:10.1103/PhysRevLett.115.072001 [arXiv:1507.03414 [hep-ex]].
  * (6) R. Aaij et al. [LHCb], Phys. Rev. Lett. 122, no.22, 222001 (2019) doi:10.1103/PhysRevLett.122.222001 [arXiv:1904.03947 [hep-ex]].
  * (7) R. Aaij et al. [LHCb], Sci. Bull. 66, 1278-1287 (2021) doi:10.1016/j.scib.2021.02.030 [arXiv:2012.10380 [hep-ex]].
  * (8) R. Aaij et al. [LHCb], Phys. Rev. Lett. 128, no.6, 062001 (2022) doi:10.1103/PhysRevLett.128.062001 [arXiv:2108.04720 [hep-ex]].
  * (9) R. Aaij et al. [LHCb], Nature Phys. 18, no.7, 751-754 (2022) doi:10.1038/s41567-022-01614-y [arXiv:2109.01038 [hep-ex]].
  * (10) R. Aaij et al. [LHCb], Nature Commun. 13, no.1, 3351 (2022) doi:10.1038/s41467-022-30206-w [arXiv:2109.01056 [hep-ex]].
  * (11) R. Aaij et al. [LHCb], Sci. Bull. 65, no.23, 1983-1993 (2020) doi:10.1016/j.scib.2020.08.032 [arXiv:2006.16957 [hep-ex]].
  * (12) H. X. Chen, W. Chen, X. Liu and S. L. Zhu, Phys. Rept. 639, 1-121 (2016) doi:10.1016/j.physrep.2016.05.004 [arXiv:1601.02092 [hep-ph]].
  * (13) A. Esposito, A. Pilloni and A. D. Polosa, Phys. Rept. 668, 1-97 (2017) doi:10.1016/j.physrep.2016.11.002 [arXiv:1611.07920 [hep-ph]].
  * (14) F. K. Guo, C. Hanhart, U. G. Meißner, Q. Wang, Q. Zhao and B. S. Zou, Rev. Mod. Phys. 90, no.1, 015004 (2018) [erratum: Rev. Mod. Phys. 94, no.2, 029901 (2022)] doi:10.1103/RevModPhys.90.015004 [arXiv:1705.00141 [hep-ph]].
  * (15) Y. R. Liu, H. X. Chen, W. Chen, X. Liu and S. L. Zhu, Prog. Part. Nucl. Phys. 107, 237-320 (2019) doi:10.1016/j.ppnp.2019.04.003 [arXiv:1903.11976 [hep-ph]].
  * (16) N. Brambilla, S. Eidelman, C. Hanhart, A. Nefediev, C. P. Shen, C. E. Thomas, A. Vairo and C. Z. Yuan, Phys. Rept. 873, 1-154 (2020) doi:10.1016/j.physrep.2020.05.001 [arXiv:1907.07583 [hep-ex]].
  * (17) H. X. Chen, W. Chen, X. Liu, Y. R. Liu and S. L. Zhu, Rept. Prog. Phys. 86, no.2, 026201 (2023) doi:10.1088/1361-6633/aca3b6 [arXiv:2204.02649 [hep-ph]].
  * (18) L. Meng, B. Wang, G. J. Wang and S. L. Zhu, [arXiv:2204.08716 [hep-ph]].
  * (19) R. Aaij et al. [LHCb], Phys. Rev. D 95, no.1, 012002 (2017) doi:10.1103/PhysRevD.95.012002 [arXiv:1606.07898 [hep-ex]].
  * (20) R. Aaij et al. [LHCb], Phys. Rev. Lett. 118, no.2, 022003 (2017) doi:10.1103/PhysRevLett.118.022003 [arXiv:1606.07895 [hep-ex]].
  * (21) R. Aaij et al. [LHCb], Phys. Rev. Lett. 127, no.8, 082001 (2021) doi:10.1103/PhysRevLett.127.082001 [arXiv:2103.01803 [hep-ex]].
  * (22) T. Aaltonen et al. [CDF], Phys. Rev. Lett. 102, 242002 (2009) doi:10.1103/PhysRevLett.102.242002 [arXiv:0903.2229 [hep-ex]].
  * (23) T. Aaltonen et al. [CDF], Mod. Phys. Lett. A 32, no.26, 1750139 (2017) doi:10.1142/S0217732317501395 [arXiv:1101.6058 [hep-ex]].
  * (24) S. Chatrchyan et al. [CMS], Phys. Lett. B 734, 261-281 (2014) doi:10.1016/j.physletb.2014.05.055 [arXiv:1309.6920 [hep-ex]].
  * (25) V. M. Abazov et al. [D0], Phys. Rev. D 89, no.1, 012004 (2014) doi:10.1103/PhysRevD.89.012004 [arXiv:1309.6580 [hep-ex]].
  * (26) J. P. Lees et al. [BaBar], Phys. Rev. D 91, no.1, 012003 (2015) doi:10.1103/PhysRevD.91.012003 [arXiv:1407.7244 [hep-ex]].
  * (27) W. Chen and S. L. Zhu, Phys. Rev. D 83, 034010 (2011) doi:10.1103/PhysRevD.83.034010 [arXiv:1010.3397 [hep-ph]].
  * (28) F. Stancu, J. Phys. G 37, 075017 (2010) [erratum: J. Phys. G 46, no.1, 019501 (2019)] doi:10.1088/0954-3899/37/7/075017 [arXiv:0906.2485 [hep-ph]].
  * (29) H. X. Chen, E. L. Cui, W. Chen, X. Liu and S. L. Zhu, Eur. Phys. J. C 77, no.3, 160 (2017) doi:10.1140/epjc/s10052-017-4737-5 [arXiv:1606.03179 [hep-ph]].
  * (30) D. Ebert, R. N. Faustov and V. O. Galkin, Eur. Phys. J. C 58, 399-405 (2008) doi:10.1140/epjc/s10052-008-0754-8 [arXiv:0808.3912 [hep-ph]].
  * (31) Q. F. Lü and Y. B. Dong, Phys. Rev. D 94, no.7, 074007 (2016) doi:10.1103/PhysRevD.94.074007 [arXiv:1607.05570 [hep-ph]].
  * (32) J. Wu, Y. R. Liu, K. Chen, X. Liu and S. L. Zhu, Phys. Rev. D 94, no.9, 094031 (2016) doi:10.1103/PhysRevD.94.094031 [arXiv:1608.07900 [hep-ph]].
  * (33) X. Liu, H. Huang, J. Ping, D. Chen and X. Zhu, Eur. Phys. J. C 81, no.10, 950 (2021) doi:10.1140/epjc/s10052-021-09752-y [arXiv:2103.12425 [hep-ph]].
  * (34) C. Deng, J. Ping, H. Huang and F. Wang, Phys. Rev. D 98, no.1, 014026 (2018) doi:10.1103/PhysRevD.98.014026 [arXiv:1801.00164 [hep-ph]].
  * (35) Z. G. Wang, Eur. Phys. J. C 77, no.2, 78 (2017) doi:10.1140/epjc/s10052-017-4640-0 [arXiv:1606.05872 [hep-ph]].
  * (36) Z. G. Wang, Adv. High Energy Phys. 2021, 4426163 (2021) doi:10.1155/2021/4426163 [arXiv:2103.04236 [hep-ph]].
  * (37) L. J. Reinders, H. Rubinstein and S. Yazaki, Phys. Rept. 127, 1 (1985) doi:10.1016/0370-1573(85)90065-1
  * (38) M. A. Shifman, A. I. Vainshtein and V. I. Zakharov, Nucl. Phys. B 147, 385-447 (1979) doi:10.1016/0550-3213(79)90022-1
  * (39) A. Türkan, J. Y. Süngü and E. V. Veliev, [arXiv:2103.05515 [hep-ph]].
  * (40) Z. Y. Yang and W. Chen, [arXiv:2206.06051 [hep-ph]].
  * (41) X. Liu and S. L. Zhu, Phys. Rev. D 80, 017502 (2009) [erratum: Phys. Rev. D 85, 019902 (2012)] doi:10.1103/PhysRevD.85.019902 [arXiv:0903.2529 [hep-ph]].
  * (42) N. Mahajan, Phys. Lett. B 679, 228-230 (2009) doi:10.1016/j.physletb.2009.07.043 [arXiv:0903.3107 [hep-ph]].
  * (43) Z. G. Wang, Eur. Phys. J. C 63, 115-122 (2009) doi:10.1140/epjc/s10052-009-1097-9 [arXiv:0903.5200 [hep-ph]].
  * (44) Z. G. Wang, Z. C. Liu and X. H. Zhang, Eur. Phys. J. C 64, 373-386 (2009) doi:10.1140/epjc/s10052-009-1156-2 [arXiv:0907.1467 [hep-ph]].
  * (45) R. M. Albuquerque, M. E. Bracco and M. Nielsen, Phys. Lett. B 678, 186-190 (2009) doi:10.1016/j.physletb.2009.06.022 [arXiv:0903.5540 [hep-ph]].
  * (46) T. Branz, T. Gutsche and V. E. Lyubovitskij, Phys. Rev. D 80, 054019 (2009) doi:10.1103/PhysRevD.80.054019 [arXiv:0903.5424 [hep-ph]].
  * (47) G. J. Ding, Eur. Phys. J. C 64, 297-308 (2009) doi:10.1140/epjc/s10052-009-1146-4 [arXiv:0904.1782 [hep-ph]].
  * (48) J. R. Zhang and M. Q. Huang, J. Phys. G 37, 025005 (2010) doi:10.1088/0954-3899/37/2/025005 [arXiv:0905.4178 [hep-ph]].
  * (49) L. L. Shen, X. L. Chen, Z. G. Luo, P. Z. Huang, S. L. Zhu, P. F. Yu and X. Liu, Eur. Phys. J. C 70, 183-217 (2010) doi:10.1140/epjc/s10052-010-1441-0 [arXiv:1005.0994 [hep-ph]].
  * (50) X. Liu, Z. G. Luo and S. L. Zhu, Phys. Lett. B 699, 341-344 (2011) [erratum: Phys. Lett. B 707, 577 (2012)] doi:10.1016/j.physletb.2011.12.019 [arXiv:1011.1045 [hep-ph]].
  * (51) S. I. Finazzo, M. Nielsen and X. Liu, Phys. Lett. B 701, 101-106 (2011) doi:10.1016/j.physletb.2011.05.042 [arXiv:1102.2347 [hep-ph]].
  * (52) J. He and X. Liu, Eur. Phys. J. C 72, 1986 (2012) doi:10.1140/epjc/s10052-012-1986-1 [arXiv:1102.1127 [hep-ph]].
  * (53) C. Hidalgo-Duque, J. Nieves and M. P. Valderrama, Phys. Rev. D 87, no.7, 076006 (2013) doi:10.1103/PhysRevD.87.076006 [arXiv:1210.5431 [hep-ph]].
  * (54) Z. G. Wang, Eur. Phys. J. C 74, no.7, 2963 (2014) doi:10.1140/epjc/s10052-014-2963-7 [arXiv:1403.0810 [hep-ph]].
  * (55) L. Ma, Z. F. Sun, X. H. Liu, W. Z. Deng, X. Liu and S. L. Zhu, Phys. Rev. D 90, no.3, 034020 (2014) doi:10.1103/PhysRevD.90.034020 [arXiv:1403.7907 [hep-ph]].
  * (56) L. Ma, X. H. Liu, X. Liu and S. L. Zhu, Phys. Rev. D 91, no.3, 034032 (2015) doi:10.1103/PhysRevD.91.034032 [arXiv:1406.6879 [hep-ph]].
  * (57) M. Karliner and J. L. Rosner, Nucl. Phys. A 954, 365-370 (2016) doi:10.1016/j.nuclphysa.2016.03.057 [arXiv:1601.00565 [hep-ph]].
  * (58) A. Martínez Torres, K. P. Khemchandani, J. M. Dias, F. S. Navarra and M. Nielsen, Nucl. Phys. A 966, 135-157 (2017) doi:10.1016/j.nuclphysa.2017.06.022 [arXiv:1606.07505 [hep-ph]].
  * (59) H. Mutuk, Eur. Phys. J. C 82, no.12, 1142 (2022) doi:10.1140/epjc/s10052-022-11120-3 [arXiv:2211.14836 [hep-ph]].
  * (60) N. V. Drenska, R. Faccini and A. D. Polosa, Phys. Rev. D 79, 077502 (2009) doi:10.1103/PhysRevD.79.077502 [arXiv:0902.2803 [hep-ph]].
  * (61) D. Ebert, R. N. Faustov and V. O. Galkin, Phys. Part. Nucl. 41, 931-934 (2010) doi:10.1134/S1063779610060225
  * (62) J. Vijande and A. Valcarce, Phys. Lett. B 736, 325-328 (2014) doi:10.1016/j.physletb.2014.07.048 [arXiv:1507.08045 [hep-ph]].
  * (63) S. Patel, M. Shah and P. C. Vinodkumar, Eur. Phys. J. A 50, 131 (2014) doi:10.1140/epja/i2014-14131-9 [arXiv:1402.3974 [hep-ph]].
  * (64) V. V. Anisovich, M. A. Matveev, A. V. Sarantsev and A. N. Semenova, Int. J. Mod. Phys. A 30, no.32, 1550186 (2015) doi:10.1142/S0217751X15501869 [arXiv:1507.07232 [hep-ph]].
  * (65) Z. g. Wang and Y. f. Tian, Int. J. Mod. Phys. A 30, 1550004 (2015) doi:10.1142/S0217751X15500049 [arXiv:1502.04619 [hep-ph]].
  * (66) P. Zhou, C. R. Deng and J. L. Ping, Chin. Phys. Lett. 32, no.10, 101201 (2015) doi:10.1088/0256-307X/32/10/101201
  * (67) M. Padmanath, C. B. Lang and S. Prelovsek, Phys. Rev. D 92, no.3, 034501 (2015) doi:10.1103/PhysRevD.92.034501 [arXiv:1503.03257 [hep-lat]].
  * (68) R. F. Lebed and A. D. Polosa, Phys. Rev. D 93, no.9, 094024 (2016) doi:10.1103/PhysRevD.93.094024 [arXiv:1602.08421 [hep-ph]].
  * (69) E. van Beveren and G. Rupp, [arXiv:0906.2278 [hep-ph]].
  * (70) X. Luo and S. X. Nakamura, Phys. Rev. D 107, no.1, 1 (2023) doi:10.1103/PhysRevD.107.L011504 [arXiv:2207.12875 [hep-ph]].
  * (71) S. X. Nakamura, Phys. Lett. B 834, 137486 (2022) doi:10.1016/j.physletb.2022.137486 [arXiv:2111.05115 [hep-ph]].
  * (72) R. Molina and E. Oset, Phys. Rev. D 80, 114013 (2009) doi:10.1103/PhysRevD.80.114013 [arXiv:0907.3043 [hep-ph]].
  * (73) T. Branz, R. Molina and E. Oset, Phys. Rev. D 83, 114015 (2011) doi:10.1103/PhysRevD.83.114015 [arXiv:1010.0587 [hep-ph]].
  * (74) X. Liu, Phys. Lett. B 680, 137-140 (2009) doi:10.1016/j.physletb.2009.08.049 [arXiv:0904.0136 [hep-ph]].
  * (75) H. Zhang, J. Liao, E. Wang, Q. Wang and H. Xing, Phys. Rev. Lett. 126, no.1, 012301 (2021) [arXiv:2004.00024 [hep-ph]].
  * (76) Z. W. Lin, C. M. Ko, B. A. Li, B. Zhang and S. Pal, Phys. Rev. C 72, 064901 (2005) doi:10.1103/PhysRevC.72.064901 [arXiv:nucl-th/0411110 [nucl-th]].
  * (77) C. Zhang, L. Zheng, F. Liu, S. Shi and Z. W. Lin, Phys. Rev. C 99, no.6, 064906 (2019) doi:10.1103/PhysRevC.99.064906 [arXiv:1903.03292 [nucl-th]].
  * (78) S. Acharya et al. [ALICE], JHEP 01, 174 (2022) doi:10.1007/JHEP01(2022)174 [arXiv:2110.09420 [nucl-ex]].
  * (79) Z. W. Lin, Phys. Rev. C 90, no.1, 014904 (2014) doi:10.1103/PhysRevC.90.014904 [arXiv:1403.6321 [nucl-th]].
  * (80) S. Acharya et al. [ALICE], Phys. Lett. B 827, 136986 (2022) doi:10.1016/j.physletb.2022.136986 [arXiv:2110.10006 [nucl-ex]].
  * (81) A. M. Sirunyan et al. [CMS], Phys. Rev. Lett. 128, no.3, 032001 (2022) doi:10.1103/PhysRevLett.128.032001 [arXiv:2102.13048 [hep-ex]].
  * (82) S. Chatrchyan et al. [CMS], JHEP 08, 141 (2011) [arXiv:1107.4800 [nucl-ex]].
  * (83) E. Abbas et al. [ALICE], Phys. Lett. B 726, 610-622 (2013) [arXiv:1304.0347 [nucl-ex]].


L. Angela Mihai, School of Mathematics, Cardiff University, Cardiff, UK.

# A theoretical model for power generation via liquid crystal elastomers

L. Angela Mihai11affiliationmark:  11affiliationmark:  School of Mathematics,
Cardiff University, Cardiff, UK<EMAIL_ADDRESS>
###### Abstract

Motivated by the need for new materials and green energy production and
conversion processes, a class of mathematical models for liquid crystal
elastomers integrated within a theoretical charge pump electrical circuit is
considered. The charge pump generates higher voltage from a lower voltage
supplied by a battery and recharges the battery by harnessing the chemical and
mechanical properties of liquid crystal elastomers transitioning from the
nematic to isotropic phase when illuminated or heated. For the material
constitutive model, purely elastic and neoclassical-type strain energy
densities, applicable to a wide range of monodomain nematic elastomers, are
combined with a modified Maier-Saupe mean field model for photo-thermal
responses. By varying the model parameters of the elastic and neoclassical
terms, it is found that liquid crystal elastomers are more effective than
rubber when used as dielectric material within a charge pump capacitor.

###### keywords:

liquid crystal elastomers, actuation, parallel plate capacitor, nonlinear
elasticity, mathematical modelling

##  1 Introduction

To end the use of fossil fuels, more materials that enable green energy
production and conversion processes are sought and developed Gallo:2016:etal ;
Renewables:2022 . In particular, flexible energy harvesters made of rubber-
like materials demonstrate great potential for generating low carbon renewable
energy in emerging technologies Lu:2020:etal ; Moretti:2020:etal .

This paper considers a liquid crystal elastomer (LCE) incorporated in a
theoretical charge pump electrical circuit. Charge pumps are convenient and
economical devices that use capacitors to generate higher voltages from a
lower voltage supplied by a source battery. The simplest capacitor consists of
two parallel plate electrical conductors separated by air or an insulating
material known as the dielectric. The plates are connected to two terminals,
which can be wired into an electric circuit. When the performance of a
capacitor changes by altering the distance between plates or the amount of
plate surface area, a variable capacitor is obtained.

Figure 1: Schematic of parallel plate capacitor with LCE dielectric between
two compliant electrodes. Figure 2: Natural deformation of LCE dielectric
under the nematic-isotropic phase transition caused by temperature changes.
The thickness $d=\lambda d_{0}$ and surface area $A$ of the LCE change, while
the volume remains constant, $\Omega=A_{0}d_{0}=Ad$. In the nematic phase, the
director n is aligned parallel to the surface, in the first Cartesian
direction, while in the isotropic phase, there is no director alignment.

LCEs are top candidates for dielectric material because they are capable of
large strain deformations which are reversible and repeatable under external
stimuli, such as heat or light deJeu:2012 ; Warner:2007:WT . This is due to
their unique molecular architecture combining the flexibility of polymeric
networks with liquid crystal self-organisation deGennes:1993:dGP . A capacitor
with LCE dielectric between two compliant electrodes is presented
schematically in Figure 1. An example of LCE nematic-isotropic phase
transition under thermal stimuli is depicted in Figure 2. Light-induced shape
changes for LCEs containing photoisomerizing dye molecules can be represented
in a similar manner. In Hiscock:2011:HWPM , a hypothetical charge pump which
converts solar heat into DC electricity was proposed where the LCE was
described by the neoclassical model Bladon:1994:BTW ; Warner:1988:WGV ;
Warner:1991:WW .

To model the LCE, here, the purely elastic and neoclassical-type strain energy
densities are combined with a modified Maier-Saupe mean field model
Bai:2020:BB ; Corbett:2006:CW ; Corbett:2008:CW ; Maier:1959:MS  for
optothermal stimulation. The resulting composite model is applicable to a wide
range of nematic elastomers and can be reduced to either the neo-Hookean model
for rubber Treloar:1944  or the neoclassical model for LCEs. By varying the
model parameters of the elastic and neoclassical terms, it is shown that LCEs
are more effective than rubber when used as dielectric material within a
charge pump capacitor.

##  2 The charge pump circuit

The electrical energy potential stored by a capacitor, known as capacitance,
is equal to $C=q/V$, where $q$ is the magnitude of the charge stored when the
voltage across the capacitor is equal to $V$, and is measured in farads (F).
Capacitance depends on both the geometry and the materials that the capacitor
is made of.

Figure 3: Charge pump electrical circuit with a single-cell supply battery of
lower voltage $V_{1}$ and a variable capacitor with capacitance $C$ that
attains a higher voltage $V_{2}$ Hiscock:2011:HWPM .

Figure 4: Charge pump electrical circuit with a single-cell supply battery of
lower voltage $V_{1}$ and a variable capacitor with capacitance $C$ that
attains a higher voltage $V_{2}$ and recharges the supply battery
Hiscock:2011:HWPM .  Figure 5: Operating cycle of variable capacitor with the
dielectric made of LCE generating a higher voltage $V_{2}$ from a lower
voltage $V_{1}$ Hiscock:2011:HWPM .

Figure 4 represents schematically a charge pump electrical circuit where the
voltage $V_{1}$ supplied by an external source (battery) is raised to voltage
$V_{2}$ by using a variable capacitor with capacitance $C$. The diodes prevent
backflow of charge and act as voltage activated switches. The alternative
circuit illustrated in Figure 4 allows for the supply battery to be recharged.
Figure 5 shows the operating cycle of a variable capacitor with LCE dielectric
containing the following states Hiscock:2011:HWPM :

  * (A)

The LCE is at the lowest temperature corresponding to the nematic state,
$T_{n}=T_{A}$, the input voltage from a supply battery is $V_{1}$ and the
capacitor is charged to the initial charge $q_{1}$. At this state, the
capacitance is equal to $C_{n}=C_{A}=q_{1}/V_{1}$;

  * (B)

The temperature rises to $T_{B}>T_{A}$, so the capacitance decreases to
$C_{B}<C_{A}$, while the charge remains equal to $q_{1}$. Thus the voltage
across the capacitor increases to $V_{2}=q_{1}/C_{B}>q_{1}/C_{A}=V_{1}$;

  * (C)

The temperature continues to increase to the isotropic state,
$T_{i}=T_{C}>T_{B}$, hence the capacitance further decrease to
$C_{i}=C_{C}<C_{B}$, but the voltage remains equal to $V_{2}$, so the charge
decreases to $q_{2}=C_{C}V_{2}<C_{B}V_{2}=q_{1}$;

  * (D)

The temperature drops to $T_{D}<T_{i}$, while the charge remains equal to
$q_{2}$ and the capacitance increases to
$C_{D}=q_{2}/V_{1}>q_{2}/V_{2}=C_{C}$. Upon further cooling to the initial
temperature $T_{n}$, the battery with voltage $V_{1}$ charges the capacitor to
the initial charge state $q_{1}$ and the cycle can be repeated.

During one cycle, the external source provides an electrical energy
$W_{in}=V_{1}\left(q_{1}-q_{2}\right)$ and produces
$W_{out}=V_{2}\left(q_{1}-q_{2}\right)$. Thus the net output generated by this
cycle is equal to

$W=W_{out}-W_{in}=\left(V_{2}-V_{1}\right)\left(q_{1}-q_{2}\right)=-C_{n}V_{1}^{2}\left(\frac{C_{i}}{C_{B}}-1\right)\left(\frac{C_{n}}{C_{B}}-1\right).$
(1)

By defining the capacitance ratio

$\xi=\frac{C_{n}}{C_{i}}>1,$ (2)

the maximum generated output per cycle is equal to

$W_{m}(\xi)=C_{n}V_{1}^{2}\frac{(\xi-1)^{2}}{4\xi}.$ (3)

This is attained when

$C_{B}=\frac{2C_{n}C_{i}}{C_{n}+C_{i}},$ (4)

or equivalently, when

$\frac{V_{2}}{V_{1}}=\frac{\xi+1}{2}.$ (5)

For the LCE transitioning from a nematic to an isotropic state and vice versa,
changes in light instead of temperature can be used.

##  3 The LCE strain-energy function

To describe the incompressible nematic LCE, we consider the following elastic
strain-energy function Conti:2002a:CdSD ; Fried:2004:FS ; Fried:2005:FS ;
Fried:2006:FS ; Lee:2021:LB ; Mihai:2020a:MG ; Mihai:2021a:MG ,

$\mathcal{W}^{(el)}(\textbf{F},\textbf{Q},\textbf{n})=\mathcal{W}^{(1)}\left(\textbf{F}\right)+\mathcal{W}^{(2)}\left(\textbf{G}^{-1}\textbf{F}\textbf{G}_{0}\right),$
(6)

where F denotes the deformation gradient from the reference cross-linking
state, such that $\det\textbf{F}=1$, and n is a unit vector for the localised
direction of uniaxial nematic alignment in the present configuration, and is
known as the director. The first term on the right-hand side of equation (6)
represents the strain-energy density associated with the overall macroscopic
deformation, and the second term is the strain-energy density of the polymer
microstructure. In the second term, $\textbf{G}_{0}$ and G denote the natural
(or spontaneous) deformation tensor in the reference and current
configuration, respectively Finkelmann:2001:FGW ; Warner:2007:WT .

Assuming that the LCE is intrinsically uniaxial, the natural deformation
tensor takes the form

$\textbf{G}=a^{-1/6}\textbf{I}+\left(a^{1/3}-a^{-1/6}\right)\textbf{n}\otimes\textbf{n},$
(7)

where $\textbf{I}=\mathrm{diag}(1,1,1)$ is the identity tensor and

$a=\frac{1+2Q}{1-Q}$ (8)

denotes the natural shape parameter, with $Q$ representing the scalar uniaxial
order of the liquid crystal mesogens ($Q=1$ corresponds to perfect nematic
order, while $Q=0$ is for the case when the mesogens are randomly oriented).
In the reference configuration, G is replaced by $\textbf{G}_{0}$, with
$\textbf{n}_{0}$, $a_{0}$ and $Q_{0}$ instead of n, $a$ and $Q$, respectively.

We set the nematic director in the reference and current configuration as
$\textbf{n}_{0}=\textbf{n}=[1,0,0]^{\mathrm{T}}$ and define the associated
natural deformation tensors as follows,

$\textbf{G}_{0}=\mathrm{diag}\left(a_{0}^{1/3},a_{0}^{-1/6},a_{0}^{-1/6}\right),\qquad\textbf{G}=\mathrm{diag}\left(a^{1/3},a^{-1/6},a^{-1/6}\right).$
(9)

The components of the LCE model are defined by the neo-Hookean-type strain-
energy density function Treloar:1944 ,

$\mathcal{W}^{(1)}(\lambda_{1},\lambda_{2},\lambda_{3})=\frac{\mu^{(1)}}{2}\left(\lambda_{1}^{2}+\lambda_{2}^{2}+\lambda_{3}^{2}\right),$
(10)

where $\mu^{(1)}>0$ is a constant independent of the deformation and
$\\{\lambda_{1}^{2},\lambda_{2}^{2},\lambda_{3}^{2}\\}$ are the eigenvalues of
the tensor $\textbf{F}^{T}\textbf{F}$, such that
$\lambda_{1}\lambda_{2}\lambda_{3}=1$, and

$\mathcal{W}^{(2)}(\alpha_{1},\alpha_{2},\alpha_{3})=\frac{\mu^{(2)}}{2}\left(\alpha_{1}^{2}+\alpha_{2}^{2}+\alpha_{3}^{2}\right),$
(11)

where $\mu^{(2)}>0$ is a constant independent of the deformation and
$\\{\alpha_{1}^{2},\alpha_{2}^{2},\alpha_{3}^{2}\\}$ are the eigenvalues of
the elastic Cauchy-Green tensor $\textbf{A}^{T}\textbf{A}$, with the local
elastic deformation tensor
$\textbf{A}=\textbf{G}^{-1}\textbf{F}\textbf{G}_{0}$, such that
$\alpha_{1}\alpha_{2}\alpha_{3}=1$.

The composite model defined by equation (6) then takes the form

$\mathcal{W}^{(el)}=\frac{\mu^{(1)}}{2}\left(\lambda_{1}^{2}+\lambda_{2}^{2}+\lambda_{3}^{2}\right)+\frac{\mu^{(2)}}{2}\left(\alpha_{1}^{2}+\alpha_{2}^{2}+\alpha_{3}^{2}\right).$
(12)

The strain-energy function described by equation (12) reduces to the
neoclassical model for LCEs when $\mu^{(1)}=0$ Bladon:1994:BTW ;
Warner:1988:WGV ; Warner:1991:WW  and to the neo-Hookean model for rubber when
$\mu^{(2)}=0$ Treloar:1944 .

###  3.1 Photo-thermal responses

When azobenzene mesogens are embedded in the nematic elastomeric network,
then, if photons are absorbed, the so-called Weigert effect Ebralidze:1995 ;
Kakichashvili:1974  occurs where the dye molecules change from straight
_trans_ \- to bent _cis_ -isomers, causing a reduction in the nematic order.
To account for photo-thermal deformation of a dielectric LCE, we adopt the
following modified Maier-Saupe mean field model Bai:2020:BB ; Corbett:2006:CW
; Corbett:2008:CW ; Maier:1959:MS ,

$\mathcal{W}^{(ms)}(Q,\textbf{n})=\mu_{n}\left(1-c(Q,\textbf{n})\right)\left[g^{-1}(Q)Q-\log
Z(Q)-\frac{\tilde{J}}{2}Q^{2}\left(1-c(Q,\textbf{n})\right)\right],$ (13)

where $Q$ and n are defined as before, $\mu_{n}=N_{n}kT$, with $N_{n}$ the
total number of mesogens per unit volume and $kT$ the temperature per unit of
energy, $\tilde{J}$ represents the average interaction between two mesogens in
the unit of energy, and $c(\cdot,\cdot)$ denotes the fractional number of
_cis_ molecules.

  * •

If the light is _polarised_ and $\alpha$ denotes the angle between nematic
director n and the light polarisation, then

$c(Q,\textbf{n})=f\frac{I\left[1+Q\left(3\cos^{2}\alpha-1\right)\right]}{3+I\left[1+Q\left(3\cos^{2}\alpha-1\right)\right]},$
(14)

where $f$ is the fraction of photo-active mesogens and $I$ is the non-
dimensional homogeneous light intensity.

  * •

If the light is _unpolarised_ and $\alpha$ is the angle between nematic
director n and the light beam direction, then

$c(Q,\textbf{n})=f\frac{I\left[1-(Q/2)\left(3\cos^{2}\alpha-1\right)\right]}{3+I\left[1-(Q/2)\left(3\cos^{2}\alpha-1\right)\right]}.$
(15)

The expressions for the functions $g(\cdot)$ and $Z(\cdot)$ are, respectively,

$g(x)=-\frac{1}{2}-\frac{1}{2x}+\frac{1}{2x}\sqrt{\frac{3x}{2}}\frac{\exp\left(3x/2\right)}{\int_{0}^{\sqrt{\left(3x/2\right)}}\exp\left(y^{2}\right)\text{d}y}$
(16)

and

$Z(Q)=\frac{\exp\left(g^{-1}(Q)\right)}{1+g^{-1}(Q)\left(1+2Q\right)}.$ (17)

When there is no light, $I=0$ and the energy function defined by (13) reduces
to

$\mathcal{W}^{(ms)}(Q,\textbf{n})=\mu_{n}\left(g^{-1}(Q)Q-\log
Z(Q)-\frac{\tilde{J}}{2}Q^{2}\right).$ (18)

In addition, the ratio between photoactive and non-photoactive mesogens can
also be taken into account Corbett:2006:CW ; Corbett:2008:CW ;
Goriely:2022:GMM .

##  4 Energy conversion

In our analysis, at state (A), the LCE dielectric is considered either in its
natural configuration or pre-stretched perpendicular or parallel to the
nematic director, so that the surface area is increased and the distance
between plates is reduced, hence the initial capacitance increases. In all
these cases, the LCE is then actuated by illumination or heating.

Assuming $\mu\ll\mu_{n}$, the energy functions $\mathcal{W}^{(el)}$ and
$\mathcal{W}^{(ms)}$ described by equations (12) and (13), respectively, can
be treated separately, and the equilibrium scalar order parameter $Q$ obtained
by minimising the function $\mathcal{W}^{(ms)}$. Experiments on photo-thermal
shape changes in LCEs were reported in Finkelmann:2001:FNPW ; Yu:2003:YN . A
general theoretical model for photo-mechanical responses in nematic-elastic
rods is presented in Goriely:2022:GMM . Various light-induced mechanical
effects were reviewed in Ambulo:2020:etal ; McCracken:2020:etal ; Warner:2020
; Wen:2020:etal .

Similarly, when heat instead of light is absorbed Hiscock:2011:HWPM , it is
assumed for simplicity that the uniaxial order parameter $Q$ can be determined
by minimising the following Landau-de Gennes approximation of the nematic
energy density

$\mathcal{W}^{(lc)}(Q)=\frac{\texttt{A}}{2}Q^{2}-\frac{\texttt{B}}{3}Q^{3}+\frac{\texttt{C}}{4}Q^{4},$
(19)

where $\texttt{A},\texttt{B},\texttt{C}$ are material constants, with
$\texttt{A}=\texttt{A}(T)$ depending on temperature (Warner:2007:WT, , p. 15).
For nematic LCEs, the contribution of the above function to the total strain-
energy density including both the isotropic elastic energy and the nematic
energy functions was originally analysed in Finkelmann:2001:FGW  and more
recently in Mihai:2021:MWGG .

###  4.1 Natural deformation

When the LCE is in its natural configuration at state (A), as the capacitor is
connected to the source battery, the total strain energy of the system takes
the form (see also Hiscock:2011:HWPM )

$\mathcal{W}(\lambda_{1},\lambda)=\frac{\mu^{(1)}}{2}\left(\lambda_{1}^{2}+\lambda_{1}^{-2}\lambda^{-2}+\lambda^{2}\right)+\frac{\mu^{(2)}}{2}a^{1/3}a_{0}^{-1/3}\left(\lambda_{1}^{2}a^{-1}a_{0}+\lambda_{1}^{-2}\lambda^{-2}+\lambda^{2}\right)-\frac{CV^{2}}{2\Omega},$
(20)

where the last term represents the electrical energy per unit volume,
$\Omega=Ad$ is the volume of the dielectric, with $A$ the surface area and
$d=\lambda d_{0}$ the distance between the conductive plates, $V$ is the
voltage across the capacitor, and $C$ is the capacitance given by

$C=\frac{\varepsilon_{\perp}\varepsilon_{0\perp}A}{d}=\frac{\varepsilon_{\perp}\varepsilon_{0\perp}\Omega}{d_{0}^{2}\lambda^{2}},$
(21)

where $\varepsilon_{0\perp}$ is the permittivity for the perfectly nematic
phase and $\varepsilon_{\perp}$ is the relative permittivity when the director
is perpendicular to the electric field. The $\perp$ notation stands for the
electric field being applied perpendicular to the nematic director.

We denote by $V_{m}$ the voltage where the strain energy is comparable to the
stored elastic energy, satisfying the ratio

$V_{m}^{2}=\frac{\mu d_{0}^{2}}{\varepsilon_{0\perp}}.$ (22)

Minimising the strain energy described by equation (20) with respect to
$\lambda_{1}$ gives

$\begin{split}&\mathcal{W}(\lambda)=\frac{\mu^{(1)}}{2}\left[\left(\frac{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}\right)^{1/2}\lambda^{-1}+\left(\frac{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}\right)^{-1/2}\lambda^{-1}+\lambda^{2}\right]\\\
&+\frac{\mu^{(2)}}{2}\frac{a^{1/3}}{a_{0}^{1/3}}\left[\left(\frac{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}\right)^{1/2}\lambda^{-1}a^{-1}a_{0}+\left(\frac{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}\right)^{-1/2}\lambda^{-1}+\lambda^{2}\right]\\\
&-\frac{\mu}{2}\varepsilon_{\perp}V^{2}V_{m}^{-2}\lambda^{-2}.\end{split}$
(23)

At the initial state (A), there is no light, i.e., $I=0$, and minimising the
energy function defined by equation (18) with respect to $Q$ yields the
optimal value $Q_{0}$. At this state also, $V=V_{1}$ and $a=a_{0}$, hence the
strain energy given by equation (23) becomes

$\mathcal{W}_{n}(\lambda)=\frac{\mu}{2}\left(\lambda^{2}+2\lambda^{-1}\right)-\frac{\mu}{2}v\lambda^{-2},$
(24)

where

$v=\frac{\varepsilon_{1}V_{1}^{2}}{V_{m}^{2}}$ (25)

denotes the operating voltage.

By solving for $\lambda=\lambda_{n}$ and $v=v_{n}$ the following system of
equations,

$\frac{\partial\mathcal{W}_{n}}{\partial\lambda}=0,\qquad\frac{\partial^{2}\mathcal{W}_{n}}{\partial^{2}\lambda}=0,$
(26)

we obtain

$\lambda_{n}=\frac{1}{4^{1/3}}\qquad\mbox{and}\qquad
v_{n}=\frac{3\lambda_{n}}{4}=\frac{3}{4^{4/3}}.$ (27)

When light intensity increases so that phase transition to the isotropic state
occurs, such that $Q=0$, the order parameter is reduced and the capacitance of
the system decreases. At state (C), where $V=V_{2}$, the strain energy takes
the form

$\begin{split}&\mathcal{W}_{i}(\lambda)=\frac{\mu^{(1)}}{2}\left[\left(\frac{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}\right)^{1/2}\lambda^{-1}+\left(\frac{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}\right)^{-1/2}\lambda^{-1}+\lambda^{2}\right]\\\
&+\frac{\mu^{(2)}}{2}\frac{a^{1/3}}{a_{0}^{1/3}}\left[\left(\frac{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}\right)^{1/2}\lambda^{-1}a^{-1}a_{0}+\left(\frac{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}\right)^{-1/2}\lambda^{-1}+\lambda^{2}\right]\\\
&-\frac{\mu}{2}\varepsilon_{2}\varepsilon_{1}^{-1}V_{2}^{2}V_{1}^{-2}v\lambda^{-2}.\end{split}$
(28)

Then solving for $\lambda=\lambda_{i}$ and $v=v_{i}$ the following system of
equations,

$\frac{\partial\mathcal{W}_{i}}{\partial\lambda}=0,\qquad\frac{\partial^{2}\mathcal{W}_{i}}{\partial^{2}\lambda}=0,$
(29)

yields

$\lambda_{i}=\frac{1}{4^{1/3}}\frac{\left(\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}\right)^{1/6}}{\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)^{1/6}}$
(30)

and

$\begin{split}v_{i}&=3\left(\lambda_{i}\right)^{4}\frac{\varepsilon_{1}}{\varepsilon_{2}}\frac{V_{1}^{2}}{V_{2}^{2}}\frac{\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)}{\mu^{(1)}/\mu^{(2)}+1}\\\
&=\frac{3}{4^{4/3}}\frac{\varepsilon_{1}}{\varepsilon_{2}}\frac{V_{1}^{2}}{V_{2}^{2}}\frac{\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)^{1/3}\left(\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}\right)^{2/3}}{\mu^{(1)}/\mu^{(2)}+1}.\end{split}$
(31)

From equations (2) and (21) we find the capacitance ratio

$\xi=\frac{C_{n}}{C_{i}}=\frac{\varepsilon_{1}\lambda_{i}^{2}}{\varepsilon_{2}\lambda_{n}^{2}},$
(32)

and from equations (3), (22) and (25) we obtain

$W_{m}(\xi)=C_{n}V_{1}^{2}\frac{(\xi-1)^{2}}{4\xi}=\mu\Omega\frac{v(\xi-1)^{2}}{4\xi\lambda_{n}^{2}}$
(33)

as the optimal generated output per cycle.

###  4.2 The effect of pre-stretching perpendicular to the director

Next, we consider the LCE dielectric to be pre-stretched in the second
direction, i.e., perpendicular to the director, at state (A), with a
prescribed stretch ratio $\lambda_{0}>1$ Corbett:2009a:CW ; Corbett:2009b:CW ;
Hiscock:2011:HWPM ; Kofod:2008 . As pre-stretching increases the area of the
dielectric and reduces the distance between plates, the amount of charge that
can be taken from the battery increases. However, two types of instability may
occur in this case, namely shear striping and wrinkling Corbett:2009a:CW ;
Corbett:2009b:CW . The formation of shear stripes caused by director rotation
in elongated nematic LCEs is well understood Finkelmann:1997:FKTW ;
Higaki:2013:HTU ; Kundler:1995:KF ; Kundler:1998:KF ; Petelin:2009:PC ;
Petelin:2010:PC ; Talroze:1999:etal ; Zubarev:1999:etal  and has been modelled
extensively Conti:2002a:CdSD ; DeSimone:2000:dSD ; DeSimone:2009:dST ;
Fried:2004:FS ; Fried:2005:FS ; Fried:2006:FS ; Goriely:2021:GM ;
Kundler:1995:KF ; Mihai:2022 ; Mihai:2020a:MG ; Mihai:2021a:MG ;
Mihai:2023:MRGMG . Wrinkling in compressed LCEs was examined theoretically in
Goriely:2021:GM ; Krieger:2019:KD ; Soni:2016:SPP . In pre-stretched LCEs,
wrinkles can form if the voltage is too high, due to the so-called
electrostrictive effect observed when charging the electrodes, since the
applied Maxwell stress causes contraction in the field direction and
elongation in the perpendicular directions. We assume that the input voltage
is below but close to the critical magnitude causing wrinkling, and that any
reorientation of the nematic director that might occur is reverted (see also
Appendix A).

In this case, at the state (A), where $V=V_{1}$ and $a=a_{0}$, the strain
energy given by (23) becomes

$\mathcal{W}_{n}(\lambda,\lambda_{0})=\frac{\mu}{2}\left(\lambda^{-2}\lambda_{0}^{-2}+\lambda_{0}^{2}+\lambda^{2}\right)-\frac{\mu}{2}v\lambda^{-2}.$
(34)

Solving for $\lambda=\lambda_{n}$ the equation

$\frac{\partial\mathcal{W}_{n}}{\partial\lambda}=0$ (35)

gives

$\lambda_{n}=\left(\lambda_{0}^{-2}-v\right)^{1/4}.$ (36)

In addition, solving for $v=v_{1}$ the following system of equations,

$\frac{\partial\mathcal{W}_{n}}{\partial\lambda}=0,\qquad\frac{\partial\mathcal{W}_{n}}{\partial\lambda_{0}}=0,$
(37)

gives the wrinkling voltage

$v_{1}=\lambda_{0}^{-2}-\lambda_{0}^{-8}.$ (38)

Note that wrinkling instability occurs when the stress in the second (pre-
stretched) direction becomes zero, i.e.,
$P_{2}=\partial\mathcal{W}_{n}/\partial\lambda_{0}=0$.

At state (C), where $V=V_{2}$, the strain energy is equal to

$\begin{split}\mathcal{W}_{i}(\lambda,\lambda_{0})&=\frac{\mu^{(1)}}{2}\left(\lambda^{-2}\lambda_{0}^{-2}+\lambda_{0}^{2}+\lambda^{2}\right)+\frac{\mu^{(2)}}{2}a^{1/3}a_{0}^{-1/3}\left(\lambda^{-2}\lambda_{0}^{-2}a^{-1}a_{0}+\lambda_{0}^{2}+\lambda^{2}\right)\\\
&-\frac{\mu}{2}\varepsilon_{2}\varepsilon_{1}^{-1}V_{2}^{2}V_{1}^{-2}v\lambda^{-2}.\end{split}$
(39)

Solving for $\lambda=\lambda_{i}$ the equation

$\frac{\partial\mathcal{W}_{i}}{\partial\lambda}=0$ (40)

gives

$\lambda_{i}=\left[\frac{\lambda_{0}^{-2}\left(\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}\right)-\left(\mu^{(1)}/\mu^{(2)}+1\right)\varepsilon_{2}\varepsilon_{1}^{-1}V_{2}^{2}V_{1}^{-2}v}{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}\right]^{1/4}.$
(41)

Then solving for $v=v_{2}$ the system of equations

$\frac{\partial\mathcal{W}_{i}}{\partial\lambda}=0,\qquad\frac{\partial\mathcal{W}_{i}}{\partial\lambda_{0}}=0$
(42)

gives the wrinkling voltage

$\begin{split}v_{2}&=\frac{\varepsilon_{1}}{\varepsilon_{2}}\frac{V_{1}^{2}}{V_{2}^{2}}\left\\{\frac{\mu^{(1)}/\mu^{(2)}}{\mu^{(1)}/\mu^{(2)}+1}\left[\lambda_{0}^{-2}-\lambda_{0}^{-8}\left(\frac{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}\right)^{2}\right]\right.\\\
&+\left.\frac{a^{1/3}/a_{0}^{1/3}}{\mu^{(1)}/\mu^{(2)}+1}\left[\lambda_{0}^{-2}a_{0}/a-\lambda_{0}^{-8}\left(\frac{\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}}{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}\right)^{2}\right]\right\\}.\end{split}$
(43)

In this case also, wrinkling is attained when the stress in the second
direction is equal to zero, i.e.,
$P_{2}=\partial\mathcal{W}_{i}/\partial\lambda_{0}=0$.

From equations (2), (32), (36), and (41), we obtain the operating voltage

$v=\lambda_{0}^{-2}\frac{\xi^{2}\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)-\varepsilon_{1}^{2}\varepsilon_{2}^{-2}\left(\mu^{(1)}/\mu^{(2)}+a_{0}^{2/3}/a^{2/3}\right)}{\xi^{2}\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)-\varepsilon_{1}\varepsilon_{2}^{-1}\left(\mu^{(1)}/\mu^{(2)}+1\right)\left(\xi+1\right)^{2}/4},$
(44)

which is a nonlinear function of the pre-stretch ratio $\lambda_{0}$ and the
capacitance ratio $\xi$.

###  4.3 The effect of pre-stretching parallel to the director

We also consider the case when the LCE is pre-stretched parallel to the
nematic director at state (A), with a prescribed stretch ratio
$\lambda_{0}>1$. In this case, at the state $(A)$, where $V=V_{1}$ and
$a=a_{0}$, the strain energy given by equation (23) becomes

$\mathcal{W}_{n}(\lambda,\lambda_{0})=\frac{\mu}{2}\left(\lambda_{0}^{2}+\lambda^{-2}\lambda_{0}^{-2}+\lambda^{2}\right)-\frac{\mu}{2}v\lambda^{-2}.$
(45)

As before, solving for $\lambda=\lambda_{n}$ the equation

$\frac{\partial\mathcal{W}_{n}}{\partial\lambda}=0$ (46)

gives

$\lambda_{n}=\left(\lambda_{0}^{-2}-v\right)^{1/4}.$ (47)

Then solving for $v=v_{1}$ the following system of equations,

$\frac{\partial\mathcal{W}_{n}}{\partial\lambda}=0,\qquad\frac{\partial\mathcal{W}_{n}}{\partial\lambda_{0}}=0,$
(48)

yields the same wrinkling voltage as that given by equation (38).

At state (C), where $V=V_{2}$, the strain energy is equal to

$\begin{split}\mathcal{W}_{i}(\lambda,\lambda_{0})&=\frac{\mu^{(1)}}{2}\left(\lambda_{0}^{2}+\lambda^{-2}\lambda_{0}^{-2}+\lambda^{2}\right)+\frac{\mu^{(2)}}{2}a^{1/3}a_{0}^{-1/3}\left(\lambda_{0}^{2}a^{-1}a_{0}+\lambda^{-2}\lambda_{0}^{-2}+\lambda^{2}\right)\\\
&-\frac{\mu}{2}\varepsilon_{2}\varepsilon_{1}^{-1}V_{2}^{2}V_{1}^{-2}v\lambda^{-2}.\end{split}$
(49)

Then solving for $\lambda=\lambda_{i}$ the equation

$\frac{\partial\mathcal{W}_{i}}{\partial\lambda}=0$ (50)

gives

$\lambda_{i}=\left[\frac{\lambda_{0}^{-2}\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)-\left(\mu^{(1)}/\mu^{(2)}+1\right)\varepsilon_{2}\varepsilon_{1}^{-1}V_{2}^{2}V_{1}^{-2}v}{\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}}\right]^{1/4}.$
(51)

Since the LCE tends to contract in the pre-stretched direction while expanding
in the thickness direction, there is no wrinkling.

By equations (2), (32), (47), and (51), we obtain

$v=\lambda_{0}^{-2}\frac{\xi^{2}\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)-\varepsilon_{1}^{2}\varepsilon_{2}^{-2}\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)}{\xi^{2}\left(\mu^{(1)}/\mu^{(2)}+a^{1/3}/a_{0}^{1/3}\right)-\varepsilon_{1}\varepsilon_{2}^{-1}\left(\mu^{(1)}/\mu^{(2)}+1\right)\left(\xi+1\right)^{2}/4},$
(52)

as a nonlinear function of $\lambda_{0}$ and $\xi$.

##  5 Numerical results

In this section, we present a set of numerical results to illustrate the
performance of the theoretical model for the LCE-based charge pump. Following
Hiscock:2011:HWPM , we choose the shear modulus
$\mu=\mu^{(1)}+\mu^{(2)}=10^{6}$ Pa and initial thickness $d_{0}=50\cdot
10^{-6}$ m for the LCE, and the dielectric constants
$\varepsilon_{0\parallel}=10$ and $\varepsilon_{0\perp}=20$.

  * •

When heat is absorbed, we take the scalar uniaxial order parameters
$Q_{0}=0.5$ and $Q=0.02$, corresponding to the nematic and isotropic states,
respectively Hiscock:2011:HWPM ;

  * •

When light is absorbed, we set $\alpha=\pi/2$, $f=1/6$, $\tilde{J}=5$ and
$\mu/\mu_{n}=0.05$ Bai:2020:BB ; Corbett:2006:CW ; Corbett:2008:CW . At the
initial state where there is no light, $I=0$ and the optimal order parameter
is $Q_{0}=0.61$, while at the isotropic state, $Q=0$ (see Figure 6).

The following parameters are then calculated:
$\varepsilon_{1}=\left[\varepsilon_{0\parallel}\left(1-Q_{0}\right)+\varepsilon_{0\perp}\left(2+Q_{0}\right)\right]/3=\bar{\varepsilon}+Q_{0}\left(\varepsilon_{0\perp}-\bar{\varepsilon}\right)$
and
$\varepsilon_{2}=\bar{\varepsilon}+Q\left(\varepsilon_{0\perp}-\bar{\varepsilon}\right)\approx\bar{\varepsilon}$,
where
$\bar{\varepsilon}=\left(2\varepsilon_{0\perp}+\varepsilon_{0\parallel}\right)/3$.

The energy needed to heat the LCE from the nematic to the isotropic phase is
assumed $H_{heat}\approx 3\cdot 10^{6}$ J/m3 (Warner:2007:WT, , Section 2.3).
When light instead of heat is absorbed, the required (solar) input is
considered $H_{light}\approx 10^{7}$ J/m3 Hiscock:2011:HWPM . The efficiency
of the system absorbing heat or light in the three cases when the elastomer is
not pre-stretched and when it is a pre-stretched either perpendicular or
parallel to the nematic director, are summarised in Figure 7.

###  5.1 Energy efficiency under natural deformation

When the LCE is in its natural configuration at initial state (A), $C=C_{n}$,
$\varepsilon_{\perp}=\varepsilon_{1}$, and $\lambda=\lambda_{n}$, while at
state (C), $C=C_{i}$, $\varepsilon_{\perp}=\varepsilon_{2}$, and
$\lambda=\lambda_{i}$. The stretch ratio $\lambda_{i}$ and the corresponding
operating voltage $v_{i}$, defined by equations (30) and (31), respectively,
are plotted as functions of the parameter ratio $\mu^{(1)}/\mu^{(2)}$ in
Figure 8. By varying the parameter ratio $\mu^{(1)}/\mu^{(2)}$, the maximum
optimal output per unit volume is shown in Figure 9. For example, if
$\mu^{(1)}/\mu^{(2)}=1$, then:

  * •

When the LCE absorbs heat, the maximum optimal output is equal to
$W_{m}/\Omega\approx 3.5\cdot 10^{4}$ J/m3 per cycle. The efficiency is
$W_{m}/\left(H_{heat}\Omega\right)\approx 0.035/3\approx 1.2\%$. The operating
voltages are $V_{1}=\left[\mu
v_{i}d_{0}^{2}/\left(\varepsilon_{1}\varepsilon_{0\perp}\right)\right]^{1/2}\approx
1.8$ kV and $V_{2}=V_{1}\left(\xi+1\right)/2\approx 2.2$ kV (see also Appendix
B).

  * •

When the LCE absorbs light, the maximum optimal output is $W_{m}/\Omega\approx
0.6\cdot 10^{5}$ J/m3 per cycle and the efficiency is
$W_{m}/\left(H_{light}\Omega\right)\approx 0.6\%$. The operating voltages are
$V_{1}\approx 1.8$ kV and $V_{2}\approx 2.3$ kV.

For these two cases, Figure 7 shows that efficiency decreases as
$\mu^{(1)}/\mu^{(2)}$ increases.

(a) (b)

Figure 6: The modified Maier-Saupe mean field energy given by equation (13) as
a function of the order parameter $Q$ when there is no light ($I=0$) or the
light is: (a) polarised or (b) unpolarised, with varying intensity $I>0$ and
$\alpha=\pi/2$, $f=1/6$, $\tilde{J}=5$. As $I$ increases, there is a
transition from the nematic phase to the isotropic phase. This transition
occurs at higher values of $I$ when the light is polarised than when it is
unpolarised. When $I=0$, the minimum energy is attained for $Q_{0}=0.61$.

(a)  (b)

Figure 7: The efficiency bound as a function of the parameter ratio
$\mu^{(1)}/\mu^{(2)}>0$ for the LCE dielectric absorbing: (a) heat or (b)
light when there is no initial pre-stretch or when there is a pre-stretch
either perpendicular or parallel to the nematic director, with ratio
$\lambda_{0}=1.25$.

(a)  (b)

Figure 8: (a) The stretch ratio $\lambda_{i}$ given by equation (30) and (b)
the voltage $v_{i}$ satisfying equation (31) as functions of the parameter
ratio $\mu^{(1)}/\mu^{(2)}$ for the LCE dielectric absorbing heat or light.

(a)  (b)

Figure 9: The scaled optimal generated output $W_{m}/(\mu\Omega)$ for
$\mu^{(1)}/\mu^{(2)}\in\\{0,1,10,100,1000\\}$ when the LCE dielectric absorbs:
(a) heat or (b) light.

(a) (b)

Figure 10: Wrinkling (a) input voltage $v_{1}$ given by equation (38) as a
function of the pre-stretch ratio $\lambda_{0}$ and (b) output voltage $v_{2}$
given by equation (43) as a function of the parameter ratio
$\mu^{(1)}/\mu^{(2)}$ when the LCE dielectric is pre-stretched perpendicular
to the director with ratio $\lambda_{0}=1.25$ and absorbs heat or light. The
maximum input wrinkling voltage is attained for $\lambda_{0}=1.2424$.

(a) (b)

Figure 11: (a) The stretch ratio $\lambda_{i}$ given by equation (41) and (b)
the capacitance ratio $\xi$ satisfying equation (44) as functions of the
parameter ratio $\mu^{(1)}/\mu^{(2)}$ when the LCE dielectric is pre-stretched
perpendicular to the director with ratio $\lambda_{0}=1.25$ and absorbs heat
or light.

(a) (b)

Figure 12: (a) The stretch ratio $\lambda_{i}$ given by equation (51) and (b)
the capacitance ratio $\xi$ satisfying equation (52) as functions of the
parameter ratio $\mu^{(1)}/\mu^{(2)}$ when the LCE dielectric is pre-stretched
parallel to the director with ratio $\lambda_{0}=1.25$ and absorbs heat or
light.

###  5.2 Energy efficiency when pre-stretching perpendicular to the director

If, at state (A), the LCE dielectric is pre-stretched perpendicular to the
director, with initial stretch ratio $\lambda_{0}>1$, then the input wrinkling
voltage, given by equation (38) is independent of $\mu^{(1)}/\mu^{(2)}$ and
attains its maximum for $\lambda_{0}=1.2424$, while the output wrinkling
voltage, given by equation (43), decreases as the parameter ratio
$\mu^{(1)}/\mu^{(2)}$ increases. The input wrinkling voltage is plotted in
Figure 10(a). To maximise the operating voltage, we choose an input close to
the corresponding wrinkling voltage. For example, when $\lambda_{0}=1.25$ and
$v=0.4$, if the parameter ratio $\mu^{(1)}/\mu^{(2)}$ varies, then the output
wrinkling voltage $v_{2}$ given by equation (43) is displayed in Figure 10(b),
and the stretch ratio $\lambda_{i}$ given by equation (41) and capacitance
ratio $\xi$ satisfying equation (44) are plotted in Figure 11. In this case,
Figure 7 suggests that efficiency decreases as $\mu^{(1)}/\mu^{(2)}$
increases. In particular, if $\mu^{(1)}/\mu^{(2)}=1$, then:

  * •

When the LCE absorbs heat, the maximum optimal output is equal to
$W_{m}/\Omega 2\cdot 10^{5}$ J/m3 per cycle. The efficiency is
$W_{m}/\left(H_{heat}\Omega\right)=0.2/3\approx 6\%$. The operating voltages
are $V_{1}\approx 1.8$ kV and $V_{2}\approx 2.4$ kV.

  * •

When the LCE absorbs light, the maximum optimal output is
$W_{m}/\Omega=3.4\cdot 10^{5}$ J/m3 per cycle. The efficiency is
$W_{m}/\left(H_{light}\Omega\right)\approx 3.4\%$. The operating voltages are
$V_{1}\approx 1.8$ kV and $V_{2}\approx 2.6$ kV.

For the numerical values of the given parameters, the auxiliary results
presented in Appendix A suggest that director rotation can be ignored when the
LCE dielectric is pre-stretched perpendicular to the director.

###  5.3 Energy efficiency when pre-stretching parallel to the director

When the LCE dielectric is pre-stretched parallel to the director, the input
wrinkling voltage is the same as that in Figure 10(a). Choosing again
$\lambda_{0}=1.25$ and $v=0.4$, when the parameter ratio $\mu^{(1)}/\mu^{(2)}$
varies, the stretch ratio $\lambda_{i}$ given by equation (51) and capacitance
ratio $\xi$ satisfying equation (52) are represented in Figure 12. We note
that $\xi=C_{n}/C_{i}>1$ for $\mu^{(1)}/\mu^{(2)}>1.6$ when the LCE is heated
and for $\mu^{(1)}/\mu^{(2)}>1.5$ when the LCE is illuminated. In this case,
the efficiency shown in Figure 7 increases with $\mu^{(1)}/\mu^{(2)}$.

##  6 Conclusion

In this study, a theoretical model is developed for a charge pump with a
parallel plate capacitor where the dielectric is made of LCE material that
naturally responds to environmental changes such as light or heat.
Specifically, heating or illuminating the LCE induces a transition from a
nematic to an isotropic state. In addition, the geometry of the dielectric
changes, and the contact area with the conducting plates and the distance
between them are altered. In the charge pump electrical circuit, at the
beginning of a reversible cycle of heating and cooling or illumination and
absence of light, first the dielectric is assumed to be in a relaxed natural
state, then it is pre-stretched so that the capacitance increases by
increasing the contact area with the plates and decreasing the distance
between them. The LCE is described by a composite strain-energy function
which, when taking its constitutive parameters to their limiting values, can
be reduced to either the purely elastic neo-Hookean model or the neoclassical
model for nematic elastomers.

From the above analysis, we infer that: (i) LCE is more efficient than rubber
when used as dielectric in a parallel plate capacitor; (ii) When the
dielectric is pre-stretched perpendicular to the director at the initial state
of the proposed cycle, the capacitor becomes more effective in raising the
voltage supplied by the source battery. To make these results analytically
tractable, in the proposed model, the coupling between elastic deformation and
photo-thermal responses was neglected. This coupling can be included for a
more accurate model. Further refinements can also take into account the ratio
between photoactive and non-photoactive mesogens when light is absorbed. Other
geometries can also be considered.

The author thanks Professor Peter Palffy-Muhoray (Advanced Materials and
Liquid Crystal Institute, Kent State University, Ohio, USA) for useful and
stimulating discussions.

##  7 Appendix A. Shear stripes formation when pre-stretching

In this appendix, we determine the stretch interval where shear stripes can
form in the nematic LCE described by the elastic strain-energy function given
by equation (12). Setting the nematic director in the relaxed and stretched
configuration, respectively, as follows,

$\textbf{n}_{0}=\left[\begin{array}[]{c}1\\\ 0\\\
0\end{array}\right],\qquad\textbf{n}=\left[\begin{array}[]{c}\cos\theta\\\
\sin\theta\\\ 0\end{array}\right],$ (A.1)

where $\theta\in[0,\pi/2]$ is the angle between n and $\textbf{n}_{0}$, the
associated natural deformation tensors, given by equation (7) are,
respectively,

$\textbf{G}_{0}=\left[\begin{array}[]{ccc}a^{1/3}&0&0\\\ 0&a^{-1/6}&0\\\
0&0&a^{-1/6}\end{array}\right]$ (A.2)

and

$\textbf{G}=\left[\begin{array}[]{ccc}a^{-1/6}+\left(a^{1/3}-a^{-1/6}\right)\cos^{2}\theta&\left(a^{1/3}-a^{-1/6}\right)\sin\theta\cos\theta&0\\\
\left(a^{1/3}-a^{-1/6}\right)\sin\theta\cos\theta&a^{-1/6}+\left(a^{1/3}-a^{-1/6}\right)\sin^{2}\theta&0\\\
0&0&a^{-1/6}\end{array}\right].$ (A.3)

To demonstrate shear-striping instability, we consider the following perturbed
deformation gradient

$\textbf{F}=\left[\begin{array}[]{ccc}\lambda^{-1/2}&0&0\\\
\varepsilon&\lambda&0\\\ 0&0&\lambda^{-1/2}\end{array}\right],$ (A.4)

where $\lambda>1$ is the stretch ratio in the direction of the applied tensile
force, and $0<\varepsilon\ll 1$ is a small shear parameter. The elastic
deformation tensor $\textbf{A}=\textbf{G}^{-1}\textbf{F}\textbf{G}_{0}$ is
then equal to

$\begin{split}\textbf{A}&=\left[\begin{array}[]{ccc}\lambda^{-1/2}\left(a^{1/2}\sin^{2}\theta+\cos^{2}\theta\right)&\lambda\left(a^{-1/2}-1\right)\sin\theta\cos\theta&0\\\
\lambda^{-1/2}\left(1-a^{1/2}\right)\sin\theta\cos\theta&\lambda\left(a^{-1/2}\sin^{2}\theta+\cos^{2}\theta\right)&0\\\
0&0&\lambda^{-1/2}\end{array}\right]\\\
&+\varepsilon\left[\begin{array}[]{ccc}\left(1-a^{1/2}\right)\sin\theta\cos\theta&0&0\\\
\left(\sin^{2}\theta+a^{1/2}\cos^{2}\theta\right)&0&0\\\
0&0&0\end{array}\right].\end{split}$ (A.5)

The eigenvalues $\\{\lambda_{1}^{2},\lambda_{2}^{2},\lambda_{3}^{2}\\}$ of the
Cauchy-Green tensor $\textbf{F}\textbf{F}^{T}$ and
$\\{\alpha_{1}^{2},\alpha_{2}^{2},\alpha_{3}^{2}\\}$ of the tensor
$\textbf{A}\textbf{A}^{T}$ satisfy the following relations, respectively,

$\lambda_{1}^{2}+\lambda_{2}^{2}+\lambda_{3}^{2}=\lambda^{2}+2\lambda^{-1}+\varepsilon^{2}$
(A.6)

and

$\begin{split}\alpha_{1}^{2}+\alpha_{2}^{2}+\alpha_{3}^{2}&=\left[\lambda^{-1/2}\left(a^{1/2}\sin^{2}\theta+\cos^{2}\theta\right)+\varepsilon\left(1-a^{1/2}\right)\sin\theta\cos\theta\right]^{2}\\\
&+\left[\lambda^{-1/2}\left(1-a^{1/2}\right)\sin\theta\cos\theta+\varepsilon\left(\sin^{2}\theta+a^{1/2}\cos^{2}\theta\right)\right]^{2}\\\
&+\left[\lambda\left(a^{-1/2}-1\right)\sin\theta\cos\theta\right]^{2}+\left[\lambda\left(a^{-1/2}\sin^{2}\theta+\cos^{2}\theta\right)\right]^{2}+\lambda^{-1}.\end{split}$
(A.7)

(a)  (b)

Figure A.1: The stretch ratio interval where shear striping occurs in the LCE
dielectric pre-stretched with ratio $\lambda_{0}$, perpendicular to the
director, while the uniaxial order parameter is equal to: (a) $Q_{0}=0.5$ or
(b) $Q_{0}=0.61$.

We define the following function,

$w(\lambda,\varepsilon,\theta)=\mathcal{W}^{(el)}(\lambda_{1},\lambda_{2},\lambda_{3},\theta,a),$
(A.8)

with
$\mathcal{W}^{(el)}(\lambda_{1},\lambda_{2},\lambda_{3},\theta,a)=W^{(el)}$
described by equation (12). Differentiating the above function with respect to
$\varepsilon$ and $\theta$, respectively, gives

$\begin{split}\frac{\partial
w(\lambda,\varepsilon,\theta)}{\partial\varepsilon}&=\mu^{(1)}\varepsilon\\\
&+\mu^{(2)}\left\\{\left[\lambda^{-1/2}\left(a^{1/2}\sin^{2}\theta+\cos^{2}\theta\right)+\varepsilon\left(1-a^{1/2}\right)\sin\theta\cos\theta\right]\left(1-a^{1/2}\right)\sin\theta\cos\theta\right.\\\
&\left.+\left[\lambda^{-1/2}\left(1-a^{1/2}\right)\sin\theta\cos\theta+\varepsilon\left(\sin^{2}\theta+a^{1/2}\cos^{2}\theta\right)\right]\left(\sin^{2}\theta+a^{1/2}\cos^{2}\theta\right)\right\\}\end{split}$
(A.9)

and

$\begin{split}\frac{\partial
w(\lambda,\varepsilon,\theta)}{\partial\theta}&=\mu^{(2)}\left\\{\left(a^{1/2}-1\right)\left[2\lambda^{-1/2}\sin\theta\cos\theta+\varepsilon\left(\sin^{2}\theta-\cos^{2}\theta\right)\right]\right.\\\
&\left.\cdot\left[\lambda^{-1/2}\left(a^{1/2}\sin^{2}\theta+\cos^{2}\theta\right)+\varepsilon\left(1-a^{1/2}\right)\sin\theta\cos\theta\right]\right.\\\
&\left.+\left(a^{1/2}-1\right)\left[\lambda^{-1/2}\left(\sin^{2}\theta-\cos^{2}\theta\right)-2\varepsilon\sin\theta\cos\theta\right]\right.\\\
&\left.\cdot\left[\lambda^{-1/2}\left(1-a^{1/2}\right)\sin\theta\cos\theta+\varepsilon\left(\sin^{2}\theta+a^{1/2}\cos^{2}\theta\right)\right]\right.\\\
&\left.+\lambda^{2}\left(a^{-1}-1\right)\sin\theta\cos\theta\right\\}.\end{split}$
(A.10)

The equilibrium solution minimises the energy, and thus satisfies the
simultaneous equations

$\frac{\partial
w(\lambda,\varepsilon,\theta)}{\partial\varepsilon}=0\qquad\mbox{and}\qquad\frac{\partial
w(\lambda,\varepsilon,\theta)}{\partial\theta}=0.$ (A.11)

At $\varepsilon=0$ and $\theta=0$, the partial derivatives defined by
equations (A.9)-(A.10) are both equal to zero. Hence, this trivial solution is
always an equilibrium state. For sufficiently small values of $\varepsilon$
and $\theta$, we can write the second order approximation

$w(\lambda,\varepsilon,\theta)\approx
w(\lambda,0,0)+\frac{1}{2}\left(\varepsilon^{2}\frac{\partial^{2}w}{\partial\varepsilon^{2}}(\lambda,0,0)+2\varepsilon\theta\frac{\partial^{2}w}{\partial\varepsilon\partial\theta}(\lambda,0,0)+\theta^{2}\frac{\partial^{2}w}{\partial\theta^{2}}(\lambda,0,0)\right),$
(A.12)

where

$\displaystyle\frac{\partial^{2}w}{\partial\varepsilon^{2}}(\lambda,0,0)=\mu^{(1)}+\mu^{(2)}a,$
(A.13)
$\displaystyle\frac{\partial^{2}w}{\partial\varepsilon\partial\theta}(\lambda,0,0)=\mu^{(2)}\lambda^{-1/2}\left(1-a\right),$
(A.14)
$\displaystyle\frac{\partial^{2}w}{\partial\theta^{2}}(\lambda,0,0)=\mu^{(2)}\left(\lambda^{2}-\lambda^{-1}a\right)\left(a^{-1}-1\right).$
(A.15)

First, we find the equilibrium value $\theta_{0}$ for $\theta$ as a function
of $\varepsilon$ by solving the second equation in (A.11). By the
approximation (A.12), the respective equation takes the form

$\varepsilon\frac{\partial^{2}w}{\partial\varepsilon\partial\theta}(\lambda,0,0)+\theta\frac{\partial^{2}w}{\partial\theta^{2}}(\lambda,0,0)=0,$
(A.16)

and implies

$\theta_{0}(\varepsilon)=-\varepsilon\frac{\partial^{2}w}{\partial\varepsilon\partial\theta}(\lambda,0,0)/\frac{\partial^{2}w}{\partial\theta^{2}}(\lambda,0,0).$
(A.17)

Next, substituting $\theta=\theta_{0}(\varepsilon)$ in (A.12) gives the
following function of $\varepsilon$,

$w(\lambda,\varepsilon,\theta_{0}(\varepsilon))-w(\lambda,0,0)\approx\frac{\varepsilon^{2}}{2}\left[\frac{\partial^{2}w}{\partial\varepsilon^{2}}(\lambda,0,0)-\left(\frac{\partial^{2}w}{\partial\varepsilon\partial\theta}(\lambda,0,0)\right)^{2}/\frac{\partial^{2}w}{\partial\theta^{2}}(\lambda,0,0)\right].$
(A.18)

Depending on whether the expression on the right-hand side in (A.18) is
positive, zero, or negative, the respective equilibrium state is stable,
neutrally stable, or unstable Mihai:2020a:MG ; Mihai:2021a:MG ;
Mihai:2023:MRGMG . We deduce that the equilibrium state with $\varepsilon=0$
and $\theta=0$ is unstable if

$a^{1/3}\left(\frac{\mu^{(1)}/\mu^{(2)}+1}{\mu^{(1)}/\mu^{(2)}+a}\right)^{1/3}<\lambda<a^{1/3}.$
(A.19)

Similarly, at $\varepsilon=0$ and $\theta=\pi/2$, both the partial derivatives
defined by (A.9)-(A.10) are equal to zero, and

$\displaystyle\frac{\partial^{2}w}{\partial\varepsilon^{2}}(\lambda,0,\pi/2)=\mu^{(1)}+\mu^{(2)},$
(A.20)
$\displaystyle\frac{\partial^{2}w}{\partial\varepsilon\partial\theta}(\lambda,0,\pi/2)=\mu^{(2)}\lambda^{-1/2}\left(a-1\right),$
(A.21)
$\displaystyle\frac{\partial^{2}w}{\partial\theta^{2}}(\lambda,0,\pi/2)=\mu^{(2)}\left(\lambda^{2}-\lambda^{-1}a\right)\left(1-a^{-1}\right).$
(A.22)

Thus the equilibrium state with $\varepsilon=0$ and $\theta=\pi/2$ is unstable
if

$a^{1/3}<\lambda<a^{1/3}\left(\frac{\mu^{(1)}/\mu^{(2)}+a}{\mu^{(1)}/\mu^{(2)}+1}\right)^{1/3}.$
(A.23)

In Figure A.1, we plot the bounds given by equations (A.19) and (A.23) for the
LCE dielectric pre-stretched by ratio $\lambda=\lambda_{0}$, perpendicular to
the director, when $a=a_{0}$ is given by equation (8) with $Q=Q_{0}$ and
$Q_{0}=0.5$ or $Q_{0}=0.61$. For example, if $\lambda_{0}=1.25$, then shear
striping cannot occur for $\mu^{(1)}/\mu^{(2)}>1.9$ when the LCE is heated and
for $\mu^{(1)}/\mu^{(2)}>1.5$ when the LCE is illuminated.

##  8 Appendix B. Step-length ratios

We remark here that, for $\mu^{(1)}=0$, the LCE model defined by equation (12)
reduces to the neoclassical model considered originally in Hiscock:2011:HWPM .
Assuming that the step lengths in equation (8) there satisfy the
incompressibility constraint $l_{\parallel}\left(l_{\perp}\right)^{2}=1$, and
similarly, $l^{0}_{\parallel}\left(l^{0}_{\perp}\right)^{2}=1$ for the relaxed
state, the step-length ratios are, respectively:

$\begin{split}p_{\parallel}&=\frac{l^{0}_{\parallel}}{l_{\parallel}}=\frac{(1+2Q)^{1/3}(1-Q)^{2/3}}{(1+2Q_{A})^{1/3}(1-Q_{A})^{2/3}}\cdot\frac{1+2Q_{A}}{1+2Q}=\frac{(1+2Q_{A})^{2/3}(1-Q)^{2/3}}{(1-Q_{A})^{2/3}(1+2Q)^{2/3}},\\\
p_{\perp}&=\frac{l^{0}_{\perp}}{l_{\perp}}=\frac{(1+2Q)^{1/3}(1-Q)^{2/3}}{(1+2Q_{A})^{1/3}(1-Q_{A})^{2/3}}\cdot\frac{1-Q_{A}}{1-Q}=\frac{(1+2Q)^{1/3}(1-Q_{A})^{1/3}}{(1-Q)^{1/3}(1+2Q_{A})^{1/3}}.\end{split}$
(B.1)

In our notation, these ratios take the equivalent form

$p_{\parallel}=\frac{a_{0}^{2/3}}{a^{2/3}}=\frac{(1+2Q_{0})^{2/3}(1-Q)^{2/3}}{(1-Q_{0})^{2/3}(1+2Q)^{2/3}},\qquad
p_{\perp}=\frac{a^{1/3}}{a_{0}^{1/3}}=\frac{(1+2Q)^{1/3}(1-Q_{0})^{1/3}}{(1-Q)^{1/3}(1+2Q_{0})^{1/3}}.$
(B.2)

Therefore, in the case of natural deformation under temperature changes with
in-plane nematic director, an upper efficiency bound of $2.7\%$ is found in
this study when $\mu^{(1)}=0$ (see Figure 7(a)) compared to $1\%$ reported in
Hiscock:2011:HWPM  .

## References

  * (1) Ambulo CP, Tasmin S, Wang S, Abdelrahman MK, Zimmern PE, Ware TH. 2020. Processing advances in liquid crystal elastomers provide a path to biomedical applications, Journal of Applied Physics 128, 140901 (doi: 10.1063/5.0021143).
  * (2) Bai R, Bhattacharya K. 2020. Photomechanical coupling in photoactive nematic elastomers, Journal of the Mechanics and Physics of Solids 144, 104115 (doi: 10.1016/j.jmps.2020.104115).
  * (3) Bladon P, Terentjev EM, Warner M. 1994. Deformation-induced orientational transitions in liquid crystal elastomers, Journal de Physique II 4, 75-91 (doi: 10.1051/jp2:1994100).
  * (4) Conti S, DeSimone A, Dolzmann G. 2002. Soft elastic response of stretched sheets of nematic elastomers: a numerical study, Journal of the Mechanics and Physics of Solids 50, 1431-1451 (doi: 10.1016/S0022-5096(01)00120-X).
  * (5) Corbett D, Warner M. 2006. Nonlinear photoresponse of disordered elastomers, Physical Review Letters 96, 237802 (doi: 10.1103/PhysRevLett.96.237802).
  * (6) Corbett D, Warner M. 2008. Polarization dependence of optically driven polydomain elastomer mechanics, Physical Review E 78, 061701 (doi: 10.1103/PhysRevE.78.061701).
  * (7) Corbett D, Warner M. 2009. Anisotropic electrostatic actuation, Journal of Physics D. Applied Physics 42, 115505 (doi: 10.1088/0022-3727/42/11/115505).
  * (8) Corbett D, Warner M. 2009. Electromechanical elongation of nematic elastomers for actuation, Sensors and Actuators A 149, 120-129 (doi: 10.1016/j.sna.2008.10.006).
  * (9) de Gennes PG, Prost J. 1993. The Physics of Liquid Crystals, 2nd ed, Clarendon Press, Oxford.
  * (10) de Jeu WH (ed). 2012. Liquid Crystal Elastomers: Materials and Applications, Springer, New York.
  * (11) DeSimone A, Dolzmann G. Material instabilities in nematic elastomers. Physica D. 2000; 136(1-2): 175-191 (doi: S0167-2789(99)00153-0).
  * (12) DeSimone A, Teresi L. Elastic energies for nematic elastomers. Eur Phys J E. 2009; 29: 191-204 (doi: 10.1140/epje/i2009-10467-9).
  * (13) Ebralidze TD. 1995. Weigert hologram, Applied Optics 34(8), 1357-1362 (doi: 10.1364/AO.34.001357).
  * (14) Finkelmann H, Greve A, Warner M. 2001. The elastic anisotropy of nematic elastomers, The European Physical Journal E 5, 281-293 (doi: 10.1007/s101890170060).
  * (15) Finkelmann H, Kundler I, Terentjev EM, Warner M. Critical stripe-domain instability of nematic elastomers. J Phys II. 1997; 7: 1059-1069 (doi: 10.1051/jp2:1997171).
  * (16) Finkelmann H, Nishikawa E, Pereita GG, Warner M. 2001. A new opto-mechanical effect in solids, Physical Review Letters 87(1), 015501 (doi:10.1103/PhysRevLett.87.015501).
  * (17) Fried E, Sellers S. 2004. Free-energy density functions for nematic elastomers, Journal of the Mechanics and Physics of Solids 52(7), 1671-1689 (doi: 10.1016/j.jmps.2003.12.005).
  * (18) Fried E, Sellers S. 2005. Orientational order and finite strain in nematic elastomers, The Journal of Chemical Physics 123(4), 043521 (doi: 10.1063/1.1979479).
  * (19) Fried E, Sellers S. 2006. Soft elasticity is not necessary for striping in nematic elastomers, Journal of Applied Physics 100, 043521 (doi: 10.1063/1.2234824).
  * (20) Gallo AB, Simões-Moreira JR, Costa HKM, Santos MM. 2016. Energy storage in the energy tansition context: A technology review, Renewable and Sustainable Energy Reviews 65, 800-822 (doi: 10.1016/j.rser.2016.07.028).
  * (21) Goriely A, Mihai LA. Liquid crystal elastomers wrinkling. Nonlinearity. 2021; 34(8): 5599-5629 (doi: 10.1088/1361-6544/ac09c1).
  * (22) Goriely A, Moulton DE, Mihai LA. 2022. A rod theory for liquid crystalline elastomers, Journal of Elasticity, (doi: 10.1007/s10659-021-09875-z).
  * (23) Higaki H, Takigawa T, Urayama K. Nonuniform and uniform deformations of stretched nematic elastomers. Macromolecules. 2013; 46: 5223-5231 (doi: 10.1021/ma400771z).
  * (24) Hiscock T, Warner M, Palffy-Muhoray P. 2011. Solar to electrical conversion via liquid crystal elastomers, Journal of Applied Physics 109, 104506 (doi: 10.1063/1.3581134).
  * (25) Kakichashvili SD. 1974. Method for phase polarization recording of holograms, Soviet Journal of Quantum Electronics - IOPscience 4(6), 795-798.
  * (26) Kofod G. 2008. The static actuation of dielectric elastomer actuators: how does pre-stretch improve actuation?, Journal of Physics D. Applied Physics 41, 215405 (doi: 10.1088/0022-3727/41/21/215405).
  * (27) Krieger MS, Dias MA. 2019. Tunable wrinkling of thin nematic liquid crystal elastomer sheets, Physical Review E 100, 022701 (doi: 10.1103/PhysRevE.100.022701).
  * (28) Kundler I, Finkelmann H. Strain-induced director reorientation in nematic liquid single crystal elastomers, Macromol Rapid Commun. 1995; 16: 679-686 (doi: 10.1002/marc.1995.030160908).
  * (29) Kundler I, Finkelmann H. Director reorientation via stripe-domains in nematic elastomers: influence of cross-link density, anisotropy of the network and smectic clusters. Macromol Chem Phys. 1998; 199: 677-686.
  * (30) Lee V, Bhattacharya K. 2021. Actuation of cylindrical nematic elastomer balloons, Journal of Applied Physics 129, 114701 (doi: 10.1063/5.0041288).
  * (31) Lu T, Ma C, Wang T. 2020. Mechanics of dielectric elastomer structures: A review, Extreme Mechanics Letters 38, 100752 (doi: 10.1016/j.eml.2020.100752).
  * (32) Maier W, Saupe A. 1959. Eine einfache molekular-statistische Theorie der nematischen kristallinflüssigen Phase. Teil l1, Zeitschrift für Naturforschung A 14, 882-889 (doi: 10.1515/zna-1959-1005).
  * (33) McCracken JM, Donovan BR, White TJ. 2020. Materials as machines, Advanced Materials 32, 1906564 (doi: 10.1002/adma.201906564).
  * (34) Mihai LA. 2022. Stochastic Elasticity: A Nondeterministic Approach to the Nonlinear Field Theory, Springer Cham, Switzerland (doi: https://doi.org/10.1007/978-3-031-06692-4).
  * (35) Mihai LA, Goriely A. 2020. Likely striping in stochastic nematic elastomers, Mathematics and Mechanics of Solids 25(10), 1851-1872 (doi: 0.1177/1081286520914958).
  * (36) Mihai LA, Goriely A. 2021. Instabilities in liquid crystal elastomers, Material Research Society (MRS) Bulletin 46 (doi: 10.1557/s43577-021-00115-2).
  * (37) Mihai LA, Wang H, Guilleminot J, Goriely A. 2021. Nematic liquid crystalline elastomers are aeolotropic materials, Proceedings of the Royal Society A 477, 20210259 (doi: 10.1098/rspa.2021.0259).
  * (38) Mihai LA, Raistrick T, Gleeson HF, Mistry D, Goriely A. 2023. A predictive theoretical model for stretch-induced instabilities in liquid crystal elastomers, Liquid Crystals (doi: 10.1080/02678292.2022.2161655).
  * (39) Moretti G, Rosset S, Vertechy R, Anderson I, Fontana M. 2020. A review of dielectric elastomer generator systems, Advanced Intelligent Systems 2(10), 2000125 (doi: 10.1002/aisy.202000125).
  * (40) Petelin A, Čopič M. Observation of a soft mode of elastic instability in liquid crystal elastomers. Phys Rev Lett. 2009; 103: 077801 (doi: 10.1103/PhysRevLett.103.077801).
  * (41) Petelin A, Čopič M. Strain dependence of the nematic fluctuation relaxation in liquid-crystal elastomerss. Phys Rev E. 2010; 82: 011703 (doi: 10.1103/PhysRevE.82.011703).
  * (42) Renewables 2022 global status report $https://www.ren21.net/gsr-2022/$.
  * (43) Soni H, Pelcovits RA, Powers TR. 2016. Wrinkling of a thin film on a nematic liquid-crystal elastomer, Physical Review E 94, 012701 (doi: 10.1103/PhysRevE.94.012701).
  * (44) Talroze RV, Zubarev ER, Kuptsov SA, Merekalov AS, Yuranova TI, Plate NA, Finkelmann H. Liquid crystal acrylate-based networks: polymer backbone-LC order interaction. React Funct Polym. 1999; 41: 1-11 (doi: 10.1016/S1381-5148(99)00032-2).
  * (45) Treloar LRG. 1944. Stress-strain data for vulcanized rubber under various types of deformation, Transactions of the Faraday Society 40, 59-70 (doi: 10.1039/TF9444000059).
  * (46) Warner M. 2020. Topographic mechanics and applications of liquid crystalline solids, Annual Review of Condensed Matter Physics 11, 125-145 (doi: 10.1146/annurev-conmatphys-031119-050738).
  * (47) Warner M, Gelling KP, Vilgis TA. 1988. Theory of nematic networks, The Journal of Chemical Physics 88, 4008-4013 (doi: 10.1063/1.453852).
  * (48) Warner M, Terentjev EM. 2007. Liquid Crystal Elastomers, paper back, Oxford University Press, Oxford, UK.
  * (49) Warner M, Wang, XJ. 1991. Elasticity and phase behavior of nematic elastomers, Macromolecules 24, 4932-4941 (doi: 10.1021/ma00017a033).
  * (50) Wen Z, Yang K, Raquez JM. 2020. A review on liquid crystal polymers in free-standing reversible shape memory materials, Molecules 25, 1241 (doi:10.3390/molecules25051241).
  * (51) Yu Y, Nakano M, Ikeda T. 2003. Directed bending of a polymer film by light, Nature 425, 145-145 (doi: 10.1038/425145a).
  * (52) Zubarev ER, Kuptsov SA, Yuranova TI, Talroze RV, Finkelmann H. Monodomain liquid crystalline networks: reorientation mechanism from uniform to stripe domains. Liq Cryst. 1999; 26: 1531-1540 (doi: 10.1080/026782999203869).