# Turbulent Prandtl number from isotropically forced turbulence

Petri J. Käpylä1<EMAIL_ADDRESS>Nishant K. Singh2  1 Institute for
Astrophysics and Geophysics, Göttingen University, Friedrich-Hund-Platz 1,
37077 Göttingen, Germany 2 Inter-University Centre for Astronomy and
Astrophysics, Post Bag 4 Ganeshkhind, Savitribai Phule Pune University Campus,
Pune 411 007, India

###### Abstract

Turbulent motions enhance the diffusion of large-scale flows and temperature
gradients. Such diffusion is often parameterized by coefficients of turbulent
viscosity ($\nu_{\rm t}$) and turbulent thermal diffusivity ($\chi_{\rm t}$)
that are analogous to their microscopic counterparts. We compute the turbulent
diffusion coefficients by imposing large-scale velocity and temperature
gradients on a turbulent flow and measuring the response of the system. We
also confirm our results using experiments where the imposed gradients are
allowed to decay. To achieve this, we use weakly compressible three-
dimensional hydrodynamic simulations of isotropically forced homogeneous
turbulence. We find that the turbulent viscosity and thermal diffusion, as
well as their ratio the turbulent Prandtl number, $\mbox{\rm Pr}_{\rm
t}=\nu_{\rm t}/\chi_{\rm t}$, approach asymptotic values at sufficiently high
Reynolds and Peclét numbers. We also do not find a significant dependence of
$\mbox{\rm Pr}_{\rm t}$ on the microscopic Prandtl number $\mbox{\rm
Pr}=\nu/\chi$. These findings are in stark contrast to results from the
$k-\epsilon$ model which suggests that $\mbox{\rm Pr}_{\rm t}$ increases
monotonically with decreasing Pr. The current results are relevant for the
ongoing debate of, for example, the nature of the turbulent flows in the very
low Pr regimes of stellar convection zones.

###### keywords:

Turbulence

##  1 Introduction

The fluids in stellar convection zones are generally characterized by a low
microscopic Prandtl number, $\mbox{\rm Pr}=\nu/\chi$, where $\nu$ is the
kinematic viscosity and $\chi$ is the thermal diffusivity (e.g. Ossendrijver,
2003; Augustson et al., 2019). Typical values in the bulk of the solar
convection zone, for example, range between $10^{-6}$ and $10^{-3}$
(Schumacher & Sreenivasan, 2020). Recently, several studies have explored the
possibility that solar convection operates at a high effective Prandtl number
regime, meaning that the turbulent Prandtl number $\mbox{\rm Pr}_{\rm t}$
exceeds unity (e.g. O’Mara et al., 2016; Bekki et al., 2017; Karak et al.,
2018), as a possible solution to the too high velocity amplitudes in
simulations in comparison to the Sun (e.g. Hanasoge et al., 2012; Schumacher &
Sreenivasan, 2020). However, few attempts have been made to actually measure
the turbulent Prandtl number from simulations. A notable exception is the
study of Pandey et al. (2021) who reported that the turbulent Prandtl number
decreases steeply as a function of the molecular Prandtl number such that
$\mbox{\rm Pr}_{\rm t}\propto\mbox{\rm Pr}^{-1}$ in simulations of standard
Boussinesq and variable heat conductivity Boussinesq convection.

Here we set out to measure turbulent viscosity and thermal diffusivity from a
simpler system of isotropically forced homogeneous turbulence. This is done by
imposing large-scale gradients of velocity and temperature (equivalently
specific entropy) and measuring the response of the system. The turbulent
diffusion coefficients are computed from the Boussinesq ansatz and an
analogous expression for the enthalpy flux. This method provides a direct
measurement of the diffusion coefficients without the need to resort to
turbulent closures. Similar methods were used recently to measure the
turbulent magnetic Prandtl number (Käpylä et al., 2020). We compare our
results with those from the widely used $k-\epsilon$ model which was also used
by Pandey et al. (2021). We show that the direct results and those from the
$k-\epsilon$ model are systematically different and that the latter yields
misleading results.

##  2 The model

We model isotropically forced, non-isothermal, turbulence in a fully periodic
cube of volume $(2\pi)^{3}$. We solve the equations of fully compressible
hydrodynamics

$\displaystyle\frac{D\ln\rho}{Dt}$ $\displaystyle=$
$\displaystyle-\bm{\nabla}\bm{\cdot}{\bm{u}},$ (1)
$\displaystyle\frac{D{\bm{u}}}{Dt}$ $\displaystyle=$
$\displaystyle-\frac{1}{\rho}(\bm{\nabla}p-\bm{\nabla}\bm{\cdot}2\nu\rho\bm{\mathsf{S}})+{\bm{f}}-\frac{1}{\tau}({\bm{u}}-\overline{\mbox{\boldmath$u$}}{}_{0}),$
(2) $\displaystyle T\frac{Ds}{Dt}$ $\displaystyle=$
$\displaystyle-\frac{1}{\rho}\left(\bm{\nabla}\bm{\cdot}\mbox{\boldmath${\cal
F}$}{}_{\rm rad}-{\cal
C}\right)+2\nu\bm{\mathsf{S}}^{2}-\frac{T}{\tau}(s-\overline{s}_{0}),$ (3)

where $D/Dt=\partial/\partial t+{\bm{u}}\cdot\bm{\nabla}$ is the advective
derivative, $\rho$ is the density, ${\bm{u}}$ is the velocity, $p$ is the
pressure, $\nu$ is the kinematic viscosity, $\bm{\mathsf{S}}$ is the traceless
rate-of-strain tensor with

$\displaystyle\mathsf{S}_{ij}={\textstyle{1\over
2}}(u_{i,j}+u_{j,i})-{\textstyle{1\over
3}}\delta_{ij}\bm{\nabla}\bm{\cdot}{\bm{u}},$ (4)

${\bm{f}}$ is the external forcing, $\tau$ is a relaxation timescale, and
$\overline{\mbox{\boldmath$u$}}{}_{0}$ is the target mean velocity profile.
Furthermore, $T$ is the temperature, $s$ is the specific entropy,
$\mbox{\boldmath${\cal F}$}{}_{\rm rad}$ is the radiative flux, ${\cal C}$ is
a cooling term, and $\overline{s}_{0}$ is the target mean specific entropy
profile. Radiation is modeled via the diffusion approximation, with the
radiative flux given by

$\displaystyle\mbox{\boldmath${\cal F}$}{}_{\rm
rad}=-c_{P}\rho\chi\bm{\nabla}T,$ (5)

where $c_{P}$ is the specific heat in constant pressure and $\chi$ is the
thermal diffusivity. The ideal gas equation of state $p=(c_{P}-c_{V})\rho
T={\cal R}\rho T$ is assumed, where ${\cal R}$ is the gas constant, and
$c_{V}$ is the specific heat capacity at constant volume. In the presence of
an imposed large-scale flow, viscous dissipation of kinetic energy acts as a
source for thermal energy and leads to a linear increase of the temperature.
Additional volumetric cooling is applied to counter this with

$\displaystyle{\cal C}({\bm{x}})=\rho c_{P}\frac{T({\bm{x}})-\langle
T_{0}\rangle}{\tau_{\rm cool}},$ (6)

where $\langle T_{0}\rangle$ is the volume-averaged initial temperature and
$\tau_{\rm cool}$ is a cooling timescale. We use $\tau=\tau_{\rm
cool}=(c_{s0}k_{1})^{-1}$, where $c_{s0}$ is the initial uniform value of the
sound speed and $k_{1}$ is the wavenumber corresponding to the box scale.

The external forcing is given by (see Brandenburg, 2001)

$\displaystyle{\bm{f}}=\Real\\{N(t){\bm{f}}_{{\bm{k}}(t)}\exp[{\rm
i}{\bm{k}}(t)\bm{\cdot}{\bm{x}}-{\rm i}\phi(t)]\\},$ (7)

where ${\bm{x}}$ is the position vector and $N(t)=f_{0}c_{\rm s}(k(t)c_{\rm
s}/\delta t)^{1/2}$ is a normalization factor where $f_{0}$ is the forcing
amplitude, $k=|{\bm{k}}|$, $\delta t$ is the length of the time step, and
$-\pi<\phi(t)<\pi$ is a random delta-correlated phase. The vector
${\bm{f}}_{\bm{k}}$ describes non-helical transverse waves, and is given by

$\displaystyle{\bm{f}}_{\bm{k}}=\frac{{\bm{k}}\times\hat{\bm{e}}}{\sqrt{{\bm{k}}^{2}-({\bm{k}}\bm{\cdot}\hat{\bm{e}})^{2}}},$
(8)

where $\hat{\bm{e}}$ is an arbitrary unit vector, and where the wavenumber
${\bm{k}}$ is randomly chosen. The target profiles of mean velocity and
specific entropy are given by

$\displaystyle\overline{\mbox{\boldmath$u$}}{}_{0}$ $\displaystyle=$
$\displaystyle u_{0}\sin(k_{1}z){\hat{\mbox{\boldmath$e$}}{}}_{y},$ (9)
$\displaystyle\overline{s}_{0}$ $\displaystyle=$ $\displaystyle
s_{0}\sin(k_{1}z).$ (10)

In addition to the physical diffusion, the advective terms in (1) to (3) are
implemented in terms of fifth-order upwinding derivatives with sixth-order
hyperdiffusive corrections and flow-dependent diffusion coefficients; see
Appendix B of Dobler et al. (2006).

The Pencil Code (Pencil Code Collaboration et al.,
2021)111http://github.com/pencil-code, which uses high-order finite
differences for spatial and temporal discretisation, was used to produce the
numerical simulations.

###  2.1 Units, system parameters, and diagnostics

The equations are non-dimensionalized by choosing the units

$\displaystyle[x]=k_{1}^{-1},\ \ [\rho]=\rho_{0},\ \ [u]=c_{\rm s0},\ \
[s]=c_{P},$ (11)

where $\rho_{0}$ is the initial uniform density and $c_{\rm
s0}=\sqrt{\gamma{\cal R}T_{0}}$ is the sound speed corresponding to the
initial temperature $T_{0}$. The level of velocity fluctuations is determined
by the forcing amplitude $f_{0}$ along with the kinematic viscosity. A key
system parameter is the ratio of kinematic viscosity and thermal diffusion or
the Prandtl number

$\displaystyle\mbox{\rm Pr}=\frac{\nu}{\chi},$ (12)

which is varied between $0.01$ and $10$ in the present study. The Reynolds and
Péclet numbers quantify the level of turbulence of the flows:

$\displaystyle\mbox{\rm Re}=\frac{u_{\rm rms}}{\nu k_{\rm f}},\ \ \ \mbox{\rm
Pe}=\mbox{\rm Pr}\mbox{\rm Re}=\frac{u_{\rm rms}}{\chi k_{\rm f}},$ (13)

where $u_{\rm rms}=\sqrt{\langle({\bm{u}}-\overline{\bm{u}}_{0})^{2}\rangle}$
is the volume-averaged fluctuating rms-velocity and $k_{\rm f}$ is the average
forcing wavenumber characterizing the energy injection scale. The latter is
chosen from a uniformly distributed narrow range in the vicinity of $5k_{1}$.
The imposed gradients of large-scale flow and entropy are quantified by

$\displaystyle\mbox{\rm Ma}_{s}=\frac{u_{0}}{c_{\rm s0}},\ \ \mbox{\rm
Ma}_{g}=\frac{[(\gamma-1)s_{0}T_{0}]^{1/2}}{c_{\rm s0}},$ (14)

where $\mbox{\rm Ma}_{s}$ is the Mach number of the mean flow and
$k_{U}=k_{s}=k_{1}$. The Mach number of the turbulent flow is given by

$\displaystyle\mbox{\rm Ma}=\frac{u_{\rm rms}}{\langle c_{\rm s}\rangle},$
(15)

where $\langle c_{\rm s}\rangle$ is the volume-averaged speed of sound.

Mean values are taken to be horizontal averages denoted by overbars, that is

$\displaystyle\overline{f}=\frac{1}{(2\pi)^{2}}\int_{x}\int_{y}f({\bm{x}})dxdy.$
(16)

Often an additional time average over the statistically steady part of the
simulation is taken. Volume averages are denoted by angle brackets
$\langle.\rangle$ apart from the rms-values which are always assumed to be
volume-averaged unless otherwise stated. Errors were estimated by dividing the
time series in three parts and averaging over each subinterval. The greatest
deviation from the average computed over the whole time series was taken as
the error estimate.

Set | Pr | Re | Pe | Ma | $\mbox{\rm Ma}_{s}$ | $\mbox{\rm Ma}_{g}$ | # runs
---|---|---|---|---|---|---|---
i001 | $0.01$ | $27\ldots 779$ | $0.27\ldots 7.8$ | $0.069\ldots 0.079$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $10$
i002 | $0.02$ | $27\ldots 779$ | $0.54\ldots 16$ | $0.069\ldots 0.080$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $12$
i005 | $0.05$ | $27\ldots 780$ | $1.4\ldots 39$ | $0.069\ldots 0.080$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $12$
i010 | $0.1$ | $12\ldots 800$ | $1.2\ldots 80$ | $0.060\ldots 0.082$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $30$
i020 | $0.2$ | $27\ldots 781$ | $5.4\ldots 156$ | $0.069\ldots 0.080$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $10$
i025 | $0.25$ | $397$ | $99$ | $0.081$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $2$
i050 | $0.5$ | $27\ldots 781$ | $14\ldots 390$ | $0.069\ldots 0.081$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $12$
i075 | $0.75$ | $399$ | $301$ | $0.081$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $2$
i100 | $1.0$ | $27\ldots 1582$ | $27\ldots 1582$ | $0.069\ldots 0.081$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $14$
i200 | $2.0$ | $22$ | $44$ | $0.068$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $2$
i500 | $5.0$ | $22.2$ | $111$ | $0.068$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $2$
i1000 | $10.0$ | $22.3$ | $223$ | $0.068$ | $0.01\ldots 0.03$ | $0.07\ldots 0.1$ | $2$
du001 | $0.01$ | $153$ | $1.5$ | $0.078$ | $0.01$ | $-$ | $10$
du010 | $0.1$ | $153\ldots 391$ | $15\ldots 39$ | $0.078$ | $0.01\ldots 0.03$ | $-$ | $11$
du100a | $1.0$ | $154$ | $154$ | $0.079$ | $0.01$ | $-$ | $10$
du100b | $1.0$ | $154$ | $154$ | $0.079$ | $0.03$ | $-$ | $10$
ds001 | $0.01$ | $153$ | $1.5$ | $0.078$ | $-$ | $0.1$ | $10$
ds010 | $0.1$ | $153\ldots 391$ | $15\ldots 39$ | $0.078$ | $-$ | $0.07\ldots 0.1$ | $11$
ds100 | $1.0$ | $153$ | $153$ | $0.078$ | $-$ | $0.1$ | $10$
Table 1: Summary of runs. Runs with imposed velocity or specific entropy
gradients are denoted with prefix i, whereas decay experiments of velocity
(specific entropy) are identified by prefix du (ds). Grid resolutions range
between $144^{3}$ and $1152^{3}$.

##  3 Results

The simulations discussed in the present study are listed in table 1.

###  3.1 Turbulent viscosity and heat diffusion from imposed flow and entropy
methods

We measure the turbulent viscosity and thermal diffusivity in two ways. First,
we impose sinusoidal large-scale profiles of velocity (9) or entropy (10). The
response of the system are non-zero Reynolds stress and vertical enthalpy flux
profiles that are parameterized with gradient diffusion terms (e.g. Rüdiger,
1989)

$\displaystyle{F^{\rm enth}_{z}}(z)=c_{P}\overline{(\rho
u_{z})^{\prime}T^{\prime}}\approx
c_{P}\overline{\rho}\overline{u_{z}^{\prime}T^{\prime}}=-\chi_{\rm
t}\overline{\rho}\overline{T}\frac{\partial\overline{s}}{\partial z},$ (17)

and

$\displaystyle{R_{yz}}(z)=\overline{u_{y}^{\prime}u_{z}^{\prime}}=-\nu_{\rm
t}\frac{\partial\overline{u}_{y}}{\partial z},$ (18)

where primes denote fluctuations from the mean, e.g.,
${\bm{u}}^{\prime}={\bm{u}}-\overline{\bm{u}}$. The Mach number in the current
simulations is always less than 0.1. Therefore we neglect density-dependent
terms in our analysis because they scale with $\mbox{\rm Ma}^{2}$.

Figure 1: Normalized turbulent viscosity $\tilde{\nu_{\rm t}}=\nu_{\rm
t}/\nu_{\rm t0}$ (squares) and heat diffusivity $\tilde{\chi_{\rm
t}}=\chi_{\rm t}/\chi_{\rm t0}$ (circles) as functions of Reynolds and Péclet
numbers. The crosses ($\times$) and pluses ($+$) indicate results from decay
experiments. The colours of the symbols indicate the microscopic Prandtl
number as shown by the colourbar. The dotted horizontal lines show fit to the
data for $\mbox{\rm Pe},\mbox{\rm Re}>10$ and a line proportional to Pe is
shown for low Pe.

The coefficients $\chi_{\rm t}$ and $\nu_{\rm t}$ are assumed to be scalars
and were obtained from linear fits between time-averaged ${F^{\rm enth}_{z}}$
and $\overline{\rho}\overline{T}\partial_{z}\overline{s}$ and between
${R_{yz}}$ and $\partial_{z}\overline{u}_{y}$, respectively. Results from our
simulations are shown in figure 1. We normalize $\nu_{\rm t}$ and $\chi_{\rm
t}$ by

$\displaystyle\nu_{\rm t0}=\chi_{\rm t0}={\textstyle{1\over 3}}u_{\rm
rms}k_{\rm f}^{-1},$ (19)

which is an order of magnitude estimate for the turbulent diffusion
coefficients. We note that in the parameter regimes studied here, the
estimates $\nu_{\rm t0}$ and $\chi_{\rm t0}$ are very similar in all of our
runs. Our results show that for low Péclet numbers the turbulent heat
diffusion increases in proportion to Pe for $\mbox{\rm Pe}\lesssim 1$. This is
consistent with earlier numerical results for turbulent viscosity (e.g. Käpylä
et al., 2020), magnetic diffusivity (e.g. Sur et al., 2008), and passive
scalar diffusion (Brandenburg et al., 2009), and with corresponding analytic
results in the diffusion dominated ($\mbox{\rm Pe}\ll 1$) regime. For
sufficiently large Pe, $\tilde{\chi_{\rm t}}$ tends to a constant value. The
turbulent viscosity is also roughly constant in the parameter space covered
here. For low fluid Reynolds numbers $\nu_{\rm t}$ is proportional to Re as
has been shown in Käpylä et al. (2020). However, we do not cover this
parameter regime with the current simulations.

###  3.2 Turbulent viscosity and heat diffusion from decay experiments

Figure 2: Panels (a) and (b): $\overline{u}_{y}(t,z)$ and $\overline{s}(t,z)$
normalized by $c_{s}$ and $c_{P}$, respectively, from decay experiments with
$\mbox{\rm Pr}=1$ and $\mbox{\rm Re}=157$. Red vertical lines denote end times
of exponential fits. Panels (c) and (d): temporal decays of $k_{z}/k_{1}=1$
mode of $\overline{u}_{y}$ and $\overline{s}$, respectively; black line shows
the progenitor run, and the red/gray lines indicate the decaying runs; red
part is used to fit exponential decay; blue dotted lines show the exponential
fit.

Decay experiments were made as an independent way to measure the turbulent
viscosity and heat diffusion. Snapshots from the imposed velocity/entropy
gradient runs were used as initial conditions and the relaxation terms of the
rhs of the Navier–Stokes and entropy equations were deactivated. The large-
scale velocity and entropy profiles in such runs decay due to the combined
effect of molecular and turbulent diffusion. To measure the decay rate we
monitored the amplitude of the $k=k_{1}$ components of $\overline{u}_{y}$ and
$\overline{s}$. Exponential decay laws

$\displaystyle\overline{u}_{y}(t,k_{1})=\overline{u}_{y}(t_{0},k_{1})e^{-(\nu_{\rm
t}+\nu)t},\ \ \overline{s}(t,k_{1})=\overline{s}(t_{0},k_{1})e^{-(\chi_{\rm
t}+\chi)t},$ (20)

were then fitted to the numerical data. Representative examples from decay
experiments of large-scale velocity and entropy are shown in figure 2. The
upper panels (a) and (b) show $\overline{u}_{y}(z,t)$ and
$\overline{s}_{y}(z,t)$ from typical decay experiments. The $k_{1}$ components
of these fields decay exponentially when the forcing is turned off; see panels
(c) and (d) of figure 2. Ultimately the amplitude of the $k_{1}$ mode
decreases sufficiently such that it cannot be distinguished from the
background turbulence. The time it takes to reach this state varies and
depends on the initial amplitudes $u_{0}$ and $s_{0}$. However, at the same
time these amplitudes need to be kept as low as possible to avoid non-linear
effects becoming important (see, e.g. Käpylä et al., 2020). This is
particularly important for the velocity field due to which the range from
which turbulent viscosity can be estimated is limited which necessitates
running several experiments with different snapshots as initial conditions to
reach converged values for $\nu_{\rm t}$ and $\chi_{\rm t}$.

Due to this, only a limited subset of the parameter range covered by the
imposed cases were repeated with decay experiments. We used ten snapshots from
each run for the decay experiments. The separation between the snapshots is
roughly $\Delta t=40u_{\rm rms}k_{\rm f}$ such that the realizations can be
considered uncorrelated. Results from the decay experiments are shown in
figure 1 with crosses ($\chi_{\rm t}$) and pluses ($\nu_{\rm t}$). We find
that the results from the decay experiments are consistent with those from the
imposed flow and entropy gradient methods.

###  3.3 The $k-\epsilon$ model

Figure 3: Turbulent Prandtl number $\mbox{\rm Pr}_{\rm t}^{(k-\epsilon)}$
according to (23) as a function of Péclet number. The colour of the symbols
denotes the molecular Prandtl number as indicated by the colour bar. Inset:
$\mbox{\rm Pr}_{\rm t}^{(k-\epsilon)}$ versus Pr from runs with $\mbox{\rm
Pe}>20$.

To facilitate a comparison with Pandey et al. (2021) we use the expressions of
$\nu_{\rm t}$ and $\chi_{\rm t}$ derived under the $k-\epsilon$ model with

$\displaystyle\nu_{\rm t}^{(k-\epsilon)}$ $\displaystyle=$ $\displaystyle
c^{\prime}_{\nu}k_{u}^{2}/\epsilon_{K},$ (21) $\displaystyle\chi_{\rm
t}^{(k-\epsilon)}$ $\displaystyle=$ $\displaystyle
c^{\prime}_{\chi}k_{u}k_{T}/\epsilon_{T},$ (22)

where $k_{u}=\langle{\bm{u}}^{\prime 2}\rangle/2$ is the turbulent kinetic
energy, $k_{T}=\langle T^{\prime 2}\rangle$ is the variance of the temperature
fluctuations, $c^{\prime}_{\nu}$ and $c^{\prime}_{\chi}$ are assumed to be
universal constants222Primes indicate the constancy of $c^{\prime}_{\nu}$ and
$c^{\prime}_{\chi}$, but see section 3.4 where this assumption is lifted..
Viscous and thermal dissipation rates are defined as
$\epsilon_{K}=2\nu[\langle\bm{\mathsf{S}}^{2}\rangle-\langle\bm{\mathsf{S}}_{0}\rangle^{2}]$
and
$\epsilon_{T}=\chi\langle(\bm{\nabla}T^{\prime})^{2}\rangle=\chi[\langle(\bm{\nabla}T)^{2}\rangle-\langle(\bm{\nabla}\overline{T})^{2}\rangle]$,
respectively, where we have removed contributions from the mean flow and the
mean entropy; $\bm{\mathsf{S}}_{0}$ denotes the traceless rate-of-strain
tensor as defined in (4) but with $\overline{\bm{u}}_{0}$ instead of
${\bm{u}}$. Pandey et al. (2021) computed $\mbox{\rm Pr}_{\rm t}$ using the
$k-\epsilon$ model by fixing the ratio of
$c^{\prime}_{\nu}/c^{\prime}_{\chi}$, which yields

$\displaystyle\mbox{\rm Pr}_{\rm t}^{(k-\epsilon)}=\frac{\nu_{\rm
t}^{(k-\epsilon)}}{\chi_{\rm
t}^{(k-\epsilon)}}=\frac{c^{\prime}_{\nu}}{c^{\prime}_{\chi}}\frac{k_{u}\epsilon_{T}}{k_{T}\epsilon_{u}}.$
(23)

For simplicity, we assume $c^{\prime}_{\nu}/c^{\prime}_{\chi}=1$ in this
subsection. The results are shown in figure 3. We find that taking the ratio
$c^{\prime}_{\nu}/c^{\prime}_{\chi}$ to be a constant leads to results where
$\mbox{\rm Pr}_{\rm t}^{(k-\epsilon)}$ increases monotonically with decreasing
Pr when Pe is larger than about 20; see the inset in figure 3 which reveals a
dependence of $\mbox{\rm Pr}^{-0.25}$. Qualitatively, this result is in
agreement with the one in Pandey et al. (2021). However, we would like to note
here that a strong assumption was made to reach this conclusion, namely, that
the ratio $c^{\prime}_{\nu}/c^{\prime}_{\chi}$ is fixed, and that
$c^{\prime}_{\nu}$ and $c^{\prime}_{\chi}$ are universal constants,
independent of control parameters such as Re and Pe. Henceforth, we relax
these assumptions, and also omit primes from the coefficients $c_{\nu}$ and
$c_{\chi}$.

###  3.4 Relaxing the assumption that $c_{\nu}$ and $c_{\chi}$ are constants

It is reasonable to assume that for sufficiently large Reynolds and Péclet
numbers $k_{u}$ and $k_{T}$ tend to constant values. Furthermore, there is
evidence from numerical simulations that $\epsilon_{K}$ also tends to a non-
zero constant value for large Reynolds numbers. Similar evidence for
$\epsilon_{T}$ has not been presented. Therefore it is not clear whether the
assumption of universality of $c_{\nu}$ and $c_{\chi}$ is valid. This is
particularly important for numerical simulations such as those in the current
study where the Reynolds and Péclet numbers are still modest. Since we have
independently measured $\nu_{\rm t}$ and $\chi_{\rm t}$ using the imposed flow
and entropy method (section 3.1) and from decay experiments (section 3.2), we
can estimate $c_{\nu}$ and $c_{\chi}$ using:

$\displaystyle c_{\nu}$ $\displaystyle=$ $\displaystyle\nu_{\rm
t}/(k_{u}^{2}/\epsilon_{K}),$ (24) $\displaystyle c_{\chi}$ $\displaystyle=$
$\displaystyle\chi_{\rm t}/(k_{u}k_{T}/\epsilon_{T}),$ (25)

where $\nu_{\rm t}$ and $\chi_{\rm t}$ are the ones obtained above in section
3.1 with the imposed field method. The results are shown in figure 4. Our
results indicate that $c_{\nu}$ and $c_{\chi}$ are highly variable and that
they depend not only on Re and Pe but also on Pr. Furthermore, for
sufficiently large Reynolds and Péclet numbers, $c_{\nu}$ and $c_{\chi}$ show
decreasing trends proportional to roughly $-0.25$ power of Re and Pe,
respectively. This shows that any estimate of $\nu_{\rm t}$ or $\chi_{\rm t}$
with the $k-\epsilon$ model in the parameter regime studied here would require
prior knowledge of $c_{\nu}$ and $c_{\chi}$ for the particular parameters (Re,
Pe, Pr) of that system.

Figure 4: Similar to figure 1 but for $c_{\nu}(\mbox{\rm Re})$ and
$c_{\kappa}(\mbox{\rm Pe})$ from (24) and (25). The colours again indicate the
microscopic Prandtl number. The grey dotted lines indicate fits to the data
for $(\mbox{\rm Pe},\mbox{\rm Re}>10)$ and a line proportional to Pe is shown
for low Pe. Figure 5: Turbulent Prandtl number $\mbox{\rm Pr}_{\rm t}=\nu_{\rm
t}/\chi_{\rm t}$ as a function of Péclet number. The colour of the symbols
denotes the molecular Prandtl number as indicated by the colour bar. The
crosses ($\times$) show results from decay experiments. Linear and power law
fits to data for $\mbox{\rm Pe}>10$ are shown by the dashed and dotted lines,
respectively, and a line proportional to $\mbox{\rm Pe}^{-1}$ is shown for low
Pe.

###  3.5 Turbulent Prandtl number

Our results for the turbulent Prandtl number

$\displaystyle\mbox{\rm Pr}_{\rm t}=\frac{\nu_{\rm t}}{\chi_{\rm t}}$ (26)

are shown in figure 5 with $\nu_{t}$ and $\chi_{t}$ as discussed in sections
3.1 and 3.2. We find that for $\mbox{\rm Pe}\lesssim 1$ the turbulent Prandtl
number is roughly inversely proportional to Pe for low molecular Prandtl
number. We have not computed the turbulent Prandtl number for cases where both
Re and Pe are smaller than unity. For sufficiently high Péclet number the
turbulent Prandtl number tends to a constant value which is close to 0.7. This
is in accordance with theoretical estimates in that $\mbox{\rm Pr}_{\rm t}$ is
somewhat smaller than unity. For example, Rüdiger (1989) derived $\mbox{\rm
Pr}_{\rm t}=2/5$ using first-order smoothing approximation. The turbulent
Prandtl number plays an important role also in the atmospheric boundary layer
where several methods yield values of the order of unity (Li, 2019, and
references therein).

That, the turbulent Prandtl number $\mbox{\rm Pr}_{\rm t}$ reaches a constant
value at sufficiently large Pe, independent of Pr, is in stark contrast to the
results obtained from the $k-\epsilon$ model with a fixed $c_{\nu}/c_{\chi}$;
compare figures 3 and 5. Now we make an attempt to understand the reason for
this discrepancy. From figure 4 we note the following approximate scaling
relations at sufficiently large Re and Pe: $c_{\nu}\propto\mbox{\rm
Re}^{-0.25}$ and $c_{\chi}\propto\mbox{\rm Pe}^{-0.25}$, suggesting thus that
the ratio $c_{\nu}/c_{\chi}$ scales with the Prandtl number as $\mbox{\rm
Pr}^{+0.25}$. With this, if we let
$c^{\prime}_{\nu}/c^{\prime}_{\chi}\propto\mbox{\rm Pr}^{+0.25}$ in 23,
instead of a fixed ratio, and note from the inset of figure 3 that the factor
$k_{u}\epsilon_{T}/k_{T}\epsilon_{u}\propto\mbox{\rm Pr}^{-0.25}$, we would
obtain from 23 that $\mbox{\rm Pr}_{\rm t}^{(k-\epsilon)}$ becomes independent
of Pr, agreeing thus qualitatively with our results as shown in figure 5.
Therefore we conclude that the results from the $k-\epsilon$ model with a
fixed value for the ratio $c_{\nu}/c_{\chi}$ are unreliable and that the
strong dependence of $\mbox{\rm Pr}_{\rm t}$ on Pr found in Pandey et al.
(2021) is due to the restrictive assumption in their model.

##  4 Conclusions

Using simulations of weakly compressible isotropically forced turbulence with
imposed large-scale gradients of velocity and temperature, and corresponding
decay experiments, we find that the turbulent Prandtl number $\mbox{\rm
Pr}_{\rm t}$ is roughly 0.7 and independent of the microscopic Prandtl number
Pr provided that the Peclét number is higher than about ten. This is in stark
contrast from the recent results of Pandey et al. (2021) who found that
$\mbox{\rm Pr}_{\rm t}\propto\mbox{\rm Pr}^{-1}$ from non-Boussinesq
simulations of convection. Although the physical setups are quite different,
we were able to qualitatively reproduce their finding under the strong
assumption that $c_{\nu}/c_{\chi}$ is fixed, and note that the method by which
the turbulent viscosity and thermal diffusivity were obtained in Pandey et al.
(2021) produce unreliable results even in the simpler cases considered here.
Relaxing their assumption of universality of $c_{\nu}$ and $c_{\chi}$, we find
that these depend not only on Re and Pe, respectively, but also on Pr. This
allows us to understand the reason for the discrepancy.

[Acknowledgments]We thank S. Sridhar for his comments on an earlier draft of
the manuscript. NS thanks the hospitality provided by RRI, Bangalore where
parts of this work were done. We gratefully acknowledge the use of high
performance computing facilities at HLRN in Göttingen and Berlin, and IUCAA,
Pune.

[Funding]This work was supported by the Deutsche Forschungsgemeinschaft
Heisenberg programme (P.J.K., grant number KA4825/4-1).

[Declaration of interests]The authors report no conflict of interest.

[Author ORCID]P. J. Käpylä, https://orcid.org/0000-0001-9619-0053; N. K.
Singh, https://orcid.org/0000-0001-6097-688X.

[Author contributions]Both authors contributed equally to running the
simulations, analyzing data, reaching conclusions, and in writing the paper.

## References

  * Augustson et al. (2019) Augustson, K. C., Brun, A. S. & Toomre, J. 2019 Rossby and Magnetic Prandtl Number Scaling of Stellar Dynamos. ApJ 876 (1), 83.
  * Bekki et al. (2017) Bekki, Y., Hotta, H. & Yokoyama, T. 2017 Convective Velocity Suppression via the Enhancement of the Subadiabatic Layer: Role of the Effective Prandtl Number. ApJ 851, 74, arXiv: 1711.05960.
  * Brandenburg (2001) Brandenburg, A. 2001 The Inverse Cascade and Nonlinear Alpha-Effect in Simulations of Isotropic Helical Hydromagnetic Turbulence. ApJ 550, 824–840, arXiv: arXiv:astro-ph/0006186.
  * Brandenburg et al. (2009) Brandenburg, A., Svedin, A. & Vasil, G. M. 2009 Turbulent diffusion with rotation or magnetic fields. MNRAS 395, 1599–1606, arXiv: 0901.2112.
  * Dobler et al. (2006) Dobler, W., Stix, M. & Brandenburg, A. 2006 Magnetic Field Generation in Fully Convective Rotating Spheres. ApJ 638, 336–347, arXiv: arXiv:astro-ph/0410645.
  * Hanasoge et al. (2012) Hanasoge, S. M., Duvall, T. L. & Sreenivasan, K. R. 2012 Anomalously weak solar convection. Proc. Natl. Acad. Sci. 109, 11928–11932, arXiv: 1206.3173.
  * Käpylä et al. (2020) Käpylä, P. J., Rheinhardt, M., Brandenburg, A. & Käpylä, M. J. 2020 Turbulent viscosity and magnetic Prandtl number from simulations of isotropically forced turbulence. A&A 636, A93, arXiv: 1901.00787.
  * Karak et al. (2018) Karak, B. B., Miesch, M. & Bekki, Y. 2018 Consequences of high effective Prandtl number on solar differential rotation and convective velocity. Physics of Fluids 30 (4), 046602, arXiv: 1801.00560.
  * Li (2019) Li, Dan 2019 Turbulent prandtl number in the atmospheric boundary layer - where are we now? Atmospheric Research 216, 86–105.
  * O’Mara et al. (2016) O’Mara, B., Miesch, M. S., Featherstone, N. A. & Augustson, K. C. 2016 Velocity amplitudes in global convection simulations: The role of the Prandtl number and near-surface driving. Adv. Space Res. 58, 1475–1489, arXiv: 1603.06107.
  * Ossendrijver (2003) Ossendrijver, M. 2003 The solar dynamo. A&A Rev. 11, 287–367.
  * Pandey et al. (2021) Pandey, Ambrish, Schumacher, Jörg & Sreenivasan, Katepalli R. 2021 Non-Boussinesq convection at low Prandtl numbers relevant to the Sun. Physical Review Fluids 6 (10), 100503, arXiv: 2108.04472.
  * Pencil Code Collaboration et al. (2021) Pencil Code Collaboration, Brandenburg, Axel, Johansen, Anders, Bourdin, Philippe, Dobler, Wolfgang, Lyra, Wladimir, Rheinhardt, Matthias, Bingert, Sven, Haugen, Nils, Mee, Antony, Gent, Frederick, Babkovskaia, Natalia, Yang, Chao-Chin, Heinemann, Tobias, Dintrans, Boris, Mitra, Dhrubaditya, Candelaresi, Simon, Warnecke, Jörn, Käpylä, Petri, Schreiber, Andreas, Chatterjee, Piyali, Käpylä, Maarit, Li, Xiang-Yu, Krüger, Jonas, Aarnes, Jørgen, Sarson, Graeme, Oishi, Jeffrey, Schober, Jennifer, Plasson, Raphaël, Sandin, Christer, Karchniwy, Ewa, Rodrigues, Luiz, Hubbard, Alexander, Guerrero, Gustavo, Snodin, Andrew, Losada, Illa, Pekkilä, Johannes & Qian, Chengeng 2021 The Pencil Code, a modular MPI code for partial differential equations and particles: multipurpose and multiuser-maintained. The Journal of Open Source Software 6 (58), 2807.
  * Rüdiger (1989) Rüdiger, G. 1989 Differential Rotation and Stellar Convection. Sun and Solar-type Stars. Berlin: Akademie Verlag.
  * Schumacher & Sreenivasan (2020) Schumacher, Jörg & Sreenivasan, Katepalli R. 2020 Colloquium: Unusual dynamics of convection in the Sun. Reviews of Modern Physics 92 (4), 041001.
  * Sur et al. (2008) Sur, S., Brandenburg, A. & Subramanian, K. 2008 Kinematic $\alpha$-effect in isotropic turbulence simulations. MNRAS 385, L15–L19, arXiv: 0711.3789.


REDISCUSSION OF ECLIPSING BINARIES. PAPER XVII.

THE F-TYPE TWIN SYSTEM CW ERIDANI

By Stephen Overall and John Southworth

Astrophysics Group, Keele University, Staffordshire, ST5 5BG, UK

CW Eri is a detached eclipsing binary system of two F-type stars with an
orbital period of 2.728 d. Light curves from two sectors of observations with
the Transiting Exoplanet Survey Satellite (TESS) and previously published
radial velocity data are analysed to determine the system’s physical
properties to high precision. We find the masses of the two stars to be
$1.568\pm 0.016~{}{\rm M}_{\odot}$ and $1.314\pm 0.010~{}{\rm M}_{\odot}$, the
radii to be $2.105\pm 0.007~{}{\rm R}_{\odot}$ and $1.481\pm 0.005~{}{\rm
R}_{\odot}$ and the system’s orbit to have an eccentricity of $0.0131\pm
0.0007$. The quality of the TESS photometry allows the definition of a new
high-precision orbital ephemeris, however no evidence of pulsation is found.
We derive a distance to the system of $191.7\pm 3.8\mathrm{\,pc}$, a value
consistent with the Gaia DR3 parallax which yields a distance of
$187.9^{+0.6}_{-0.9}$ pc. The measured parameters of both stellar components
are found to be in agreement with theoretical predictions for a solar chemical
composition and an age of 1.7 Gyr.

## Introduction

Detached eclipsing binaries (dEBs) are a vital source of stellar parameters as
they allow direct measurement of the component stars’ physical properties when
combining light curves and radial velocity (RV) observations 1; 2; 3. Detached
systems are particularly useful as, in the absence of mass transfer, the
components are representative of single stars and are therefore an invaluable
source of data for testing and refining stellar evolution models 4; 5.

The volume and quality of light curve data has increased enormously in recent
years 6, especially from space-based exoplanet surveys such as CoRoT 7 and
NASA’s _Kepler_ 8, K2 9 and TESS 10 (Transiting Exoplanet Survey Satellite)
missions. This work is one of a series where we revisit known dEBs in order to
refine their characterisation with the benefit of this new era of photometry.
Here we analyse CW Eridani using TESS light curves alongside previously
published RV data.

## The dEB CW Eridani

HD 19115 was categorised as photometrically variable in 1967 by Strohmeier and
Ott 11. Popper 12 reported that the spectrum was double-lined and it was given
the designation CW Eri by Kukarkin et al. 13. Chen 14 reported on _UBV_
photometric observations made at the Rosemary Hill Observatory between 1970
and 1972, publishing an ephemeris and relative sizes for the components.
Further photometric observations were made by Mauder and Ammann 15 and, with
the addition of spectroscopic observations made available to them by Popper,
they recorded masses and radii for both components to 2–3% confidence and
assigned a spectral type of F0.

Popper and Dumont 16 included CW Eri in their program of _UBV_ photometric
observations at the Palomar and Kitt Peak observatories with the _B_ \- and
_V_ -band magnitudes given in Table I being recorded over 11 nights.
Brancewicz & Dworak 17 included it in their catalogue of eclipsing binaries
where they used numerical methods to characterise the system, giving a
spectral type of F0+ and physical parameters to $\sim$5% confidence.

The most complete investigations were carried out by Popper 12; 18 based on
spectrograms taken at the Lick Observatory between 1967 and 1974 along with
photometry from Chen 14 and Mauder & Ammann 15. Popper determined the spectral
types of the components as F1 and F4, gave their masses to 2% confidence, and
determined their radii to 2.5% and 4.5% for the primary and secondary
components respectively.

Outside automated surveys, in the years since Popper few observations have
been made. Wolf & Kern 19 recorded three observations as part of their
photometric survey of the southern hemisphere, giving a _V_ -band magnitude
ranging from $8.39$ at quadrature to $8.90$ during primary eclipse. Perry &
Christodoulou 20 included it in their _uvby $\beta$_ interstellar reddening
survey of the southern hemisphere. Nordström et al. 21 made three
spectroscopic observations as part of their RV survey of early F-type dwarfs.

Table I: Basic information on CW Eri.  Property | Value | Reference
---|---|---
Right ascension (J2000) | 03:03:59.95 | 22
Declination (J2000) |  $-$17:44:16.06 | 22
Henry Draper designation | HD 19915 | 23
Hipparcos designation | HIP 14273 | 24
Tycho designation | TYC 5868-881-1 | 25
Gaia DR3 designation | 5152756553745197952 | 22
Gaia DR3 parallax |  $5.2380\pm 0.0198$ mas | 22
TESS Input Catalog designation | TIC 98853987 | 26
$B$ magnitude | $8.79\pm 0.07$ | 16
$V$ magnitude | $8.43\pm 0.07$ | 16
$G$ magnitude | $8.306\pm 0.003$ | 22
$J$ magnitude | $7.799\pm 0.020$ | 27
$H$ magnitude | $7.659\pm 0.034$ | 27
$K_{s}$ magnitude | $7.626\pm 0.023$ | 27
Spectral type | F2 V | 28

Table I shows basic information for CW Eri. The _B_ and _V_ magnitudes are
those recorded by Popper and Dumont 16. These were explicitly based on
observations made outside of an eclipse and have since been widely used. The
_J_ , _H_ and _$K_{s}$_ magnitudes are those reported by 2MASS from
observations made at JD $2\,451\,052.9027\pm 30\sec$. At this time the system
will have been within a secondary eclipse so these will be below the system’s
maximum brightness. The spectral type of F2 V is given by Houk & Smith-Moore
28 as part of the Michigan Catalogue of HD Stars, Vol 4.

## Observational material

FIG. 1:  TESS short-cadence SAP photometry of CW Eri from sectors 4 (top) and
31 (bottom). The flux measurements have been converted to magnitude units then
rectified to zero magnitude by the subtraction of low-order polynomials.

CW Eri was observed twice by the TESS mission 10, first in sector 4 from
2018/10/19 to 2018/11/14 and again in sector 31 from 2020/10/22 to 2020/11/18,
each in short cadence mode with a $120$ s sampling rate. Both sectors show
light curves covering a period of approximately 25 days with a break near the
midpoint for data download. Unambiguous primary and secondary eclipses are
seen in addition to a sinusoidal variation resulting from the ellipsoidal
effect (Fig. 1).

The TESS time series data for the two sectors were downloaded from the MAST
archive***Mikulski Archive for Space Telescopes,
https://mast.stsci.edu/portal/Mashup/Clients/Mast/Portal.html and subsequently
processed using the lightkurve 29 and astropy 30 Python packages. These data
consist of simple aperture photometry (SAP) and pre-search data conditioning
SAP (PDCSAP) flux measurements 31. We based our analysis on the SAP data as it
is well-behaved whereas extraneous variability was seen in the PDCSAP data
from sector 4. Data points with no flux value recorded (NaN) and those with a
non-zero QUALITY flag were cut, as were those within a distorted secondary
eclipse within sector 4 from BJD $2\,458\,420.0$ to $2\,458\,423.0$. A total
of $13\,841$ data points from sector 4 and $16\,671$ from sector 31 were
considered for subsequent analysis.

The Gaia DR3 database†††https://vizier.cds.unistra.fr/viz-
bin/VizieR-3?-source=I/355/gaiadr3 was queried for potential sources of third
light within 2 arcmin of CW Eri. Six of the seven objects found are at least
10 mag fainter than CW Eri in the _G_ -band so contribute negligible light.
The remaining object, TYC 5868-428-1, has a _G_ -band magnitude of 11.053 mag
with the resulting flux ratio of 0.080 being adopted as the initial value of
the fitted third light parameter in the following analysis.

FIG. 2:  Best fit to the full TESS sector 4 light curve of CW Eri using
jktebop. The primary and secondary eclipse of the first half-sector are shown
to the left and those for the second half-sector to the right. The residuals
are shown on an enlarged scale in the lower panels.

## Light curve analysis

The remaining SAP flux data were converted to magnitudes then rectified to
zero and detrended by fitting and subtracting a quadratic polynomial across
the whole of each sector. This was refined after initial attempts at fitting,
with the best results achieved by subtracting a second quadratic fit from
those data in sector 4 following the mid-sector break. The resulting light
curves, shown in Fig. 1, consist of four isolated half-sectors over a time
interval of $\sim$759 d. We adopt the standard definition of the primary
eclipse as being the deeper of the two which occurs when the larger and
brighter component, which we label star A, is eclipsed by the smaller star B.

FIG. 3:  Same as Fig. 2 but for TESS data from sector 31.

The data were fitted using version 43 of the
jktebop‡‡‡http://www.astro.keele.ac.uk/jkt/codes/jktebop.html code 32; 33 with
a total $30\,512$ datapoints fitted as the four separate half-sectors. Each
light curve was fitted for the orbital period ($P$) and the time of mid-
primary eclipse ($T_{0}$) with our reference time being the primary eclipse
closest to the midpoint of the data. Also fitted were the sum ($r_{\rm
A}+r_{\rm B}$) and ratio ($k=r_{\rm B}/r_{\rm A}$) of the fractional radii,
the orbital inclination ($i$), the orbital eccentricity ($e$) and argument of
periastron ($\omega$) through the Poincaré elements $e\cos{\omega}$ and
$e\sin{\omega}$, the stars’ central surface brightness ratio ($J$), the amount
of third light ($L_{3}$) and each star’s reflected light.

Table II:  The fitted parameters of CW Eri for each of the four TESS half-sector light curves using the jktebop code. The uncertainties are 1$\sigma$ values derived from either Monte Carlo or residual-permutation simulations. For each parameter the uncertainties given are from the method yielding the larger weighted mean uncertainty across the half-sectors. 2400000 has been subtracted from the times to save space.  Parameter | Sector 4.1 | Sector 4.2 | Sector 31.1 | Sector 31.2
---|---|---|---|---
Fitted parameters: |  |  |  |
Primary eclipse time (BJDTDB) | $58415.482929\pm 0.000022$ | $58431.853144\pm 0.000015$ | $59152.142894\pm 0.000015$ | $59163.056373\pm 0.000020$
Orbital period (d) | $2.7283891\pm 0.0000157$ | $2.7283696\pm 0.0000092$ | $2.7283737\pm 0.0000065$ | $2.7283805\pm 0.0000117$
Orbital inclination (∘) | $86.366\pm 0.037$ | $86.373\pm 0.033$ | $86.412\pm 0.022$ | $86.313\pm 0.040$
Sum of the fractional radii | $0.30659\pm 0.00021$ | $0.30676\pm 0.00019$ | $0.30651\pm 0.00013$ | $0.30682\pm 0.00023$
Ratio of the radii | $0.7042\pm 0.0017$ | $0.7033\pm 0.0011$ | $0.7048\pm 0.0010$ | $0.7026\pm 0.0012$
Central surface brightness ratio | $0.9262\pm 0.0065$ | $0.9203\pm 0.0061$ | $0.9309\pm 0.0041$ | $0.9203\pm 0.0068$
Third light | $-0.0025\pm 0.0019$ | $-0.0014\pm 0.0017$ | $0.0024\pm 0.0012$ | $-0.0036\pm 0.0021$
LD $c$ coefficient of star A | $0.592\pm 0.027$ | $0.622\pm 0.025$ | $0.573\pm 0.017$ | $0.611\pm 0.028$
LD $c$ coefficient of star B | $0.614\pm 0.019$ | $0.601\pm 0.018$ | $0.620\pm 0.012$ | $0.608\pm 0.021$
LD $\alpha$ coefficient of star A |  $0.4676$ (fixed)
LD $\alpha$ coefficient of star B |  $0.4967$ (fixed)
$e\cos{\omega}$ | $0.00492\pm 0.00002$ | $0.00491\pm 0.00001$ | $0.00513\pm 0.00001$ | $0.00513\pm 0.00002$
$e\sin{\omega}$ | $-0.01181\pm 0.00108$ | $-0.01290\pm 0.00069$ | $-0.01135\pm 0.00066$ | $-0.01229\pm 0.00082$
Derived parameters: |  |  |  |
Fractional radius of star A | $0.17990\pm 0.00026$ | $0.18010\pm 0.00017$ | $0.17979\pm 0.00016$ | $0.18021\pm 0.00018$
Fractional radius of star B | $0.12669\pm 0.00016$ | $0.12666\pm 0.00011$ | $0.12672\pm 0.00010$ | $0.12661\pm 0.00012$
Orbital eccentricity | $0.01279\pm 0.00099$ | $0.01381\pm 0.00064$ | $0.01245\pm 0.00060$ | $0.01332\pm 0.00075$
Light ratio $\ell_{\rm B}/\ell_{\rm A}$ | $0.4535\pm 0.0014$ | $0.4537\pm 0.0010$ | $0.4542\pm 0.0009$ | $0.4510\pm 0.0010$

We adopted the power-2 limb darkening (LD) law with TESS-specific coefficients
taken from Claret & Southworth 34. The coefficients were interpolated for star
A ($T_{\rm eff}=6840\,\mathrm{K}$ and $\log{g}=4.0$) and star B ($T_{\rm
eff}=6560\,\mathrm{K}$ and $\log{g}=4.2$), each with a solar metallicity
($\mathrm{[Fe/H]}=0.0$). For both stars, the scaling coefficient $c$ was left
free to fit and $\alpha$ was fixed.

The best fits to the light curves are shown in Figs. 2 and 3 where it can be
seen that the secondary eclipse is slightly offset from phase 0.50, confirming
a small orbital eccentricity. As F-type stars may exhibit $\gamma$ Doradus or
$\delta$ Scuti pulsations 35 the residuals of the fits were analysed with
Lomb-Scargle periodograms, but no evidence of pulsation was found.

Table III:  The adopted parameters of CW Eri derived from the four TESS half-sector light curves fitted with the jktebop code. Other than the time of primary eclipse each is the weighted mean of the corresponding fitted parameter values and 1$\sigma$ uncertainties for each fitted half-sector given in Table. II.  Parameter | Value
---|---
Fitted parameters: |
Time of primary eclipse (BJDTDB) | $2458415.482929\pm 0.000022$
Orbital period (d) | $2.7283751\pm 0.0000068$
Orbital inclination (∘) | $86.381\pm 0.042$
Sum of the fractional radii | $0.30662\pm 0.00015$
Ratio of the radii | $0.7037\pm 0.0011$
Central surface brightness ratio | $0.9262\pm 0.0057$
Third light | $-0.0002\pm 0.0030$
LD $c$ coefficient of star A | $0.593\pm 0.024$
LD $c$ coefficient of star B | $0.613\pm 0.009$
LD $\alpha$ coefficient of star A |  $0.4676$ (fixed)
LD $\alpha$ coefficient of star B |  $0.4967$ (fixed)
$e\cos{\omega}$ | $0.00502\pm 0.00013$
$e\sin{\omega}$ | $-0.01210\pm 0.00076$
Derived parameters: |
Fractional radius of star-A | $0.18000\pm 0.00020$
Fractional radius of star-B | $0.12667\pm 0.00005$
Orbital eccentricity | $0.01310\pm 0.00067$
Light ratio $\ell_{\rm B}/\ell_{\rm A}$ | $0.4532\pm 0.0015$

The final values and uncertainties for the fitted parameters of each half-
sector were separately determined with 10,000 Monte Carlo (MC) simulations 32
and a residual permutation (RP) algorithm 36, as implemented by jktebop tasks
8 and 9 respectively. The latter method successively shifts the best-fit
residuals along the light curve until they are cycled back to their initial
position. With each shift a new fit is made and the final distribution of each
fitted parameter gives an estimate of its uncertainty. While the MC
simulations are sensitive to Poisson noise, the RP algorithm is additionally
sensitive to correlated noise 36. The fitted parameters for each half-sector
are given in Table II with the uncertainties being the 1$\sigma$ values of
either the MC or RP simulations. The selection of uncertainties for each
parameter is based on the method which yields the larger weighted mean
errorbar. The adopted parameters for CW Eri, as given in Table III, are the
weighted mean and uncertainty of the corresponding fitted parameters across
the four half-sectors.

## Orbital ephemeris

Table IV: Times of published mid-eclipse for CW Eri and their residuals versus the fitted ephemeris. Orbital | Eclipse time | Uncertainty | Residual | Reference
---|---|---|---|---
cycle | (BJDTDB) | (d) | (d) |
$-4130.5$ | 2441230. | 842503 |  0. | 000400 |  $-$0. | 000439 | 14
$-4121.0$ | 2441256. | 762604 | 0. | 000400 | 0. | 000144 | 14
$-4117.0$ | 2441267. | 676804 | 0. | 000400 | 0. | 000863 | 14
$-4110.0$ | 2441286. | 774304 | 0. | 000400 |  $-$0. | 000228 | 14
$-4106.0$ | 2441297. | 687604 | 0. | 000400 |  $-$0. | 000409 | 14
$-4100.5$ | 2441312. | 692704 | 0. | 000400 |  $-$0. | 001345 | 14
$-3986.5$ | 2441623. | 728610 | 0. | 000400 | 0. | 000352 | 14
$-3981.0$ | 2441638. | 734809 | 0. | 000400 | 0. | 000515 | 14
$-3967.5$ | 2441675. | 567908 | 0. | 000400 | 0. | 000616 | 14
$-3966.0$ | 2441679. | 659808 | 0. | 000400 |  $-$0. | 000039 | 14
$2168.0$ | 2458415. | 482929 | 0. | 000109 | 0. | 000003 | This work
$2174.0$ | 2458431. | 853144 | 0. | 000075 |  $-$0. | 000004 | This work
$2438.0$ | 2459152. | 142894 | 0. | 000074 | 0. | 000002 | This work
$2442.0$ | 2459163. | 056373 | 0. | 000099 |  $-$0. | 000000 | This work

With the light curve analysis yielding consistent orbital parameters we sought
to derive a high-precision orbital ephemeris for the system. In order to base
this on the longest possible dataset, historical minima timing data for CW Eri
were obtained from the TIming DAtabase at Krakow (TIDAK) team 37. While the
majority of minima timings were given without an uncertainty, all included a
weight value between 1 and 10. Where missing, estimated uncertainties were
generated by scaling a base estimate of $0.004$ by the reciprocal of the
observation’s weight. To these data were added the primary epoch and period
from fitting each of the four half-sectors of TESS data with their
uncertainties scaled up by a factor of 5 to cover any scatter.

The existing TIDAK ephemeris§§§https://www.as.up.krakow.pl/minicalc/ERICW.HTM
was used to calculate cycle numbers and assign minima types (primary or
secondary) to the eclipses after which linear, quadratic and cubic polynomials
were fitted to reveal trends in the timings. Initial attempts at fitting the
data revealed excessive scatter from a number of sources and the final fitted
ephemeris is based only on the TESS observations and those from TIDAK with a
weight of 10. With the fitting complete the quadratic and cubic fits were
discounted, as they were poorly constrained by the data, and the following
linear ephemeris was adopted:

$\mbox{Min~{}I}={\rm BJD}_{\rm TDB}~{}2452500.37624(69)+2.72837024(27)E$ (1)

with $E$ being the cycle number since the reference time and the bracketed
values being the uncertainties in the last digit of the preceding values. The
final eclipse timing data used in this analysis are given in Table IV and the
residuals of the linear fit are shown in Fig. 4.

FIG. 4:  Observed minus calculated ($O-C$) diagram of the times of primary
minimum versus the fitted linear ephemeris. Timings from the TESS data are
shown with with filled circles and those from the literature are shown as open
squares where uncertainties have been estimated. The shaded areas indicate the
$1\sigma$ uncertainty in the ephemeris determined from these data.

## Radial velocities

FIG. 5:  RVs of CW Eri measured by Popper 12 (filled circles for star A and
open circles for star B) compared to the best-fitting spectroscopic orbits
from jktebop (solid curves). The residuals are given in the lower panels
separately for the two components.

The RV measurements originally published by Popper 12 were reanalysed. The
observations were made between 1967 and 1974 at the Lick Observatory and
consist of 38 RVs for star A and 35 for star B. Popper’s data gives RVs to one
decimal place and HJD timestamps to three decimal places and, in the absence
of uncertainties, we applied equal weighting to all measurements. The RVs were
analysed with jktebop based on the ephemeris and orbital parameters derived
from the photometric fitting with the uncertainties of the fitted results
determined using Monte Carlo simulations (see Paper VI, ref. 38).

Initial fitting was carried out with fixed values for $T_{0}$ and $P$ which
yielded results very similar to Popper’s with slightly worse rms residuals for
star B. Given the low temporal resolution of the observations, we investigated
whether allowing these parameters to be varied when fitting the RV orbits
would yield an improved fit. It was found that allowing $T_{0}$ to vary
yielded a demonstrable improvement in the fitted RV orbits with lower
uncertainties and rms residuals; these are the results reported here.

The fitted orbits are shown in Fig. 5. Parameters for the spectroscopic orbits
are given in Table V which shows them to be in good agreement with Popper 12
while having lower uncertainties and residuals. Nordström et al. 21 give only
an overall systemic velocity of $V_{\gamma}=37.16\pm
2.94\,\mathrm{km\,s^{-1}}$, based on three observations, which also agrees
well with our findings. Few works have published any further RV data on CW
Eri, with Duflot et al. 39 giving $V_{\gamma}=36.4\,\mathrm{km\,s^{-1}}$ in
their Wilson-Evans-Batten catalogue and Gontcherov 40 giving a value
$V_{\gamma}=36.8\pm 2.1\,\mathrm{km\,s^{-1}}$, potentially based on the values
published by Nordström and Duflot, with both showing some overlap with our
individual RVs.

Table V:  Spectroscopic orbits for CW Eri from the literature and from the reanalysis of the RVs in the current work. All quantities are in km s-1. Source | $K_{\rm A}$ | $K_{\rm B}$ | ${V_{\gamma}}_{\rm,A}$ | ${V_{\gamma}}_{\rm,B}$ |  $rms$ residual
---|---|---|---|---|---
Popper 12 | 98.9 | 118.0 | 36.4 | 35.7 | 1.70, 2.80
|  $\pm$ 0.3 |  $\pm$ 0.6 |  $\pm$ 0.3 |  $\pm$ 0.5 |
This work | 98.7 | 117.7 | 36.1 | 36.2 | 1.55, 2.55
|  $\pm$ 0.3 |  $\pm$ 0.5 |  $\pm$ 0.3 |  $\pm$ 0.4 |

## Physical properties of CW Eri

The physical properties of CW Eri were calculated based on the parameters
derived from the light curves (Table III), the RV fitting (Table V) and the
new ephemeris calculated above. The uncertainties for $r_{\rm A}$ and $r_{\rm
B}$ were increased to 0.2% following the recommendation from Maxted et al. 41.
Effective temperature values were taken from Popper 12 where a value for both
components has been given with accompanying uncertainty. The jktabsdim code 42
was used to calculate the system’s properties given in Table VI with
uncertainties propagated using a perturbation approach. Standard formulae 43
and the reference solar values from the IAU 44 were used.

The results show that the masses and radii are determined to a precision of
better than 1.0%, meeting the criteria for inclusion in the Detached Eclipsing
Binary Catalogue (DEBCat¶¶¶https://www.astro.keele.ac.uk/jkt/debcat/ Ref. 3).
The mass measurements are in agreement with the original values published by
Popper 12, as expected as they are based on the same RV data. The measured
radii are consistent with those from Popper ($2.08\pm 0.05$ and $1.56\pm
0.07$$~{}{\rm R}_{\odot}$) but are much more precise due to the availability
of the TESS photometry.

We determined the distance to the system based on the $B$ and $V$ apparent
magnitudes from Popper & Dumont 16 and those in the $J$, $H$ and $K_{s}$-bands
from 2MASS 27 (Table I). The 2MASS magnitudes are based on observations made
during a secondary eclipse which were corrected for by subtracting the fitted
light curve model at the corresponding phase to find revised values of
$J=7.658\pm 0.022$, $H=7.518\pm 0.035$ and $K_{s}=7.485\pm 0.025$ mag. We
searched for reliable observations made in the Cousins $R$ and $I$-band but
found none. An interstellar extinction value of $E(B\\!-\\!V)=0.013\pm 0.015$
was adopted from the stilism tool ∥∥∥https://stilism.obspm.fr and bolometric
corrections from Girardi et al. 45 were used.

Table VI: Physical properties of CW Eri defined using the nominal solar units given by IAU 2015 Resolution B3 (ref. 44).  Parameter | Star A | Star B
---|---|---
Mass ratio | $0.8385\pm 0.0047$
Semimajor axis of relative orbit ($\mathcal{R}^{\rm N}_{\odot}$) | $11.694\pm 0.034$
Mass ($\mathcal{M}^{\rm N}_{\odot}$) | 1.568 $\pm$ | 0.016 | 1.314 $\pm$ | 0.010
Radius ($\mathcal{R}^{\rm N}_{\odot}$) | 2.1048 $\pm$ | 0.0074 | 1.4812 $\pm$ | 0.0052
Surface gravity ($\log$[cgs]) | 3.9869 $\pm$ | 0.0026 | 4.2156 $\pm$ | 0.0022
Density ($~{}\rho_{\odot}$) | 0.1681 $\pm$ | 0.0011 | 0.4045 $\pm$ | 0.0027
Effective temperature (K) | 6839 $\pm$ | 87 | 6561 $\pm$ | 98
Luminosity $\log(L/\hbox{$\mathcal{L}^{\rm N}_{\odot}$})$ | 0.941 $\pm$ | 0.022 | 0.564 $\pm$ | 0.026
$M_{\rm bol}$ (mag) | 2.387 $\pm$ | 0.056 | 3.330 $\pm$ | 0.065
Distance (pc) | $191.7\pm 3.8$

The Gaia DR3 22 parallax yields a distance of $187.9^{+0.6}_{-0.9}$ pc for CW
Eri, with the renormalized unit weight error (RUWE) of 1.029 indicating a
reliable astrometric solution. This is in agreement with distances based on
the bolometric corrections of Girardi et al. 45 when applied to the $B$, $V$
and $J$-band magnitudes, with the $J$-band yielding the best match at
$188.4\pm 3.0$ pc. A similar pattern is seen when using the passband-specific
surface brightness-$T_{\rm eff}$ relations of Kervella et al. 46. Both methods
yield slightly larger distances than Gaia using the $H$ and $K_{s}$-band
magnitudes, however the weighted mean of all of the derived distances is
$191.7\pm 3.8$ pc which is in agreement with the Gaia value.

## Comparison with theoretical models

FIG. 6:  Comparison between the theoretical predictions of parsec models 47
and the measured properties of CW Eri presented here for stellar mass versus
radius, $\mathrm{T_{eff}}$ and luminosity. The ages and metal abundances of
the chosen theoretical models are given in the legend within the lower right
quadrant.

To test our results, the measured properties of the components of CW Eri were
compared with predictions of parsec theoretical stellar evolutionary models 47
in plots of mass versus radius, $T_{\rm eff}$ and luminosity. The best
agreement was found for models based on a solar composition (fractional metal
abundance of $Z=0.017$) and an age of 1.7 Gyr. This gives a very good fit to
star A with star B appearing slightly larger and more luminous than the model
predictions. Choosing a model with lower metallicity gives a closer fit to
star B’s radius and luminosity but at the expense of star A which now appears
slightly smaller and less massive than predicted and both components are
cooler than the model. The converse is found when higher-metallicity models
are used with the $T_{\rm eff}$ being most sensitive to any change. Fig. 6
shows models ranging from slightly sub-solar ($Z=0.014$) through solar to
slightly super-solar metal abundance ($Z=0.020$).

This choice of age and metallicity is further supported with a Hertzsprung-
Russell diagram showing the evolutionary tracks of parsec model stars of
$Z=0.017$ and masses 1.3, 1.45 and 1.6$~{}{\rm M}_{\odot}$ (Fig. 7). Both
components of CW Eri are consistent with model stars of similar mass which
have evolved away from the ZAMS line into the upper half of the main sequence.
Comparisons were also made with equivalent MIST 48; 49 models and evolutionary
tracks, with broadly similar results except that the star A is hotter than
predicted.

While the chosen parsec model has a good fit to the mass, radii and $T_{\rm
eff}$ of the components, we note that the metallicity is in disagreement with
published values. Perry & Christodoulou 20 give $[\mathrm{Fe/H}]=-0.32$ in
their $uvby\beta$ survey of southern hemisphere A and F-type stars. A value of
$[\mathrm{Fe/H}]=-0.39$ was published in the Geneva-Copenhagen Catalogue by
Holmberg et al. 50 and subsequently recalculated as $-0.26$ by Casagrande et
al. 51. A plausible answer to this discrepancy is that the metallicity was
calculated assuming the photometry of the system represents that of a single
star rather than the combined light of two stars of different colours.

FIG. 7:  Hertzsprung-Russell diagram showing the components of CW Eri (filled
circle for star A and open circle for star B) and selected predictions from
the parsec models 47 (dotted lines). The zero-age main sequence is shown with
a long-dashed line, and an ioschrone for an age of 1.7 Gyr with a short-dashed
line. Models for 1.30, 1.45 and 1.60$~{}{\rm M}_{\odot}$ are shown (labelled)
with a metal abundance of $Z=0.017$.

## Summary

CW Eri is a dEB consisting of a pair of F-type stars that has remained largely
ignored since it was last studied by Popper in 1983 12. We have revisited the
system, making use of two sectors of TESS photometry, and determined its
photometric parameters to high precision. By combining these results with
Popper’s original RVs we refined the spectroscopic orbits and subsequently
obtained high-quality measurements of the physical properties of the system.
The residuals were analysed for evidence of pulsations with none being found.
By combining eclipse timings from the four fitted TESS half-sectors with
archival eclipse timing data we defined a new high-precision orbital
ephemeris.

The properties of both stars were found to be consistent with parsec models
for a solar metallicity and an age of 1.7 Gyr. The evolutionary tracks show
the stars to be in the second half of their main-sequence lifetime with the
more massive star A having evolved farther from the ZAMS. With two stars of
well-constrained properties and age this system lends itself to the
calibration of future stellar models, a role which could be further enhanced
by the analysis of follow-up spectroscopy to better constrain their
atmospheric characteristics.

## Acknowledgements

We gratefully acknowledge the financial support of the Science and
Technologies Facilities Council (STFC) in the form of a PhD studentship. This
paper includes data collected by the TESS mission and obtained from the MAST
data archive at the Space Telescope Science Institute (STScI). Funding for the
TESS mission is provided by NASA’s Science Mission Directorate. STScI is
operated by the Association of Universities for Research in Astronomy, Inc.,
under NASA contract NAS 5–26555. The following resources were used in the
course of this work: the NASA Astrophysics Data System; the SIMBAD database
operated at CDS, Strasbourg, France; and the ar$\chi$iv scientific paper
preprint service operated by Cornell University. Finally, we thank Jerzy
Kreiner and Bartek Zakrzewski for providing TIDAK ephemeris data for CW Eri
via private communication.

## References

  * 1 J. Andersen, A&ARv, 3, 91, 1991.
  * 2 G. Torres, J. Andersen & A. Giménez, A&ARv, 18, 67, 2010.
  * 3 J. Southworth, in Living Together: Planets, Host Stars and Binaries (S. M. Rucinski, G. Torres & M. Zejda, eds.), 2015, Astronomical Society of the Pacific Conference Series, vol. 496, p. 321\.
  * 4 A. Claret & G. Torres, A&A, 592, A15, 2016.
  * 5 A. Tkachenko et al., A&A, 637, A60, 2020.
  * 6 J. Southworth, Universe, 7, 369, 2021.
  * 7 M. Deleuil et al., A&A, 619, A97, 2018.
  * 8 W. J. Borucki, Reports on Progress in Physics, 79, 036901, 2016\.
  * 9 S. B. Howell et al., PASP, 126, 398, 2014.
  * 10 G. R. Ricker et al., Journal of Astronomical Telescopes, Instruments, and Systems, 1, 014003, 2015.
  * 11 W. Strohmeier & H. Ott, IBVS, 195, 1, 1967.
  * 12 D. M. Popper, AJ, 88, 1242, 1983.
  * 13 B. V. Kukarkin et al., IBVS, 717, 1, 1972.
  * 14 K. Y. Chen, Acta Astron, 25, 89, 1975.
  * 15 H. Mauder & M. Ammann, Mitteilungen der Astronomischen Gesellschaft Hamburg, 38, 231, 1976.
  * 16 D. M. Popper & P. J. Dumont, AJ, 82, 216, 1977.
  * 17 H. K. Brancewicz & T. Z. Dworak, Acta Astron, 30, 501, 1980.
  * 18 D. M. Popper, ARA&A, 18, 115, 1980.
  * 19 G. W. Wolf & J. T. Kern, ApJS, 52, 429, 1983.
  * 20 C. L. Perry & D. M. Christodoulou, PASP, 108, 772, 1996.
  * 21 B. Nordström et al., A&AS, 126, 21, 1997.
  * 22 Gaia Collaboration et al., A&A, 649, A1, 2021.
  * 23 A. J. Cannon & E. C. Pickering, Annals of Harvard College Observatory, 95, 1, 1920.
  * 24 ESA (ed.), The HIPPARCOS and TYCHO catalogues. Astrometric and photometric star catalogues derived from the ESA HIPPARCOS Space Astrometry Mission, ESA Special Publication, vol. 1200, 1997.
  * 25 E. Høg et al., A&A, 355, L27, 2000.
  * 26 K. G. Stassun et al., AJ, 158, 138, 2019.
  * 27 R. M. Cutri et al., 2MASS All Sky Catalog of point sources. (The IRSA 2MASS All-Sky Point Source Catalog, NASA/IPAC Infrared Science Archive., Caltech, US), 2003.
  * 28 N. Houk & M. Smith-Moore, Michigan Catalogue of Two-dimensional Spectral Types for the HD Stars. Volume 4 (University of Michigan), 1988.
  * 29 Lightkurve Collaboration et al., ‘(Lightkurve: Kepler and TESS Time Series Analysis in Python)’, Astrophysics Source Code Library, 2018.
  * 30 Astropy Collaboration et al., ApJ, 935, 167, 2022.
  * 31 J. M. Jenkins et al., in Software and Cyberinfrastructure for Astronomy IV (G. Chiozzi & J. C. Guzman, eds.), 2016, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, vol. 9913, p. 99133E.
  * 32 J. Southworth, P. F. L. Maxted & B. Smalley, MNRAS, 351, 1277, 2004.
  * 33 J. Southworth, A&A, 557, A119, 2013.
  * 34 A. Claret & J. Southworth, A&A, 664, A128, 2022.
  * 35 K. Uytterhoeven et al., A&A, 534, A125, 2011.
  * 36 J. Southworth, MNRAS, 386, 1644, 2008.
  * 37 J. M. Kreiner, Acta Astron, 54, 207, 2004.
  * 38 J. Southworth, The Observatory, 141, 234, 2021.
  * 39 M. Duflot, P. Figon & N. Meyssonnier, A&AS, 114, 269, 1995.
  * 40 G. A. Gontcharov, AL, 32, 759, 2006.
  * 41 P. F. L. Maxted et al., MNRAS, 498, 332, 2020.
  * 42 J. Southworth, P. F. L. Maxted & B. Smalley, A&A, 429, 645, 2005.
  * 43 R. W. Hilditch, An Introduction to Close Binary Stars (Cambridge University Press, Cambridge, UK), 2001.
  * 44 A. Prša et al., AJ, 152, 41, 2016.
  * 45 L. Girardi et al., A&A, 391, 195, 2002.
  * 46 P. Kervella et al., A&A, 426, 297, 2004.
  * 47 A. Bressan et al., MNRAS, 427, 127, 2012.
  * 48 A. Dotter, ApJS, 222, 8, 2016.
  * 49 J. Choi et al., ApJ, 823, 102, 2016.
  * 50 J. Holmberg, B. Nordström & J. Andersen, A&A, 475, 519, 2007\.
  * 51 L. Casagrande et al., A&A, 530, A138, 2011.