Mixed virtual element methods for a stress-velocity-rotation formulation in viscoelasticity

TL;DR

The paper addresses robust numerical approximation of viscoelastic dynamics modeled by a Zener (standard linear solid) rheology using a mixed virtual element method with weakly imposed stress symmetry. It combines a polygonal-mesh compatible H(div)-conforming stress space with velocity and skew-symmetric rotation multipliers and employs Crank–Nicolson time stepping, supported by a weakly symmetric projection to enable fully computable schemes. The authors prove unique solvability for both semi-discrete and fully discrete problems and derive optimal a priori error estimates for all unknowns, backed by numerical tests on various polygonal meshes that confirm the theoretical rates and robustness, including near-incompressible regimes and hanging-node configurations. This approach yields a flexible, locking-free framework for simulating viscoelastic materials on general meshes with rigorous error control and demonstrated practical effectiveness.

Abstract

In this paper we propose a new mixed virtual element formulation for the numerical approximation of viscoelasticity equations with weakly imposed stress symmetry. The governing equations use the Zener model and are expressed in terms of the principal unknowns of additively decomposed stress into elastic and internal viscoelastic contributions, while the rotation tensor and velocity act as Lagrange multipliers. The time discretisation uses Crank--Nicolson's scheme. We demonstrate the unique solvability of both semi-discrete and fully-discrete problems by leveraging the properties of suitable local projectors. Moreover, we establish optimal a priori error estimates for all variables that appear in the mixed formulation. To validate our theoretical findings, we present several representative numerical examples that also highlight the features of the proposed formulation.

Figures (8)

• Figure 1: Schematic representation of the rheology for the classical Zener (or standard linear solid) model in viscoelasticity. $\mathcal{C}_0,\mathcal{C}_1$ are the elasticity tensors associated with the first and second spring units, respectively, while $\mathcal{C}'$ is the tensor associated with the dashpot.
• Figure 1: Schematics of the local degrees of freedom for each row of the stress variable, focusing on polynomial degree $k=2$.
• Figure 1: Sample of coarse meshes of the two types (squares - left and hexagonal - right) employed for the convergence tests.
• Figure 2: Error history for the system variables in the $L^2$-norm using the manufactured solution $\boldsymbol{u}(x,y,t)=(e^{-y}\cos(t)\sin(x),e^{t+x})^{\tt t}$, with polynomial degree $k=1$ and on square meshes (a). Using the manufactured solution $\boldsymbol{u}(x,y,t)=(t^2x(1-x)y(1-y),0)^{\tt t}$ with polynomial degree $k=1$ and on square meshes (b), on hexagonal meshes (c) on square meshes (d), with degree $k=2$ on hexagonal meshes (e) with $k=2$ on square meshes (f), with $k=3$ on square meshes (g), and with $k=3$ on hexagonal meshes (h). The solid lines indicate fitted convergence rate slopes.
• Figure 3: Convergence plots for a nearly incompressible material with Poisson's ration $0.49$ (a) and $0.4999$ (b).

Theorems & Definitions (12)

• remark thmcounterremark
• lemma thmcounterlemma
• theorem 3.1
• remark thmcounterremark
• lemma thmcounterlemma
• lemma thmcounterlemma
• remark thmcounterremark
• theorem 4.1
• theorem 4.2
• lemma thmcounterlemma