Quintessence Dark Energy: Non-Linear Gravity

An analytically solvable non-linear f(Q) gravity model with bulk-viscous quintessence dark energy in an anisotropic Bianchi-I universe is shown to match cosmic-chronometer and Pantheon supernova data while highlighting dual-use technological opportunities for aerospace, quantum systems and healthcare.

The LupoToro Group and its Research Teams are committed to exploring advanced cosmological models with potential dual-use applications in both defense and private-sector innovation. In particular, theories of modified gravity and exotic energy components are of strategic interest for their implications in space-time engineering, advanced propulsion systems that could leverage gravitational phenomena, and non-local quantum interactions potentially mediated by cosmological fields. The following analyses quintessence dark energy model and respective findings in a non-linear f(Q) gravity theory with bulk viscosity, and does refer to and reference both internal and external findings and theories.

This model is notable for its alternative formulation of gravity (using non-metricity Q in the Teleparallel Palatini formalism) and for its inclusion of a viscous fluid simulating dark energy, which together yield a late-time accelerating universe without a cosmological constant. LupoToro Group’s strategic interest in such a model lies in its potential to deepen our understanding of how altering the geometric properties of space-time (through non-metricity rather than curvature) might permit novel manipulations of gravity and energy. An accelerated expansion driven by a quintessence field (a dynamic scalar field, as realised in this f(Q) model, could inform theoretical frameworks for controlling gravitational fields or inertia, possibly informing future propulsion technologies. Moreover, studying non-linear f(Q)gravity aligns with the Group’s interest in non-local interactions: the model’s modifications to general relativity (GR) challenge the notion of gravity as purely local geometry, raising the prospect that coupling between matter and the non-metricity field might exhibit subtler, long-range effects that could inspire new approaches to quantum entanglement or information transfer. This document maintains technical details (field equations, solutions, and observational tests) alongside commentary on their significance and potential dual-use implications.

Note: Mathematical formula and explanations herein will be denoted in text-base format; this format can be extrapolated and inserted into computer-based mathematical software, for automated calculations (testing, with full workings).

Abstract

This paper examines a cosmological model of a quintessence-type dark energy scenario within a non-linear f(Q) gravity framework, incorporating bulk viscosity in an anisotropic spacetime. A modified set of gravitational field equations is derived from a generalized action (where Q denotes the non-metricity scalar in the symmetric teleparallel formulation). An analytical solution for the cosmic scale factor is obtained under the assumption of a locally rotationally symmetric (LRS) Bianchi type-I universe with shear proportional to expansion. Key model parameters are constrained using Markov Chain Monte Carlo (MCMC) analysis on observational data from cosmic chronometers (CC) and the Pantheon Type Ia supernova sample. Using the best-fit parameters, we evaluate the cosmological dynamics: the Hubble expansion rate, the deceleration parameter, and the effective equation of state of the viscous fluid (including an anisotropic skewness parameter). The results indicate a perpetually accelerating universe without a cosmological constant, consistent with a quintessence dark energy behavior. The model yields a present-day Hubble constant of about 68 km/s/Mpc and a cosmic age of approximately 13.8 Gyr, in agreement with independent observations. Statefinder diagnostics and the Om parameter assessment suggest that the solution asymptotically approaches the standard ΛCDM model at late times while presently exhibiting mild deviations characteristic of quintessence. We discuss the implications of these findings for cosmological parameter estimation and outline how the underlying methods and insights can translate into dual-use technological innovations in aerospace, quantum computing, and healthcare domains.

Introduction

Observations in the late 1990s revealed that the universe’s expansion is accelerating, implying the dominance of an enigmatic energy component dubbed dark energy Dark energy is characterized by a large negative pressure, often quantified by its equation-of-state (EoS) parameter $\omega$ (the ratio of pressure to energy density). Empirical evidence suggests that $\omega < -\tfrac{1}{3}$ is required for cosmic acceleration. The nature of dark energy can be hypothesized through different ranges of $\omega$: a quintessence field has $-1 < \omega < -\tfrac{1}{3}$, whereas a phantom field has $\omega < -1$. Current measurements (e.g. WMAP9 and Planck satellite results) constrain the dark energy EoS around $\omega \approx -1$ within observational uncertainties, consistent with a cosmological constant (Λ) or a near-quintessence behavior.

These findings, however, also highlight tensions in our understanding of gravity. The success of General Relativity (GR) on solar-system scales does not straightforwardly extend to explaining late-time cosmic acceleration without introducing a cosmological constant or exotic fields. This has motivated numerous modified gravity theories that extend or replace GR. The prototypical extension is f(R) gravity, wherein the Einstein-Hilbert action’s Ricci scalar R is generalized to an arbitrary function f(R). Such modifications can naturally produce accelerated expansion and offer alternatives to ΛCDM, the standard cosmological model.

Among the three fundamental ways to generalize GR (incorporating curvature, torsion, or non-metricity), the symmetric teleparallel gravity approach focuses on non-metricity – essentially how the length of vectors changes under parallel transport – while setting curvature and torsion to zero. In this formulation (sometimes called Q-gravity or coincident GR), gravitation is attributed to the non-metricity scalar Qinstead of the curvature scalar R. The f(Q) gravity theory postulates a gravitational Lagrangian f(Q) that is a free function of Q. This approach has drawn recent interest: it has been shown, for example, that the successful ΛCDM cosmology can be reproduced by a simple choice like $f(Q) = Q + \text{constant}$. Subsequent studies have explored matter coupling in f(Q) gravity and provided comprehensive reviews of its theoretical implications.

In parallel, cosmologists have been investigating the influence of viscous cosmic fluids on the universe’s expansion. Viscosity in the cosmic fluid (either bulk or shear viscosity) can produce effective pressure that modifies the expansion dynamics. Bulk viscosity $\xi$, in particular, acts like an internal friction arising when the universe expands or contracts too rapidly, driving the system towards equilibrium. Prior works have considered bulk-viscous matter as a mechanism for accelerated expansion in both early and late cosmic epochs. A positive bulk viscosity coefficient leads to additional effective pressure $-3\xi H$ (with $H$ the Hubble parameter), which can damp expansion and generate entropy, consistent with the second law of thermodynamics. However, constructing consistent viscous models requires causality-preserving formulations (often using Israel-Stewart-type 2nd-order theories to avoid instantaneous propagation of perturbations).

Building on this context, the present work investigates a cosmological model where bulk viscous matter drives quintessence-like acceleration in an f(Q) gravity framework. We assume an anisotropic but homogeneous LRS Bianchi type-I spacetime as the arena for this model. This choice allows for different expansion rates along one axis versus the other two (useful for exploring potential anisotropic effects) while still being simpler than a fully general anisotropic model. By imposing a relation between shear and expansion (which effectively ties the anisotropy to the overall Hubble expansion), we obtain a tractable solution. We further assume a specific non-linear functional form for f(Q) and a phenomenological form for the bulk viscosity coefficient, then solve the modified field equations analytically.

The paper is organized as follows. Section 2 provides an overview of the f(Q) gravity framework, including the action principle and field equations in the presence of matter. Section 3 details the analytical approach: the spacetime metric, the assumptions (such as the shear-expansion proportionality), the chosen form of $f(Q)$, and the resulting exact solution for the Hubble parameter and scale factor. Section 4 describes the observational data and methods used to constrain the model parameters (using cosmic chronometers and supernova data), while Section 5 discusses the implications of these constraints for cosmological parameters and the model’s consistency with observations. Section 6 outlines the broader context and dual-use technological implications of this research, highlighting how insights from this cosmological model may inform aerospace and quantum systems as well as healthcare innovations. Finally, Section 7 presents our conclusions.

Overview of the f(Q) Gravity Framework

In f(Q) gravity, one replaces the Ricci scalar R of GR with the non-metricity scalar Q. The starting point is the action functional for the gravitational field. Figure 1 shows the action $S$ considered in this framework, which consists of the f(Q) term and the matter Lagrangian $L_m$ . This action is varied with respect to the metric (and connection) to derive the field equations:

Figure 1: Action functional in f(Q) gravity. The integration is over spacetime volume $d^4x$ with $g$ the metric determinant. The function $f(Q)$ represents a general modification of the gravity Lagrangian, and $L_m$ is the matter Lagrangian . Setting $f(Q)=Q$ recovers the GR limit (up to a boundary term), whereas non-linear choices of $f(Q)$ introduce new dynamics.

The key quantity in this theory is the non-metricity scalar $Q$, defined in terms of the non-metricity tensor $Q_{\alpha\mu\nu} \equiv \nabla_\alpha g_{\mu\nu}$ (the covariant derivative of the metric). As illustrated in Figure 2, $Q$ is obtained by contracting the non-metricity tensor with its conjugate $P^{\alpha\mu\nu}$. Here $P^{\alpha}{}{\mu\nu}$ is a rank-3 tensor introduced as the conjugate to $Q{\alpha\mu\nu}$ in the action (analogous to how the conjugate momentum is introduced in Lagrangian mechanics):

Figure 2: Definition of the non-metricity scalar $Q$ in terms of the non-metricity tensor $Q_{\alpha\mu\nu}$ and its conjugate $P^{\alpha\mu\nu}$ . This definition ensures $Q$ captures the “disformation” of spacetime due to non-metricity (with curvature and torsion set to zero). A negative sign and contraction ensure that in the limit of coincident GR ($f(Q)=Q$), $Q$ reproduces the usual relation between non-metricity and the stress-energy content.

The conjugate $P^{\alpha}{}{\mu\nu}$ is constructed from the non-metricity tensor in a way that mirrors the definition of the superpotential in teleparallel gravity. Figure 3 shows the formula for $P^{\alpha}{}{\mu\nu}$ . It is defined as a specific combination of the distortion tensor $L^{\alpha}{}{\mu\nu}$ and traces of the non-metricity, such that it yields the field equations upon variation of the action. The distortion (or disformation) tensor $L^{\alpha}{}{\mu\nu}$, given in Figure 4, represents the deviation of the affine connection from the Levi-Civita connection due to non-metricity . Together, $P^{\alpha}{}{\mu\nu}$ and $Q{\alpha\mu\nu}$ encode the geometric properties of spacetime in this theory:

Figure 3: Definition of the non-metricity conjugate tensor $P^{\alpha}{}{\mu\nu}$ . It is defined as a combination of the distortion tensor and traces of $Q{\alpha\mu\nu}$, ensuring that $P^{\alpha}{}_{\mu\nu}$ is symmetric in its last two indices (denoted by parentheses around indices $\mu$ and $\nu$). This tensor plays a role analogous to the superpotential in teleparallel gravity, and it is used to construct the gravitational field equations.

Figure 4: Definition of the distortion tensor $L^{\alpha}{}{\mu\nu}$ . This tensor measures how the affine connection (in a torsion-free, flat connection setting) deviates due to non-metricity. It is symmetric in the last two indices and is half the non-metricity tensor minus its trace part. The distortion tensor helps relate $Q{\alpha\mu\nu}$ and $P^{\alpha}{}_{\mu\nu}$. In terms of geometry, it represents the degrees of freedom of the connection associated with non-metricity.

Using the definitions above, the non-metricity scalar $Q$ can be expressed explicitly in terms of metric derivatives (partial derivatives of the metric components). Figure 5 provides an expanded form of $Q$ in terms of these derivatives . This expression is useful for practical calculations in a given spacetime since it involves the metric $g_{\mu\nu}$ and its derivatives, with no dependence on the connection (having fixed a gauge where the connection is zero, often called the coincident gauge). The form shown in Figure 5 is fully general and will later be particularized to the chosen cosmological metric:

Figure 5: Expanded expression for the non-metricity scalar $Q$ in terms of metric derivatives . Each term corresponds to contractions of the metric gradient $\nabla_{\alpha} g_{\mu\nu}$. The structure ensures that $Q$ is invariant under general coordinate transformations and encapsulates the influence of non-metricity on the gravitational dynamics.

Varying the action (Figure 1) with respect to the metric yields the modified Einstein field equations in f(Q) gravity. The resulting equations relate the generalized Einstein tensor (constructed from $f_Q \equiv df/dQ$ and $P^{\alpha}{}{\mu\nu}$) to the energy-momentum tensor of matter. Figure 6 shows the field equations in their covariant form . These equations reduce to the standard Einstein equations when $f(Q) = Q$ (so $f_Q = 1$ and higher derivatives $f{QQ}=0$). The first term in Figure 6 involves a covariant derivative of $P^{\alpha}{}{\mu\nu}$, and the second term is proportional to $g{\mu\nu} f$ (playing the role of an effective cosmological term when $f(Q)$ is not simply $Q$). The right-hand side is the matter stress-energy tensor $T_{\mu\nu}$:

Figure 6: Field equations obtained by varying the f(Q) action with respect to the metric . Here $f_Q \equiv \partial f/\partial Q$. The first term represents a divergence of the conjugate tensor $P^{\alpha}{}{\mu\nu}$ weighted by $f_Q$, the second term is a generalized “potential” term $\frac{1}{2}g{\mu\nu}f$, and the third term (proportional to $f_Q$) contains quadratic combinations of $P$ and $Q$ tensors. This equation generalizes the Einstein equation to the case of non-metricity-driven gravity.

It is often useful to raise an index in the field equations to obtain a mixed-index form. Figure 7 shows the equivalent representation of the field equations with one index raised . This form can simplify the comparison with the standard Einstein equations by yielding $T^\mu{}{\nu}$ on the right-hand side. Notably, in f(Q) gravity, the covariant derivative of the stress-energy tensor is still zero (${D}{\mu}T^\mu{}_{\nu} = 0$) as required by diffeomorphism invariance , ensuring consistency with conservation laws:

Figure 7: Equivalent form of the f(Q) field equations with mixed indices . This format emphasizes that the energy-momentum tensor $T^\mu{}{\nu}$ acts as the source for the geometry encoded in $P^\alpha{}{\mu\nu}$ and $f(Q)$. It also makes it easier to impose anisotropic conditions by comparing diagonal components $\mu=\nu$ of the field equations in non-isotropic cosmologies.

In addition to the metric variation, varying the action with respect to the affine connection (assuming independent metric and connection) yields a condition on the connection (often enforcing vanishing hypermomentum in the absence of additional fields) . In our analysis, we work in the simplest setting (the coincident gauge) where the connection is taken to be zero, and this condition is automatically satisfied.

Finally, matter enters through the stress-energy tensor $T_{\mu\nu}$ on the right-hand side of the field equations. For a given matter Lagrangian $L_m$, $T_{\mu\nu}$ is defined as shown in Figure 8 . This is the same definition as in GR: it ensures that $T_{\mu\nu}$ is symmetric (for minimally coupled matter) and covariantly conserved. In our case, $L_m$ will correspond to a viscous fluid, and we will specify the form of $T_{\mu\nu}$ appropriate for an anisotropic viscous stress:

Figure 8: Definition of the energy-momentum tensor $T_{ij}$ in terms of the matter Lagrangian $L_m$ . This variational definition ensures that $T_{ij}$ is symmetric and that its covariant divergence vanishes (when the connection field equations are satisfied). It reduces to the usual $T_{ij} = (\rho + p)u_i u_j + p,g_{ij}$ for perfect fluids when $L_m$ corresponds to a perfect fluid.

Analytical Structure and Model Assumptions

To explore cosmological solutions, we specialize to a homogeneous but anisotropic metric. We consider a locally rotationally symmetric (LRS) Bianchi type-I spacetime, which has the line element given in Figure 9 . This metric has two independent scale factors: $A(t)$ for one spatial direction (say, the $x$-axis) and $B(t)$ for the other two degenerate spatial directions ($y$ and $z$ axes). It is the simplest anisotropic generalization of the flat Friedmann-Lemaître-Robertson-Walker (FLRW) metric:

Figure 9: LRS Bianchi type-I spacetime metric adopted in this model . This metric is homogeneous and anisotropic, with $A(t)$ and $B(t)$ representing directional scale factors. For $A(t)=B(t)$, it reduces to the standard flat FLRW universe. The symmetry between $y$ and $z$ axes reflects local rotational symmetry about the $x$-axis.

Given this metric, one can compute the non-metricity scalar $Q$ using its definition (Figure 2) or the expanded formula (Figure 5). Substituting the metric (Figure 9) into the general expression, we obtain the form of $Q$ for our cosmological model. Figure 10 displays the resulting $Q$ for the LRS Bianchi-I spacetime . It depends on the time derivatives of $A$ and $B$. Notably, $Q$ is negative for expansion ($\dot{A},\dot{B}>0$), and its magnitude grows with the shear between $A$ and $B$:

Figure 10: Non-metricity scalar $Q$ specialized to the LRS Bianchi I metric . The expression is given in terms of the Hubble-like rates $\dot{A}/A$ and $\dot{B}/B$. The two terms correspond to contributions from the anisotropic expansion: the first term is twice the square of the expansion rate in the $y$–$z$ plane, and the second term involves the product of expansion rates in the $x$ and $y$–$z$ directions. This $Q$ will be used in the f(Q) field equations for our model.

Next, we specify the matter content and its stress-energy tensor. We assume a cosmic fluid that can be anisotropic in pressure and includes bulk viscosity. The energy-momentum tensor is taken diagonal in the comoving frame. Figure 11shows the general form $T^\mu{}_\nu = \mathrm{diag}[-\rho,\tilde{p}_x,\tilde{p}_y,\tilde{p}_z]$ for an anisotropic fluid . Here $\rho(t)$ is the energy density, and $\tilde{p}_x, \tilde{p}_y, \tilde{p}_z$ are effective pressures in each spatial direction. We further parametrize these pressures in terms of an equation-of-state parameter and a skewnessparameter to capture the anisotropy:

Figure 11: General form of the stress-energy tensor for an anisotropic viscous fluid in diagonal form . The energy density is $\rho$, and $\tilde{p}_x, \tilde{p}_y, \tilde{p}_z$ are pressures along the $x$, $y$, and $z$ directions. In a locally rotationally symmetric scenario, we expect $\tilde{p}_y=\tilde{p}_z$, reflecting symmetry in the $y$–$z$ plane.

We introduce an EoS parameter $\omega_v(t)$ for the viscous fluid such that the longitudinal pressure $\tilde{p}_x = \omega_v \rho$. We also define a skewness parameter $\delta_v(t)$ to characterize the deviation of pressure in the transverse directions: $\tilde{p}_y = \tilde{p}_z = (\omega_v + \delta_v)\rho$. This parametrization is inserted into the tensor as shown in Figure 12 . By construction, $\delta_v=0$ corresponds to an isotropic pressure (all directions same pressure). A positive $\delta_v$ would indicate extra pressure in the transverse directions compared to longitudinal, whereas a negative $\delta_v$ (as we will find in our case) indicates slightly lower transverse pressure:

Figure 12: Parameterization of the viscous fluid stress-energy in terms of $\omega_v$ and $\delta_v$ . We set $\tilde{p}x = \omega_v \rho$ and $\tilde{p}{y,z} = (\omega_v + \delta_v)\rho$. The skewness parameter $\delta_v$ quantifies anisotropic stress: $\delta_v=0$ recovers an isotropic fluid with equation-of-state $\omega_v$. This form is useful for solving the field equations and later interpreting the results in terms of an effective dark energy EoS.

With the above setup, we plug the metric (Figure 9) and stress tensor (Figure 12) into the f(Q) field equations (Figure 6 or 7). The field equations yield a system of differential equations relating $\dot{A},\ddot{A},\dot{B},\ddot{B}$, and the function $f(Q)$ and its derivatives. In our model, we consider a specific form of $f(Q)$ to close the system. Guided by simplicity and the desire to capture non-linear effects, we choose a quadratic form for the function.

Figure 13 shows the chosen model $f(Q) = -\alpha Q^2$ . Here $\alpha$ is a constant parameter (to be constrained by data). This form is a simple deviation from GR ($f(Q) = Q$) that still allows analytic progress. It can be viewed as the $f(Q)$ analog of a quadratic gravity theory, emphasizing strong gravity corrections when $Q$ is large (e.g., in the early universe or high curvature regimes):

Figure 13: Chosen functional form of f(Q) for this study . We take $f(Q) = -\alpha Q^2$, where $\alpha$ is an arbitrary constant. This is a simple non-linear extension of the symmetric teleparallel gravity action. The negative sign ensures that for positive $\alpha$, the correction $-\alpha Q^2$ is a small positive contribution when $Q$ is small (late universe) and grows with $Q^2$ in magnitude for earlier epochs.

To solve the modified field equations, we need to impose additional conditions since the system is underdetermined (more unknown functions than independent equations). A common technique in Bianchi models is to assume a proportionality relation between the shear and the expansion, effectively tying the anisotropy to the overall rate of expansion. We assume $\sigma \propto \theta$, where $\sigma$ is the shear scalar and $\theta$ the expansion scalar. In practice, this leads to a simple relation between the directional scale factors. Specifically, we impose $A(t) = [B(t)]^m$ for some constant $m$ (with $m\neq 1$ for anisotropy). Figure 14 shows this relation . If $m=1$, the universe is isotropic ($A=B$); for $m\neq 1$, the anisotropy is sustained in the expansion. This condition has been frequently used in anisotropic cosmologies to produce solutions that tend towards isotropy at late times if $\sigma/\theta$ decays:

Figure 14: Assumed relation between scale factors under the condition $\sigma \propto \theta$ . We set $A = B^m$, where $m$ is a constant parameter. This relation simplifies the field equations by reducing the number of independent scale factor functions. It also reflects that the anisotropy in expansion is characterized by a single parameter $m$, with $m=1$ recovering isotropy.

Given $A = B^m$, we can express both scale factors in terms of a single effective scale factor $a(t)$. We define $a(t)$ such that the spatial volume $V = A B^2 = a^3$ (hence $a$ is the geometric mean scale factor). Then $A$ and $B$ can be written in terms of $a(t)$ and $m$. Figure 15 provides these expressions . They ensure that $A B^2 = a^3$ always, and they reduce correctly in special cases (for example, if $m=0$, then $A=1$ and $B=a^{3/2}$, a highly anisotropic scenario; if $m=1$, $A=a$ and $B=a$, isotropic):

Figure 15: Solutions for $A(t)$ and $B(t)$ in terms of the average scale factor $a(t)$ . These are obtained by combining $A=B^m$ with the volume relation $V = A B^2 = a^3$. The exponents $3m/(m+2)$ and $3/(m+2)$ ensure the correct volume scaling. This parameterization will be used to rewrite the field equations in terms of $a(t)$ and its derivatives.

We now incorporate bulk viscosity into the model. Bulk viscosity enters through an effective pressure term $-3\xi(t) H$ added to the equilibrium pressure. In anisotropic contexts, it effectively contributes to each directional pressure $\tilde{p}_i$. We assume only bulk (no shear) viscosity, so the effect is isotropic in pressure but influences the dynamics through $H$. The simplest approach is to specify the bulk viscosity coefficient as a function of the Hubble parameter. Based on earlier studies , we adopt a phenomenological form: $\xi(t) = \xi_1 \dot{H} - \xi_0$, where $\xi_0$ and $\xi_1$ are constants. Figure 16 shows this chosen form of $\xi(t)$ . Such a form allows viscosity to depend on the dynamics ($\dot{H}$) and includes a constant term. It is flexible enough to fit data yet simple to handle analytically:

Figure 16: Assumed form for the bulk viscosity coefficient $\xi(t)$ . The linear combination $\xi_1 \dot{H} - \xi_0$ is chosen based on phenomenological considerations and prior usage in cosmological bulk viscous models. $\xi_0$ and $\xi_1$ are constants to be determined (with $\xi_1$ dimensionless and $\xi_0$ having units of $H$, since $\dot{H}$ has units of $H^2$). This form implies that viscosity effects grow when the Hubble rate is changing rapidly.

With these assumptions ($A=B^m$, specific $f(Q)$, and $\xi(H)$ form), the field equations can be combined and simplified into a single differential equation for $H(t)$, the average Hubble parameter. Substituting $f(Q)=-\alpha Q^2$ and the relations above, we derive an ordinary differential equation for $H(t)$. Figure 17 shows the simplified equation governing $H(t)$ . This equation is first-order (in $\dot{H}$) and nonlinear, reflecting the interplay of $f(Q)$ and viscosity. All parameters $m$, $\alpha$, $\xi_0$, $\xi_1$ appear in the coefficients. Solving this equation gives $H(t)$, from which we can integrate to get the scale factor $a(t)$:

Figure 17: Reduced differential equation for the Hubble parameter $H(t)$ after substituting the model assumptions . The terms have been combined into a form $ \dot{H} + A_1 H^2 - A_2 = 0$, where $A_1$ and $A_2$ are constants depending on $m, \alpha, \xi_0, \xi_1$. This equation is reminiscent of the Riccati equation and can be solved by standard methods or recognized as a form yielding a hyperbolic function solution.

Solving the equation in Figure 17 yields a hyperbolic (coth) solution for $H(t)$. In particular, the general solution can be written as $H(t) = k_0 \coth(k_1 t + c_0)$, with constants $k_0$, $k_1$, and an integration constant $c_0$. Figure 18presents the solution for $H(t)$ . The $\coth$ form is characteristic of models where $H$ asymptotically approaches a constant (here as $t \to \infty$, $\coth \to 1$ so $H \to k_0$). The constants $k_0$ and $k_1$ are given in terms of the original parameters ($\alpha, m, \xi_0, \xi_1$) . Physically, $k_0$ turns out to be related to the present Hubble constant, and $k_1$ controls the evolution rate of $H$:

Figure 18: Solution for the average Hubble parameter as a function of time . The functional form $H(t) = k_0 \coth(k_1 t + c_0)$ indicates that $H$ starts large at early times (coth $\approx 1/t$ for small argument) and decays to a constant $k_0$ at late times (coth $\to 1$). The constants $k_0$ and $k_1$ are positive for physically realistic scenarios (expanding universe) and are defined by combinations of $m, \alpha, \xi_0,$ and $\xi_1$ given in the text.

Integrating $H(t) = \dot{a}/a$ provides the scale factor $a(t)$. Performing the integration, we obtain $a(t)$ in terms of hyperbolic sine functions. Figure 19 shows the result for $a(t)$ . The solution features a power of $\sinh(k_1 t + c_0)$ raised to an exponent determined by the parameters. This form of $a(t)$ is consistent with a universe that undergoes a smooth transition: at early times ($t$ small), $a(t) \propto t^{1/n}$ (a decelerating power-law-like phase if $\coth \approx 1/(k_1 t)$ initially), and at late times ($t$ large), $a(t)$ grows exponentially (since $\sinh(k_1 t) \sim \frac{1}{2}e^{k_1 t}$ for large $t$, the power-law exponent times $k_1 t$ in the exponent yields accelerated expansion):

Figure 19: Solution for the scale factor $a(t)$ obtained by integrating the Hubble function . It is expressed in terms of the hyperbolic sine, indicating the presence of an early decelerating phase followed by acceleration. The exponent $\frac{36\alpha(2m+1)+\xi_1(m+2)^2}{36\alpha(2m+1)}$ is greater than 1 (for positive $\xi_1$), confirming an accelerated expansion at late times. $c_1$ is an integration constant related to the scale factor at $t=0$ (which can be set by normalization).

For comparison with observations (which are functions of redshift $z$ rather than cosmic time $t$), we invert the solution to get $H$ as a function of redshift. Using the relation $1+z = a(t_0)/a(t)$ (with $t_0$ the present time) , we derive $H(z)$. The form is algebraically cumbersome, so we do not display the full expression here; however, it can be written in terms of the parameter combination $n = \frac{72\alpha(2m+1)}{36\alpha(2m+1)+\xi_1(m+2)^2}$ for simplicity. The deceleration parameter $q(z)$ can then be obtained from $H(z)$ via $q = -1 - \frac{\dot{H}}{H^2}$, which in terms of $z$ translates to $q(z) = -1 + \frac{d\ln H}{d\ln (1+z)}$. We derive an explicit expression for $q(z)$ from our solution. Figure 20 shows the final form of the deceleration parameter as a function of redshift . The expression reveals how $q(z)$ transitions from values near 0 (no acceleration) at high $z$ to negative values at $z\approx 0$, indicating acceleration:

Figure 20: Deceleration parameter $q$ as a function of redshift in the model . It is given by a rational function of $(1+z)^n$, where $n$ encapsulates the combination of model parameters. The $-1$ term indicates the de Sitter-like asymptote at high redshift (early times $q \to -1$ in this model’s limit), and the second term shows how $q$ evolves toward that asymptote. The current values $q(0)$ will depend on $\alpha, m, \xi_1$ and are expected to be around $-0.3$ (consistent with observations of an accelerating universe ).

Finally, we derive expressions for the anisotropic pressure characteristics: the directional equation of state parameters $\omega_x(z)$ and $\omega_y(z)$ (which equal $\omega_v$ and $\omega_v + \delta_v$ respectively) and the skewness parameter $\delta_v(z)$. These are obtained by plugging the solutions for $a(t)$ (or $H(z)$) back into the formulas for $\omega_v$ and $\delta_v$ that come from the field equations (specifically, combinations of the Einstein equations along different axes). Figures 21, 22, and 23 show the resulting expressions for $\omega_x(z)$, $\omega_y(z)$, and $\delta_v(z)$ . Each of these is given as a function of redshift with parameters $m, \alpha, \xi_1$ in the combination $n$ as before. We will analyze these functions to interpret the physical properties of the viscous fluid model (for example, whether it behaves like quintessence, phantom, or standard matter at various epochs):

Figure 21: Effective equation of state parameter $\omega_x(z)$ for the viscous fluid along the $x$-direction (which we take as the longitudinal direction) . In our parameterization, $\omega_x = \omega_v$. The expression features a constant part $-\frac{3(2m+3)}{5(2m+1)}$ plus a redshift-dependent correction term proportional to $\alpha$ and $\xi_1$. This indicates that at $z=0$, $\omega_v$ will be around the constant part, and it evolves with $z$ due to the bulk viscosity and f(Q) effects.

Figure 22: Effective equation of state parameter $\omega_y(z)$ (equal to $\omega_z(z)$) for the viscous fluid along the transverse directions . $\omega_y$ differs from $\omega_x$ due to the skewness $\delta_v$. The constant part here is $-\frac{2m^2 + 8m + 5}{5(2m+1)}$, and similarly a redshift-dependent correction appears. By comparing with Figure 21, one can see the difference $\omega_y - \omega_x = \delta_v$ as a function of $z$.

Figure 23: Skewness parameter $\delta_v(z)$ describing the pressure anisotropy of the viscous fluid . The expression has a constant part $\frac{2(2 - m - m^2)}{5(2m+1)}$ plus a correction term with the same structure as in Figures 21 and 22. The sign of $\delta_v$ determines whether the fluid’s pressure is higher or lower in the transverse directions compared to the longitudinal direction. Our model yields a small negative $\delta_v$ at present, meaning slightly less pressure in y,z directions than x (though $|\delta_v| \ll |\omega_v|$ in magnitude).

With these analytical expressions derived, we are prepared to confront the model with observational data and then discuss the physical implications of the results. In the next section, we describe how we use cosmic chronometer data (which give $H(z)$ measurements) and the Pantheon Type Ia supernova dataset (which constrains the luminosity distance as a function of $z$) to fit the parameters ${H_0, \alpha, m, \xi_1, \text{and minor parameter combinations}}$ and thereby validate the model.

Empirical Constraints Using Cosmic Chronometers and Pantheon Data

To test and parametrize our model, we utilize two key observational datasets: cosmic chronometers (CC), which provide direct measurements of the Hubble parameter $H(z)$ at various redshifts, and the Pantheon Type Ia supernova compilation, which provides distance modulus (or apparent magnitude) measurements over a range of redshifts. Combining these allows simultaneous constraints on the expansion rate history and the distance-redshift relationship.

For the cosmic chronometer data, we use 31 data points of $H(z)$ covering redshifts $0.07 \le z \le 1.965$ . These measurements are obtained from differential age techniques (estimating $dz/dt$ from passively evolving galaxies) and are uncorrelated. We incorporate them via a chi-square function. Figure 24 shows the chi-square definition for the CC data fitting . In this expression, $H_{\text{obs}}(z_i)$ is the observed Hubble value at redshift $z_i$, $H_{\text{th}}(z_i)$ is the theoretical value from our model (given our parameters), and $\sigma_H(z_i)$ is the observational uncertainty. The sum runs over all CC data points:

Figure 24: Chi-square function for the cosmic chronometer $H(z)$ data . Minimizing this $\chi^2_{CC}$ yields the best-fit parameters for our model against the measured expansion rates. The CC data directly constrain combinations of parameters that affect $H(z)$, notably the Hubble constant $H_0$ and parameters like $\xi_1, m, \alpha$ which alter the shape of $H(z)$ vs $z$.

For the Pantheon supernova dataset (1048 supernovae in the redshift range $0< z \lesssim 2.3$), the model predictions are encapsulated in the luminosity distance $D_L(z)$. The apparent magnitude $m(z)$ of a supernova is related to $D_L(z)$ by the distance modulus formula. Figure 25 recaps the relation between apparent magnitude $m$, absolute magnitude $M$, and luminosity distance (in megaparsecs) . We use this to compare our model’s predicted distances to the Pantheon observations. However, since there is a degeneracy between $H_0$ and $M$ (the absolute magnitude of supernovae, which is generally treated as a nuisance parameter), we follow standard practice and combine these into a single parameter. Figure 26 shows the redefinition of the absolute magnitude $M$ into a nuisance parameter $\mathcal{M}$ that absorbs $H_0$ . This approach allows the supernova analysis to be independent of the exact value of $H_0$ and instead constrain combinations of $H_0$ and $M$:

Figure 25: Definition of apparent magnitude $m(z)$ in terms of luminosity distance $D_L$ . This is the standard relation used for Type Ia supernovae, where $M$ is the absolute magnitude (assumed standardized) and $D_L$ is in units of Mpc. The constant 25 comes from $5\log_{10}(10 \text{pc}) = -5$ added to convert into the distance modulus.

Figure 26: Redefinition of the absolute magnitude to absorb $H_0$ (the Hubble constant) . The new parameter $\mathcal{M}$ is effectively the “Hubble-constant-free” absolute magnitude. When fitting the Pantheon data, $\mathcal{M}$ is varied along with other parameters, circumventing the degeneracy between $M$ and $H_0$. Any determination of $\mathcal{M}$ can later be translated back to constraints on $H_0$ given an assumed $M$ (or vice versa).

We construct a chi-square for the Pantheon data that takes into account the full covariance matrix of the binned supernova sample. The chi-square is given by $\chi^2_P = \mathbf{V}P^i (C^{-1}){ij} \mathbf{V}P^j$, where $\mathbf{V}P^i = m{\text{obs}}(z_i) - m{\text{th}}(z_i)$ is the difference between observed and model apparent magnitude for each supernova, and $C^{-1}_{ij}$ is the inverse covariance matrix of the dataset . Figure 27 provides the compact form of this expression . We use the covariance matrix provided by the Pantheon release, which includes statistical and systematic uncertainties.

Combining both datasets, the total chi-square is simply $\chi^2_{\text{tot}} = \chi^2_{CC} + \chi^2_{P}$ (assuming the CC and supernova data are independent) . We perform a joint MCMC (Markov Chain Monte Carlo) analysis to sample the parameter space and find the best-fit values along with confidence intervals. The parameters we vary are ${H_0, \xi_1, m, \alpha, \mathcal{M}}$, and we set $\xi_0$ and $c_1$ such that the model matches the current universe (i.e., $a(t_0)=1$ at $z=0$) and ensure a good fit to data. The prior ranges for these parameters were chosen broad but physically sensible (e.g., $H_0 \in [50,100]$ km/s/Mpc, $\xi_1 \in [0,1]$, $m \in [0,2]$, $\alpha \in [0.5,1.5]$, and $\mathcal{M}$ around typical supernova values).

As a result of the MCMC, we obtain the posterior distributions and confidence contours for the parameters. Figure 28shows the marginalized 2D contour plots for the key parameters $H_0, \xi_1, m,$ and $\alpha$ using only the CC dataset . Figure 29 shows the same parameters when the CC and Pantheon data are combined . In these figures (reproduced from our analysis), the inner and outer contours correspond to 1$\sigma$ and 2$\sigma$ confidence levels, respectively. The intersection of these contours gives the best-fit values. We also compile these best-fit values in Table 1 for both the CC-only and CC+Pantheon cases. Notably, the addition of Pantheon data significantly tightens the constraint on the combination of $H_0$ and $\mathcal{M}$ (thus on $H_0$ indirectly) and also helps pin down $\alpha$ more precisely, while $m$ and $\xi_1$ remain of order unity with sizable uncertainties.

(The figures and table referenced above summarize the statistical results: Figure 28 indicates that using CC data alone, $H_0$ was found around 68.2±1.3 km/s/Mpc, $\xi_1 \approx 0.17^{+0.12}{-0.09}$, $m \approx 1.05^{+0.65}{-0.58}$, and $\alpha \approx 1.09^{+0.32}{-0.26}$ . Figure 29 (with Pantheon included) gives $H_0 \approx 68.4±1.6$ km/s/Mpc, $\xi_1 \approx 0.183^{+0.098}{-0.060}$, $m \approx 1.01±0.58$, $\alpha \approx 0.96±0.30$, and the nuisance $\mathcal{M} \approx 23.85$ . Table 1 in the paper consolidates these numbers.)

In summary, the observational constraints indicate that our model is consistent with a current Hubble constant around $68$–$69$ km/s/Mpc (within errors), a non-metricity parameter $\alpha$ of order unity, and an anisotropy parameter $m$ also close to 1 (suggesting only mild anisotropy). The bulk-viscosity parameter $\xi_1$ is constrained to about $0.18^{+0.10}_{-0.06}$, implying a modest viscosity effect. The next section interprets these results in terms of cosmological parameters and discusses how the model’s behavior matches the observed universe’s acceleration and other diagnostics.

Implications for Cosmological Parameters

With the model parameters constrained by data, we can examine the behavior of key cosmological quantities to see if the model aligns with known properties of our universe. We focus on the deceleration parameter $q(z)$, the equation-of-state parameter $\omega_v(z)$ of the viscous fluid (and its current value, which indicates the effective dark energy EoS), the skewness $\delta_v(z)$ (to quantify residual anisotropy), and the inferred current age of the universe.

Accelerating expansion: The deceleration parameter $q(z)$, derived in Figure 20 and evaluated using best-fit parameters, shows a transition from positive to negative values as the universe evolves. Specifically, at high redshift (early times), $q(z) \to -1 + \frac{36\alpha(2m+1)}{36\alpha(2m+1)+\xi_1(m+2)^2}$ as $z \to \infty$ . Plugging in the best-fit values, this asymptotic value is very close to $-1$ (for instance, using $\alpha \approx 1$ and $\xi_1 \approx 0.18$, $m\approx1$, we get $q(\infty)\approx -0.99$), indicating that the model behaves like a de Sitter expansion in the remote past. At redshift $z=0$, we find $q_0 \approx -0.32$ (using the CC+Pantheon fit) . This matches well with the current observationally inferred deceleration parameter (around $-0.5$ to $-0.2$ given Planck results and supernova data, which suggest the universe entered acceleration fairly recently). The model thus produces a continuously accelerating expansion for $z<{\sim}0.6$ and a decelerating phase ($q>0$) at earlier times. The redshift at which $q(z)$ changes sign (the acceleration onset) can be calculated from our formula; it is around $z_{\text{acc}} \approx 0.7$ for the best-fit parameters, consistent with ΛCDM expectations.

We illustrate $q(z)$ in Figure 3a (from the original paper) , which plots $q$ against $z$ for the two data-fit cases (CC only and CC+Pantheon). Both curves lie between $q=0$ and $q=-1$ for $0<z<3$, and approach $q\to -1$ as $z$ increases, reflecting the model’s approach to a steady-state accelerated expansion in both past and future. The small difference between the two curves indicates that adding the Pantheon data did not dramatically change the $q(z)$ profile, only tightened the uncertainty. The present values are $q_0 \approx -0.319$ (CC) and $-0.321$ (CC+Pantheon) , which are virtually the same given the uncertainties. This confirms that the model robustly predicts a present acceleration consistent with observations.

Effective dark energy equation of state: The viscous fluid in our model plays the role of dark energy driving acceleration. We examine the EoS parameter $\omega_v(z)$ of this fluid. In our anisotropic scenario, $\omega_v$ corresponds to the pressure/density ratio along the x-direction ($\omega_x$ in Figure 21). For the best-fit parameters, we evaluate $\omega_v$ at present. The results are $\omega_v(z=0) \approx -0.451$ (CC only fit) and $\omega_v(0) \approx -0.456$ (CC+Pantheon fit) . These values lie between $-1/3$ and $-1$, indicating a quintessence-like behavior (since they are less than $-1/3$ to cause acceleration, but still greater than $-1$). In fact, the model’s current $\omega_v$ is roughly $-0.45$, which is on the higher (less negative) side compared to the cosmological constant value of $-1$. Our model yields a dark energy component slightly more “stiff” than a cosmological constant, but still adequately negative to sustain acceleration.

Figure 4a (from the original paper) plots the evolution of $\omega_v(z)$ . It shows that $\omega_v$ becomes more negative with increasing redshift (moving from about $-0.45$ at $z=0$ down toward $-1$ as $z$ approaches 3). In fact, as $z\to \infty$, $\omega_v(z)$ approaches $-1$ in our model . This means in the early universe (deep matter-dominated era and before), the effective dark energy component acted almost like a cosmological constant (though subdominant in terms of density then). As time goes on, $\omega_v$ increases (i.e., becomes less negative), crossing the phantom divide $\omega=-1$ from below? No – it approaches from above toward $-1$ as $z$ increases in our model, which suggests that in the far past it might have been slightly > -1. Actually, looking at Figure 4a, at high z, $\omega_v$ seems to approach -1 from above (e.g., -0.99). At $z=0$, it’s -0.45. So $\omega_v(z)$ increased (became less negative) as the universe expanded, which is opposite of standard quintessence which often goes from 0 to -1. Here it was nearly -1 early, now -0.45. This hints that $\omega_v$ might have “thawed” from a cosmological constant-like state to a quintessence state. Nonetheless, at all times in the range shown, $\omega_v$ stays between -0.45 and -1, thus the model never enters a phantom regime ($\omega < -1$). This is an important consistency check, as phantom energy often leads to unphysical future singularities; our model avoids such issues by remaining quintessential.

Skewness and anisotropy: The skewness parameter $\delta_v$ quantifies any residual anisotropic stress. Our results show that $\delta_v$ is very close to 0 at present, with best-fit values $\delta_v(0) \approx -0.0065$ (CC) and $-0.0013$ (CC+Pantheon) . These are tiny numbers, indicating that the pressure in y and z directions is within less than 1% of the pressure in x direction. In other words, the viscous fluid is almost isotropic in its pressure distribution today. Over time, $\delta_v(z)$ exhibits a slight increase with redshift (becoming less negative and approaching 0 from below). Figure 4b illustrates $\delta_v(z)$ , confirming it stays in a small range $(-0.02 < \delta_v < 0)$ across the redshifts of interest. At very high redshift, $\delta_v \to 0$ (the fluid becomes isotropic in the early universe), and as $z\to -1$ (looking into the far future), $\delta_v \to -0.02$ . The negative sign means the transverse pressures are slightly lower than the longitudinal pressure, but the magnitude is so small that it would be challenging to detect. Thus, while our model allowed anisotropy, the fit suggests the universe is very close to isotropic in terms of stress, consistent with observational limits on anisotropy (e.g., CMB isotropy).

In addition to $\omega_v$ and $\delta_v$, we examine the statefinder and Om diagnostic parameters to further classify the dark energy behavior. The statefinder parameters $(r,s)$ provide a diagnostic of how the model deviates from ΛCDM: for ΛCDM, $(r,s)=(1,0)$. We computed $r(z)$ and $s(z)$ as per Figure 14 (equations 51 and 52). Using best-fit parameters, we find the current values $r_0 \approx 0.305$ and $s_0 \approx 0.283$ (for CC+Pantheon) . These values are different from the ΛCDM point (0.305 vs 1 for $r$, and 0.283 vs 0 for $s$), placing our model in the “quintessence” region of the statefinder plane. Figure 6a and 6b (statefinder evolution) and Figure 7a (statefinder plane trajectory) in the original content illustrate that as $z \to \infty$, our model’s $(r,s) \to (1,0)$, meaning it approaches the ΛCDM point in the past . This is consistent with our earlier observation that at high redshift the model behaves like a cosmological constant (since matter was dominating, the dark energy component was quasi-constant). The departure at low redshift is modest but noticeable. The Om diagnostic function $Om(z)$, given by Figure 15 (equation 53), is a simpler way to test the nature of dark energy: a negative slope of $Om(z)$ indicates quintessence, a positive slope indicates phantom, and zero slope is ΛCDM. We evaluated $Om(z)$ and found it to have a slight negative slope, as shown in Figure 8 . This conclusively identifies our model’s dark energy as quintessence-like (not phantom), and moreover, $Om(z)$ tends toward a constant at high $z$, implying the model approaches ΛCDM in the future as well (because eventually the viscous effects diminish and the model converges to a constant equation of state).

Age of the Universe: Using the derived $H(z)$, we can compute the age of the universe by integrating from $z=\infty$ (the Big Bang) to $z=0$. The formula for the age $t_0$ is given by the integral in Figure 5 (equation 48) , and its solution by Figure 5 (equation 49) . Evaluating this with best-fit parameters yields $t_0 \approx 13.82$ Gyr (using CC data alone) and $13.81$ Gyr (using CC+Pantheon) . These estimates are gratifyingly close to the current best estimates of the universe’s age from Planck (about 13.8 Gyr). The difference between the two data cases is negligible, showing the result is robust. Figure 5 in the original paper shows the evolution of the cosmic age as a function of redshift , and it aligns with the standard model expectation that at high redshift, the lookback time to the Big Bang is just a few billion years, growing to ~13.8 Gyr by $z=0$. Our model, therefore, does not suffer from any age paradox and is consistent with the ages of the oldest globular clusters, etc.

Overall, the model’s implications are as follows: It provides a consistent picture of an accelerating universe driven by an effective quintessence fluid with a (current) EoS around -0.45. It matches $H_0$ and $t_0$ observations well. It predicts slight anisotropy in pressure, but at a level far below current detection limits, effectively behaving like an isotropic dark energy. The model’s additional parameters ($\alpha, \xi_1, m$) allowed it to fit the data without a cosmological constant, thereby offering an alternative explanation for acceleration through modified gravity and viscosity. Importantly, in late times the model converges toward ΛCDM behavior (the Om diagnostic flattening out, $s \to 0, r \to 1$ as $z \to -1$) . This means that in the far future, the viscous effects diminish and the expansion tends to de Sitter-like (constant $H$ eventually). Such a feature is desirable to avoid issues like big rip singularities. The fact that our model mimics ΛCDM at both early and late times, with deviations in between, makes it a viable competitor to the cosmological constant scenario, albeit with a more complex physics interpretation (modified gravity + viscosity).

In conclusion, the constrained model is not only consistent with current key cosmological observables (expansion history, age of universe, current acceleration) but also yields insight into the possible microphysics: the requirement of $\omega_v$ close to -0.5 suggests the viscous fluid behaves differently than a pure cosmological constant, potentially pointing to a dynamic field or exotic matter as the underlying cause. The small but nonzero $\delta_v$ hints that if future observations ever detected anisotropic stresses in dark energy, models like this could accommodate it. In the next section, we speculate on how the theoretical constructs and findings of this cosmological model might translate into technological innovations across different fields.

Strategic Opportunities in Dual-Use Innovation

Beyond advancing cosmology, exploring non-linear f(Q) gravity with bulk-viscous dark energy can catalyze innovations with dual-use potential in aerospace, quantum technology, and healthcare. Here we outline several strategic opportunities where the theoretical and analytical developments from this research can be translated into practical applications:

• Precision Navigation and Timing (Aerospace): The mathematical techniques used to model anistropic expansion and to fit cosmological data (e.g., MCMC parameter estimation of the expansion history) are directly applicable to aerospace engineering challenges like satellite navigation and deep-space communication. For instance, accurately modeling the expansion of space is akin to precision timing signals for GPS or interplanetary spacecraft. The f(Q)gravity framework required dealing with complex, high-precision calculations of $H(z)$; similarly, high-fidelity gravitational models are needed to account for relativistic effects in navigation systems. Insights from the f(Q) model may inspire improved gravitational potential models for Earth and planetary orbits, enhancing the accuracy of navigation systems. Additionally, the hyperbolic functions solution ($H(t) = k_0 \coth(k_1 t + c_0)$) and its behavior could inform the design of feedback control algorithms in aerospace that require sigmoidal or asymptotic response profiles. The dual-use aspect lies in the cross-pollination: methods developed for solving and constraining cosmic evolution can be adapted to optimize trajectory planning under gravitational influences or to calibrate clocks in spacecraft (leveraging the connection between cosmic time and proper time).

• Advanced Sensors and Quantum Metrology: Our study of a non-metricity-based theory and the tiny anisotropies in stress-energy touches on the detection limits of physical effects. In doing so, it pushes the boundary of sensor precision. The requirement to detect (or constrain) $\delta_v \sim 10^{-3}$ in the cosmos encourages the development of ultra-sensitive anisotropy measurements. This translates to quantum sensor development: for example, atom interferometers and superconducting gravimeters could attempt to detect minute deviations in gravity or acceleration that mimic the effects of $\delta_v$. Moreover, the f(Q) theory’s reliance on a flat, torsion-free connection suggests a natural tie-in with quantum gravity research – particularly approaches like quantum computing simulations of gravitational systems. As an example, researchers are beginning to simulate cosmological models on quantum computers . The equations we derived (such as the first-order equation for $H(t)$ in Figure 17) could be re-framed as qubit Hamiltonians in a quantum algorithm. Mastering these simulations provides dual benefits: it advances quantum computing techniques (with applications far beyond cosmology), and it offers a new way to solve or visualize complex differential equations common in many engineering fields.

• Data Analytics and Machine Learning (Healthcare & Aerospace): The heavy use of Bayesian MCMC for parameter estimation and the handling of large datasets (Pantheon supernovae, CC data) in this study directly feed into the big-data analytics domain. Techniques honed here – such as sampling high-dimensional likelihood surfaces and marginalizing over nuisance parameters – are also essential in healthcare analytics (e.g., in medical imaging or genomic data interpretation). The ability to confidently extract a signal (like cosmic acceleration) from noisy data parallels diagnosing a condition from medical scans. Indeed, there is documented overlap in algorithms for astronomy and biomedicine . LupoToro can leverage this by transferring our cosmology-driven improvements in MCMC and machine learning to its healthcare technology division, improving pattern recognition in complex datasets (like MRI images or epidemiological models). Similarly, aerospace systems (like autonomous drones or satellites) can benefit from these improved algorithms for sensor data fusion and anomaly detection, which were tuned to pick out subtle cosmological signals.

• Materials Science and Fluid Dynamics: The concept of bulk viscosity in cosmology has a direct analog in fluid engineering. Our model treats bulk viscosity as a mechanism to generate pressure that damps expansion. In aerospace engineering (jet engines, re-entry vehicles) and even in medical science (blood flow, respiratory airflow), understanding and controlling viscous dissipation is vital. By studying a cosmological bulk-viscous fluid, we developed intuition and mathematical tools (like the form $\xi(H) = \xi_1 \dot{H} - \xi_0$) that could inspire new ways to model effective viscosities in complex systems. For instance, one could design a metamaterial or smart fluid whose viscosity varies with flow acceleration (analogous to $\dot{H}$) to optimally damp vibrations or turbulences. This would be valuable in aerospace for vibration isolation or in microfluidics for biomedical devices. The dual-use benefit here is conceptual: treating viscosity in a dynamic, time-responsive way, as we did for the cosmos, might lead to adaptive damping systems in engineering.

• Cross-Disciplinary Training and Problem-Solving: From a human capital perspective, engaging researchers in solving problems like the one in this paper – which spans general relativity, quantum theory considerations (through f(Q) ties), statistical data analysis, and fluid mechanics (viscosity) – creates a talent pool adept at crossing domains. LupoToro can position such talent to switch between fundamental research and applied R&D. For example, an analyst fluent in f(Q) gravity equations might also be the right person to tackle a problem in robotics that requires non-standard geometric approaches (like SLAM algorithms on curved manifolds). Encouraging this fluidity is a strategic asset, as it fosters innovation where others might see disciplinary silos.

In conclusion, the f(Q) gravity with bulk viscosity model is more than an abstract cosmological theory; it serves as a sandbox for high-precision computation, complex system modeling, and data analysis techniques. Each of these elements has clear analogues in aerospace, quantum tech, and healthcare – the very domains of interest to LupoToro Group. By investing in this kind of fundamental research, LupoToro seeds the development of transferable skills and technologies. The rigorous demands of matching theory with cosmological data push our capabilities in ways that directly benefit advanced engineering and analytics projects. Looking ahead, we recommend interdisciplinary workshops to formally brainstorm these connections (e.g., a LupoToro symposium on “Viscosity: From Cosmic Fluids to Biomedical Fluids” or “Non-metricity and Navigation Algorithms”), and targeted pilot projects to adapt specific cosmological analysis tools to practical problems (like applying our MCMC pipeline to optimize satellite constellation configurations). Through such initiatives, the dual-use payoff of this research can be fully realized.

Conclusion

We have presented a comprehensive analysis of a quintessence-like dark energy model within the framework of non-linear f(Q) gravity, augmented by a bulk-viscous cosmic fluid. The key achievements of this work are both theoretical and observational in nature:

• Model Formulation: Starting from the symmetric teleparallel gravity (f(Q)) approach, we derived the modified field equations for an LRS Bianchi type-I universe containing a viscous fluid. By imposing a physically motivated condition relating shear and expansion ($A = B^m$) and choosing a quadratic form for $f(Q)$, we obtained an analytically solvable system. The resulting solution for the scale factor [Figure 19] is expressed in terms of hyperbolic functions, encapsulating an early decelerating phase followed by acceleration without invoking a cosmological constant term.

• Observational Consistency: Using MCMC analysis, we constrained the model with recent $H(z)$ measurements and the Pantheon supernova dataset. The best-fit parameters ($H_0 \approx 68.4$ km/s/Mpc, $\xi_1 \approx 0.18$, $m \approx 1.0$, $\alpha \approx 0.96$, etc.) show that the model fits the expansion history nearly as well as ΛCDM. The current deceleration parameter $q_0 \approx -0.32$ and derived age $t_0 \approx 13.8$ Gyr align with observations, reinforcing the viability of the scenario. Notably, the model achieves an accelerating universe (negative $q$) in the absence of a cosmological constant, thanks to the effective negative pressure provided by the bulk-viscous fluid and the modifications of gravity.

• Physical Interpretation: The dark energy in this model arises from a combination of modified gravity and bulk viscosity. We found that the effective EoS of the viscous fluid today is $\omega_v \approx -0.45$ , placing it in the quintessence regime. The slight evolution of $\omega_v(z)$ and the small skewness $\delta_v$ indicate a dynamic dark energy that was once behaving nearly like a cosmological constant and has since “thawed” to a less negative EoS. Importantly, the model does not cross into phantom ($\omega < -1$) territory, thus avoiding pathological issues. The late-time approach of statefinder parameters $(r,s) \to (1,0)$ and $Om(z)$ flattening suggests that our model asymptotically approaches the ΛCDM behavior, a reassuring feature for its long-term consistency.

• Innovation Potential: From a LupoToro Group perspective, this research has yielded advanced analytical tools and insights with cross-domain applications. The rigorous parameter estimation techniques and handling of complex systems (gravity-fluid interaction) can be translated into improvements in aerospace guidance systems, quantum simulations, and even biomedical data analysis. The collaborative effort required to connect theory with observation in cosmology mirrors the interdisciplinary approach needed to solve cutting-edge engineering problems, exemplifying the dual-use philosophy.

In summary, our study demonstrates that a non-linear f(Q) gravity model with a bulk-viscous component can provide an elegant explanation for cosmic acceleration that is consistent with empirical data. It extends the paradigm of dark energy beyond the cosmological constant, opening a window to richer interactions between geometry and cosmic fluids. While the model introduces additional parameters, these have clear physical interpretations and are constrained to reasonable values by observations, ensuring the model’s predictive power. Future work will further test this scenario against precision cosmological probes (such as CMB anisotropies and BAO measurements) and examine its behavior in the early universe (primordial nucleosynthesis and structure formation) to ensure a comprehensive viability.

From the dual-use standpoint, we recommend fostering continued interaction between fundamental cosmology and applied sciences. The transfer of knowledge—from solving for $H(t)$ in an exotic gravity theory to improving an algorithm in a navigation system—illustrates the non-linear yet fruitful pathways through which fundamental research drives innovation. LupoToro Group Research will continue to champion such integrative research, leveraging the cosmos as both a laboratory for new physics and a catalyst for technological advancement.