# Liquid Crystals as Phase Change Materials for Thermal Stabilization

Authors: Eva Klemenčič; Mitja Slavinec
Journal: Advances in Condensed Matter Physics
 (2018)

---

## Abstract

Thermal stabilization exploiting phase change materials (PCMs) is studied theoretically and numerically. Using the heat source approach in numerical simulations, we focus on phase change temperature as a key factor in improving thermal stabilization. Our focus is to analyze possible mechanisms to tune the phase change temperature. We use thermotropic liquid crystals (LCs) as PCMs in a demonstrative system. Using the Landau-de Gennes mesoscopic approach, we show that an external electric field or appropriate nanoparticles (NPs) dispersed in LCs can be exploited to manipulate the phase change temperature.


---

## Body

## 1. Introduction

In order to achieve and maintain an optimum ambient temperature, one must regulate heat flow using different heating or cooling systems, which are usually large energy consumers. Recently, several studies [1–5] focused on the development of innovative materials to passively reduce temperature fluctuations. Phase change materials (PCMs) are suitable candidates, as they undergo a phase change during which they exchange latent heat with the surroundings.The main function of PCMs is saving excess heat that is accumulated at higher temperature conditions and releasing it at a lower temperature. Therefore, the efficiency of a system using PCMs is strongly correlated with phase change temperature. For example, if one wants to decrease room temperature oscillation, it would be optimal to have a PCM with phase change around the desired room temperature. If the room temperature exceeds the temperature of the PCM composite and reaches the phase change temperature, the PCM undergoes a phase change. During the phase change, the PCM accumulates the excess heat in form of latent heat. The opposite happens as the temperature lowers towards the phase change temperature and the latent heat is released to the surroundings.In this paper, we focus on analyzing possible mechanisms to tune the phase change temperature in order to obtain optimal thermal stabilization. We used a nematic liquid crystal (LC) system as an experimental testing ground. LCs are experimentally suitable due to their liquid character, softness, large susceptibility to relatively weak external perturbations, and optical anisotropy [6, 7]. We considered thermotropic LC in which, by decreasing the temperature from the isotropic (I) phase, the nematic (N) phase is obtained via first-order continuous symmetry breaking phase transition [8–10]. In the uniaxial bulk NLC phase, rod-like molecules exhibit long-range orientational ordering commonly presented by the nematic director n⃑. Fluctuations about n⃑ describe the uniaxial order parameter S. Note that n⃑=1 and states ±n⃑ are equivalent due to the so-called head-to-tail invariance. At the mesoscopic scale, we described nematic orientational ordering by the nematic tensor order parameterQ[11]. Our main goal was to identify different mechanisms to control and manipulate the phase change temperature. In general, the presence of nanoparticles (NPs), external electric fields, and other perturbations gives rise to elastic distortions [12–16] that can result in the temperature shift of the phase transition.First we introduce a numerical model to analyze thermal stabilization using PCMs. The model is based on the heat source approach and finite difference method. Next, we present a theoretical Landau-de Gennes model to examine the impact of the external electric field and NPs on the phase change temperature of NLC system. We present and discuss results in Section3. In the final section, we conclude our findings.


## 2. Modeling

### 2.1. Thermal Stabilization with PCM

In this section, we describe the use of PCMs for the improvement of thermal stabilization. We considered one-dimensional heat transfer through the planar surface of a two-layered composite (Figure1) with one layer integrating homogeneously distributed PCMs. Our numerical model is based on the heat source approach [17–20], in which latent heat is modeled as an additional heat source. In terms of enthalpy (h), the transient heat transfer equation in one dimension is expressed as(1)ρ∂h∂t=λ∂T2∂x2,where ρ is density and λ thermal conductivity of each layer.Figure 1
Schematic representation of two layers, one integrating PCMs.Assuming high porosity of the composite layer, the expansion of PCMs during the phase transition does not change the total volume. Additionally, we considered different physical properties of the solid(s) and liquid (l) phase of PCMs.In order to simplify their numerical simulations, other studies [17, 19, 20] determined the temperature range where phase transition occurs. In this study, we assumed that phase transition is instant at one specific critical temperature Tc. We introduced the factor r that equals zero when the PCM is in the solid phase (T<Tc) and one when the PCM is in the liquid phase (T>Tc).Using a model of effective physical quantities [18] we defined the effective specific heat capacity (ce), thermal conductivity (λe), and density (ρe) as follows:(2a)ce=αcPCMlr+cPCMs1-r+1-αcA,(2b)λe=αλPCMlr+λPCMs1-r+1-αλA,(2c)ρe=1-αρA+αρPCMlr+ρPCMs1-r,where α represents the relative amount of PCM integrated in the composite layer and subscripts A and B correspond to two different layers of studied composite (see Figure 1).We derived the heat transfer equation for the composite incorporating PCMs in terms of the enthalpyhe=ceT+αrQL as follows: (3)ρece∂T∂t+ρeαQL∂r∂t=λe∂T2∂x2,where QL is the latent heat of the integrated PCM. For numerical simulations, we used the finite difference method.


### 2.2. Phase Change Temperature Shift

Next, we theoretically analyzed possible mechanisms to tune the phase change temperature. We considered the NLC rectangular cell, where lateral confinement walls impose strong homogeneous tangential anchoring (Figure2). The characteristic confinement size R is assumed to be larger than the cell thickness h.Figure 2
NLC cell of thicknessh along the z-Cartesian coordinate with conflicting strong homogeneous tangential anchoring.Using the nematic tensor order parameter [6], we described LC orientational ordering as follows:(4)Q=∑i=13sie⃑i⊗e⃑i,where e⃑i are eigenvectors, si are eigenvalues ofQ, and the symbol ⊗ represents tensor product. In case of uniaxial ordering, it holds as follows [12]:(5)Q=Sn⃑⊗n⃑-13I,whereIis the identity matrix.Next, we expressed the free energy density as a sum of the condensation nematic termfn, elastic term fe, external field term ff, and surface term fs [12, 13]: (6a)fn=A0T-T∗2TrQ2-B3TrQ3+C4TrQ22,(6b)fe=L2∇Q2,(6c)ff=ε0Δε2E⃑·QE⃑,(6d)fs=w2TrQ-Qs,where A0, B, C are material-dependent constants and T∗ is the isotropic supercooling temperature. The elastic term evaluates deviations from a spatially homogeneous orientational ordering in terms of the elastic constantL. The external field term describes the contribution of an external electric field E⃑, where Δε stands for the electric field anisotropy constant and ε0 represents the electrical permittivity constant. Since we imposed homogeneous tangential anchoring at the lateral surfaces with the anchoring strength w,Qs describes the degree of orientational ordering enforced by the confining surface.The I-N phase change temperature is determined by material-dependent constants:(7)TIN=T∗+B24A0C.Next, assuming uniaxial ordering described by (5), we imposed a distortion to orientational ordering on a characteristic length scale R. On average, it holds as follows (see (6b)) [21, 22]:(8)fe¯~L2S2R2,where the overbar ⋯- indicates spatial averaging.Then, we expressed the average contributions of the condensation nematic term and elastic term to the average free energy densityf- as follows [23–25]:(9)f-~A0T-T∗3S¯2-B9S¯3+C9S¯4+L2S¯2R2≡A0T-T∗eff3S¯2-B9S¯3+C9S¯4,where(10)T∗eff=T∗-32LA0R2.Considering L, A0, and T∗ are material-dependent coefficients, one could vary the effective phase change temperature of the NLC system by varying the typical distortion length R. We obtained the phase change temperature shift ΔT=TIN-TIN(R) expressed as follows (see (7)) [26, 27]:(11)ΔTINR=-32LA0R2.Since the temperature shift is connected to R, we focused on imposing different distortions that can be described by a linear deformation imposed characteristic length R.Below, we present two experimentally attainable possibilities influencing the phase change temperature of NLC: the presence of an external electric field and NPs. Since the presence of NPs and the external electric field could induce biaxiality, we introduce the degree of biaxiality as follows [28]:(12)β2=1-6TrQ32TrQ23,which is zero for uniaxial states and one when the state exhibits a maximum degree of biaxiality.For further convenience, we introduce three characteristic lengths: the external field coherence lengthξf=LS/(ε0ΔεE2), the surface anchoring length de=L/w, and the bare biaxial order parameter correlation length ξb=LC/B. In addition, we introduce parametrization and scaling as described in the following [12]: (13)Qx,y=q3+q1e⃑x⊗e⃑x+q3-q1e⃑y⊗e⃑y+q2e⃑x⊗e⃑y+e⃑y⊗e⃑x-2q3e⃑z⊗e⃑z,where q1, q2, and q3 are x and y functions and e⃑3=e⃑z is an eigenvector of Q.In this scaling, we define the scaled temperature ast=(T-T∗)/(T∗∗-T∗). Here, the superheating temperature of the nematic phase is T∗∗=T∗+B2/24A0C. For numerical convenience we introduce τ=1+1-t. With this in mind, we express the equilibrium uniaxial order parameter as Seq=S∗∗τ where S∗∗=B/4C=Seq(T∗∗).Finally, we obtain bulk Euler-Lagrange equations for variation parametersq1, q2, and q3 by minimization of the free energy density [29]: (14a)Δ⊥q1ξb0R2-τ6q1+2q2q1-q123q22+q12+q32+14ξb0ξf2e⃑x·e⃑2-e⃑y·e⃑2=0,(14b)Δ⊥q2ξb0R2-τ6q2+13q12+q32-3q22-q223q22+q12+q32+112ξb0ξf2=0,(14c)Δ⊥q3ξb0h2-τ6q3+2q2q3-q323q22+q12+q32+12ξb0ξf2e⃑x·e⃑e⃑y·e⃑=0,where Δ⊥=∂2/∂x2+∂2/∂y2. Euler-Lagrange equations are solved using the standard relaxation method numerically.


## 2.1. Thermal Stabilization with PCM

In this section, we describe the use of PCMs for the improvement of thermal stabilization. We considered one-dimensional heat transfer through the planar surface of a two-layered composite (Figure1) with one layer integrating homogeneously distributed PCMs. Our numerical model is based on the heat source approach [17–20], in which latent heat is modeled as an additional heat source. In terms of enthalpy (h), the transient heat transfer equation in one dimension is expressed as(1)ρ∂h∂t=λ∂T2∂x2,where ρ is density and λ thermal conductivity of each layer.Figure 1
Schematic representation of two layers, one integrating PCMs.Assuming high porosity of the composite layer, the expansion of PCMs during the phase transition does not change the total volume. Additionally, we considered different physical properties of the solid(s) and liquid (l) phase of PCMs.In order to simplify their numerical simulations, other studies [17, 19, 20] determined the temperature range where phase transition occurs. In this study, we assumed that phase transition is instant at one specific critical temperature Tc. We introduced the factor r that equals zero when the PCM is in the solid phase (T<Tc) and one when the PCM is in the liquid phase (T>Tc).Using a model of effective physical quantities [18] we defined the effective specific heat capacity (ce), thermal conductivity (λe), and density (ρe) as follows:(2a)ce=αcPCMlr+cPCMs1-r+1-αcA,(2b)λe=αλPCMlr+λPCMs1-r+1-αλA,(2c)ρe=1-αρA+αρPCMlr+ρPCMs1-r,where α represents the relative amount of PCM integrated in the composite layer and subscripts A and B correspond to two different layers of studied composite (see Figure 1).We derived the heat transfer equation for the composite incorporating PCMs in terms of the enthalpyhe=ceT+αrQL as follows: (3)ρece∂T∂t+ρeαQL∂r∂t=λe∂T2∂x2,where QL is the latent heat of the integrated PCM. For numerical simulations, we used the finite difference method.


## 2.2. Phase Change Temperature Shift

Next, we theoretically analyzed possible mechanisms to tune the phase change temperature. We considered the NLC rectangular cell, where lateral confinement walls impose strong homogeneous tangential anchoring (Figure2). The characteristic confinement size R is assumed to be larger than the cell thickness h.Figure 2
NLC cell of thicknessh along the z-Cartesian coordinate with conflicting strong homogeneous tangential anchoring.Using the nematic tensor order parameter [6], we described LC orientational ordering as follows:(4)Q=∑i=13sie⃑i⊗e⃑i,where e⃑i are eigenvectors, si are eigenvalues ofQ, and the symbol ⊗ represents tensor product. In case of uniaxial ordering, it holds as follows [12]:(5)Q=Sn⃑⊗n⃑-13I,whereIis the identity matrix.Next, we expressed the free energy density as a sum of the condensation nematic termfn, elastic term fe, external field term ff, and surface term fs [12, 13]: (6a)fn=A0T-T∗2TrQ2-B3TrQ3+C4TrQ22,(6b)fe=L2∇Q2,(6c)ff=ε0Δε2E⃑·QE⃑,(6d)fs=w2TrQ-Qs,where A0, B, C are material-dependent constants and T∗ is the isotropic supercooling temperature. The elastic term evaluates deviations from a spatially homogeneous orientational ordering in terms of the elastic constantL. The external field term describes the contribution of an external electric field E⃑, where Δε stands for the electric field anisotropy constant and ε0 represents the electrical permittivity constant. Since we imposed homogeneous tangential anchoring at the lateral surfaces with the anchoring strength w,Qs describes the degree of orientational ordering enforced by the confining surface.The I-N phase change temperature is determined by material-dependent constants:(7)TIN=T∗+B24A0C.Next, assuming uniaxial ordering described by (5), we imposed a distortion to orientational ordering on a characteristic length scale R. On average, it holds as follows (see (6b)) [21, 22]:(8)fe¯~L2S2R2,where the overbar ⋯- indicates spatial averaging.Then, we expressed the average contributions of the condensation nematic term and elastic term to the average free energy densityf- as follows [23–25]:(9)f-~A0T-T∗3S¯2-B9S¯3+C9S¯4+L2S¯2R2≡A0T-T∗eff3S¯2-B9S¯3+C9S¯4,where(10)T∗eff=T∗-32LA0R2.Considering L, A0, and T∗ are material-dependent coefficients, one could vary the effective phase change temperature of the NLC system by varying the typical distortion length R. We obtained the phase change temperature shift ΔT=TIN-TIN(R) expressed as follows (see (7)) [26, 27]:(11)ΔTINR=-32LA0R2.Since the temperature shift is connected to R, we focused on imposing different distortions that can be described by a linear deformation imposed characteristic length R.Below, we present two experimentally attainable possibilities influencing the phase change temperature of NLC: the presence of an external electric field and NPs. Since the presence of NPs and the external electric field could induce biaxiality, we introduce the degree of biaxiality as follows [28]:(12)β2=1-6TrQ32TrQ23,which is zero for uniaxial states and one when the state exhibits a maximum degree of biaxiality.For further convenience, we introduce three characteristic lengths: the external field coherence lengthξf=LS/(ε0ΔεE2), the surface anchoring length de=L/w, and the bare biaxial order parameter correlation length ξb=LC/B. In addition, we introduce parametrization and scaling as described in the following [12]: (13)Qx,y=q3+q1e⃑x⊗e⃑x+q3-q1e⃑y⊗e⃑y+q2e⃑x⊗e⃑y+e⃑y⊗e⃑x-2q3e⃑z⊗e⃑z,where q1, q2, and q3 are x and y functions and e⃑3=e⃑z is an eigenvector of Q.In this scaling, we define the scaled temperature ast=(T-T∗)/(T∗∗-T∗). Here, the superheating temperature of the nematic phase is T∗∗=T∗+B2/24A0C. For numerical convenience we introduce τ=1+1-t. With this in mind, we express the equilibrium uniaxial order parameter as Seq=S∗∗τ where S∗∗=B/4C=Seq(T∗∗).Finally, we obtain bulk Euler-Lagrange equations for variation parametersq1, q2, and q3 by minimization of the free energy density [29]: (14a)Δ⊥q1ξb0R2-τ6q1+2q2q1-q123q22+q12+q32+14ξb0ξf2e⃑x·e⃑2-e⃑y·e⃑2=0,(14b)Δ⊥q2ξb0R2-τ6q2+13q12+q32-3q22-q223q22+q12+q32+112ξb0ξf2=0,(14c)Δ⊥q3ξb0h2-τ6q3+2q2q3-q323q22+q12+q32+12ξb0ξf2e⃑x·e⃑e⃑y·e⃑=0,where Δ⊥=∂2/∂x2+∂2/∂y2. Euler-Lagrange equations are solved using the standard relaxation method numerically.


## 3. Results and Discussion

In this section, we first analyze the numerical results of a one-dimensional heat transfer through composite wall with a focus on varying the phase change temperature. In our simulations, the composite wall (labeled GPCM) consisted of two layers: 20,0 cm thick panel at the outside and 3,0 CM thick gypsum wallboard with PCMs at the inside. The latter is also true in real building constructions, as the purpose of PCM is to stabilize the interior temperature. The outside temperature fluctuated according to typical day/night temperatures in the summertime (from 15°C to 35°C) and the room temperature (Troom) was set to a constant value, in our case to 21°C. We examined cases for three different phase change temperatures (Tc) of PCM and obtained different scenarios for the optimal thermal stabilization depending on the initial temperature (Ti) of the system.Figure3 depicts the incoming heat flux (qin) time dependency for three different cases. First, we set the initial temperature of the composite wall to the average value of the outside temperature (Ti>Troom). It is evident that the presence of PCMs affects qin as soon as the phase change temperatureTc is reached. At this point, qn is close to zero, or zero for Tc=Troom (black solid curve), until the latent heat storage capacity becomes full. By varying the phase change temperature, we show that the time period of a steady incoming heat flux is the longest when Tc is slightly above the room temperature (red dotted curve). On the other hand, in order to maintain the constant room temperature in this case, the incoming heat flux should be compensated using different cooling devices. The zero incoming heat flux is reachable only by setting the phase change temperature equal to the room temperature. Negative values of qin, for the case when Tc<Troom (blue dashed curve), correspond to a decrease of a room temperature below the desired value. Next, we analyzed cases for Ti<Troom (Figure 3(b)) and Ti=Troom (Figure 3(c)). Figures 3(b) and 3(c) confirm that the incoming heat flux can be entirely or partly stored in the form of the latent heat at the Tc. It is also evident that temperature oscillations in all three cases are reduced noticeably, especially for Tc=Troom and Tc>Troom. The numerical obtained results show that thermal stabilization could be more efficient by using a small variation of the phase change temperature. Nevertheless, it is evident that the incoming heat flux fluctuates the most for the composite with no PCMs (black dotted curve).Figure 3
The time dependency of the incoming heat flux (relative units) for three values of phase change temperatureTc=Troom (black solid curve), Tc>Troom (red dotted curve), and Tc<Troom (blue dashed curve). For comparison, the composite without PCM (black dotted curve) is included. Initial temperature is set to (a) Ti>Troom, (b) Ti<Troom, and (c) Ti=Troom.
(a) (b) (c)To better understand the efficiency of PCM composites, we numerically analyzed thermal stabilization for three different composites used in real-life applications. Figure4 shows a comparison of an alternative composite GPCM (black solid curve) with composites labeled BS1 (blue dashed curve) and BS2 (red dotted curve). BS1 and BS2 both consist of classic building materials: brick wall of thickness 50,0 cm for BS1 and 20,0 cm for BS2 and 3,0 cm thick insulation layer of Styrofoam placed at the outside of a wall. We simulated typical outside temperature fluctuations for the summertime, as described above. In this case the phase change temperature of GPCM is set to the room temperature. It is noticeable that thicker wall of BS1 lowers temperature oscillations, in comparison with BS2. This is expected due to larger thermal capacity. Furthermore, it is also evident that the alternative composite GPCM, of equal thickness as BS2, improves thermal stability even more. The incoming heat flux after 48 hours of simulations is close to 0% for GPCM and around 5% for BS1. Therefore, we can conclude that both systems BS1 and GPCM optimize the thermal stability, but by using alternative composites with PCMs construction walls can be approximately 2 times thinner.Figure 4
The time dependency of the incoming heat flux (relative units) for three real-life composites BS1 (blue dashed curve), BS2 (red dotted curve), and GPCM (black solid curve).Note that in all simulations we considered temperature oscillations between 15°C and 35°C for total time of 48 hours. In order to obtain more realistic results, additional simulations using different temperature ranges and longer simulation time are welcomed. According to our numerical outcomes (Figures3 and 4), we assume similar thermal stabilization, with the main difference in the direction of the heat flow, when simulating low temperatures. In this case, we expect that the lowest temperature oscillations would be obtained for Tc=Troom and Tc<Troom. Since even small variation in Tc affects the thermal stability, it is reasonable to explore possible mechanisms to develop PCM composites with tuneable phase change temperature.We analyzed the impact of external electric fields and NPs on the phase change temperatureTc of NLC. We considered square-shaped NP that enforces strong tangential boundary conditions at its surface in the presence and absence of an external electric field. In our study, we set the lateral confinement size R below 1 micron, comparable to ξb. Figure 5 demonstrates typical nematic configurations in the absence and presence of an external electric field. Figure 5(a) depicts the diagonal structure, which is stable for cases when R≫ξb and E=0. At the corners of the cell, opposing conditions give rise to defects. In the center of each defect, a negative uniaxial order along e⃑z is established and is surrounded by a rim of a maximal degree of biaxiality. By increasing the external electric field strength, we obtained qualitatively different configurations for R/ξf=10 (Figure 5(b)) and for R/ξf=100 (Figure 5(c)). In the latter case, the external electric field triggered a surface order-reconstruction type of structural transition. Therefore, by application of strong enough external field, one obtains the qualitatively different configuration of NLC, which would affect the phase transition temperature according to (11).Figure 5
The degree of biaxialityβ2 in the absence of NPs for R/ξb=7. Uniaxial state is presented in black; a maximum degree of biaxiality is presented in white. The diagonal structure is formed (a) in the absence of an external electric field, (b) for R/ξf=10, and (c) for R/ξf=100.
(a) (b) (c)Next, we studied the impact of square-shaped NP in the absence of an external electric field (Figure6). We placed NP at four different positions of the NLC cell: (i) at the center, (ii) at the bottom boundary, (iii) at the left boundary, and (iv) at the left bottom edge. The NP acts as a source of elastic distortions that result in locally induced biaxiality.Figure 6
Absence of external electric field,R/ξb=7. Qualitatively different configurations are obtained for NP placed at the (a) center, (b) bottom boundary, (c) left boundary, and (d) left bottom edge. The color scheme is the same as in Figure 5.
(a) (b) (c) (d)We obtained an equal configuration (but rotated) for two positions of NP: close to the left boundary and close to the bottom boundary (Figures6(b) and 6(c)). This shows that the presence of NPs in an NLC cell affects the typical distortion lengthR. Regarding (10) and (11), we conclude that NPs affect effective phase change temperatures. Note that in general several NPs could be introduced which increases complexity and richness of phenomena.Additionally, we studied the combined impact of NPs and external electric fields. The external electric field broke the symmetry of the system (see Figure7) and enabled an order-reconstruction type transition at the bottom plate.Figure 7
Configurations for combined impact of the external electric fieldR/ξf=50 and four different positions of NPs (a) at the center, (b) bottom boundary, (c) left boundary, and (d) left bottom edge, for R/ξb=7. The color scheme is the same as in Figure 5.
(a) (b) (c) (d)


## 4. Conclusions

In this study, we numerically assessed the impact of PCMs on thermal stabilization. Numerical simulations based on the heat source method confirm that PCMs integrated in the composite material reduced temperature oscillations and therefore improved thermal stabilization. We focused on different cases by varying the phase change temperature around room temperature and showed that the efficiency of the thermal stabilization depends on both the phase change temperature and the desired room temperature. Therefore, it is reasonable to develop tuneable PCM composites. Our main goal was to find and analyze possible mechanisms to control and manipulate phase change temperature. We considered NLC cells as PCMs, as they are commonly used in theoretical and experimental testing.In our theoretical section, we showed that phase change temperature can be shifted by varying the typical distortion lengthR. We then examined the impact of external electric fields and NPs on the NLC configuration and consequently on the phase change temperature. We demonstrated that we could shift the phase change temperature by changing the external electric field in the absence of NPs. When an external electric field is absent, an NP can also effectively change the typical characteristic length of the NLC system.LCs can be used as PCMs to improve thermal stabilization. Furthermore, one can even tune the phase change temperature with relatively simple mechanisms. One of the main disadvantages of using LCs for thermal stabilization is relatively high cost, corresponding to the relatively small amount of latent heat. Nevertheless, LCs have potential for future applications as PCMs, especially in the space industry.



---

*Source: 1878232-2018-05-15.xml*