Global regularity for nonlinear systems with symmetric gradients

We study global regularity of nonlinear systems of partial differential equations depending on the symmetric part of the gradient with Dirichlet boundary conditions. These systems arise from variational problems in plasticity with power growth. We cover the full range of exponents $p \in (1,\infty)$. As a novelty the degenerate case for $p>2$ is included. We present a unified approach for all exponents by showing the regularity for general systems of Orlicz growth.


Introduction
In this paper we investigate the global regularity of solutions of nonlinear elliptic systems arising from problems in plasticity, see e.g.[FS00].In particular, we are interested in systems of the form (1.1) Here, Ω ⊂ R n is a bounded domain, p ∈ (1, ∞), δ ≥ 0, u : Ω → R n is a vector field and εu := 1 2 (∇u + (∇u) T ) is the symmetric part of ∇u, also abbreviated as symmetric gradient.
Similar systems arise in the context of non-Newtonian fluids (power-law fluids), where additionally a pressure and a incompressibility constraint appears, see e.g.[MNRR96].
For δ = 0 our system is also known as the symmetric p-Laplace system.Sometimes in the literature the case δ = 0 and 1 < p < 2 is called singular and the case δ = 0 and p > 2 is called degenerate.
The system (1.1) arises as the Euler-Lagrange equation of the functional where φ is the Young function φ(t) = t 0 (δ + s) p−2 s ds.
For δ = 0 we obtain φ(t) = 1 p t p .In fact, although the system (1.1) is our main interest, we also consider systems arising from more general Young functions φ.
Our technique has the advantage that we do not have to distinguish in our proofs the cases p < 2 and p > 2 as is often the case in the literature.
In the study of our system (1.1) two quantities are of special interest The former A(εu) is sometimes called the stress tensor.The latter V (εu) appears naturally in the study of regularity of nonlinear systems with power growth.Note that A and V fulfill the relation A(εu) : εu = |V (εu)| 2 .Our goal is to establish global regularity properties for solutions u of (1.1).This is best expressed in terms of V (εu).In particular, we will show that V (εu) ∈ W 1,2 (Ω) up to the boundary.This is the natural regularity when applying the technique of difference quotients which formally corresponds to testing the equation with ∆u.
We start by giving an overview of already existing results.Let us start with the simpler case of the standard p-Laplace system, which depends on the gradient ∇u rather than the symmetric gradient εu.There the proof of V (∇u) ∈ W 1,2 loc (Ω) is rather classical and can be found e.g. in [Mor08,Theorem 1.11.1],[Uhl77, Lemma 3.1] and [AF89, Proposition 2.4] or in [DE08,Theorem 11] for Orlicz growth.Global regularity V (∇u) ∈ W 1,2 (Ω) have been established in [CD89].
Let us now discuss the known results on regularity for our system (1.1) with symmetric gradients.Interior regularity V (εu) ∈ W 1,2 loc (Ω) for A(εu) roughly of the form A(εu) = (1 + |εu|) p−2 εu with p ∈ (1, 2) has been shown in [Ser92].The technique using difference quotients is similar to the one for the case of dependence on ∇u.Only in the case of global regularity V (εu) ∈ W 1,2 (Ω) the question becomes significantly more difficult.The first result on global regularity is due to Seregin and Shilkin [SS97] who proved the global regularity V (εu) ∈ W 1,2 (Ω) for the case of a flat boundary and for operators of the form A(εu) := (δ + |εu|) p−2 εu, where 1 < p ≤ 2 and δ > 0. This result was extended by Berselli and Růžička [BR20] to non-flat boundaries ∂Ω of class C 2,1 with 1 < p ≤ 2 and δ ≥ 0. The case p ≥ 2 but with δ > 0 (non-degenerate case) has been proved in [BR22], where also the parabolic situation is investigated.The estimates depend strongly on δ > 0 and deteriorate for δ → 0. In this paper we will prove V (εu) ∈ W 1,2 (Ω) for all 1 < p < ∞ and δ ≥ 0. To include all the possible variants of growth conditions we work in the even more general setting of Orlicz growth.For the special case of p-growth our main result is the following Theorem 1.1.The general case can be found in Theorem 2.4.
Our assumption of a C 2,1 -regular boundary is the same that is made in [BR20].It is a technical assumption that we need because of the way we locally flatten the boundary.In particular, we straighten the boundary locally only in one direction.A more flexible method would allow slightly less regular boundaries.
Also our assumptions on f are in some cases more stringent than necessary.This is due to the fact that φ is a general N-function and does not need to have p growth for one specific p.If it was, e.g., of p growth with 1 < p ≤ 2, then f ∈ L p ′ (Ω) would suffice.However, for p > 2 it is not sufficient to assume f ∈ L p ′ , see [BS18, Section 5].It might be possible to reduce f ∈ W 1,p ′ for p > 2 to W s,p ′ for some s ∈ ( p−2 p , 1) as done in [BS18].The condition that f vanishes at the boundary is needed in order to perform integration by parts without getting any boundary terms and is important for p ≥ 2.
Let us make a few additional remarks where the regularity theory of the symmetric p-Laplace system lacks behind the one for the p-Laplace system involving ∇u instead of εu.Higher Hölder regularity u ∈ C 1,α , in other words V (∇u) ∈ C 0,β , has been shown in [Uhl77; AF89] for the p-Laplace system and in [DSV09] for the corresponding Orlicz growth.Such regularity up to the boundary has been proved in [CD89] even for parabolic problems.Such C 1,α -regularity results are completely open for the symmetric p-Laplacian for n ≥ 3. The case n = 2 has been solved in [KMS99;DKS14], where even the p-Stokes system is covered.
Another aspect of regularity is the one of the stress A(∇u) in the case of the p- Ω) has been proved in [BCDS22].The maximal regularity A(∇u) ∈ W 1,2 (Ω) for the p-Laplace system with f ∈ L 2 (Ω) and p > 2(2− √ 2) = 1.171532 . . .has been shown in [CM19;BCDM21].The corresponding result A(εu) ∈ W 1,2 (Ω) or even A(εu) ∈ W 1,2 loc (Ω) is unsolved for the symmetric p-Laplace system.The only known result on regularity transfer from the data to the stress can be found in [DK13].There it is shown that if F ∈ L q (Ω) with p ≤ q < 3p and n = 3, then A(εu) ∈ L q loc (Ω).This is even shown for the p-Stokes system.Let us make a few final remarks on the related p-Stokes system, where an additional pressure and incompressibility constraint appear.Here V (εu) ∈ W 1,2 (Ω) is only known for systems with quadratic growth.In the language of N-functions this means φ ′′ (t) ≂ 1, see e.g.[FS00].For 1 < p ≤ 2 and δ > 0 one has V (εu) ∈ W 1,q (Ω) for some q < 2, see [BR17].
The rest of the paper is organized as follows.In Section 2 we introduce uniformly convex N-functions and formulate our main result in the language of Orlicz spaces.In Section 3 we derive Korn's and Caccioppoli's inequalities in our setting.In Section 4 we prove the main result.We do this by showing it first for N-functions with quadratic growth (Section 4.2) and then use an approximation argument to prove the general case (Section 4.3).In the Appendix A we collect some elementary properties of uniformly convex N-functions that are needed throughout the proofs.

Objective and main result
Our main focus is the study of system (1.1).This fits nicely in the more general framework of systems with Orlicz growth.In this section we introduce the corresponding partial differential system and the necessary background on Orlicz spaces and N-functions.
2.1.Uniformly convex Orlicz spaces.We begin with a bit of standard notation.With a ≲ b we mean that there exists a constant c > 0 independent of all important quantities such that a ≤ cb.We use the notation a ∨ b := max(a, b) and a ∧ b := min(a, b).For a measurable set U ⊂ R n and f ∈ L 1 loc (R n ) we denote by |U | the Lebesgue measure L n (U ) and define For a ball B ⊂ R n let 2B denote the ball with the same center but twice the diameter.Usually, we denote the radius of a ball B by r.
The systems that we study require that φ fulfills the following condition of uniform convexity, see [DFTW20].
We call p − the lower and p + the upper index (of uniform convexity) of φ.
Remark 2.3.An example for a non uniformly convex φ is given by where µ ∈ (1, ∞).In this case we have p + = 2 but p − = 1.The interior regularity of the corresponding Euler-Lagrange equation with the symmetric gradient (and Uniform convexity can also be defined for slightly less regular N-functions with φ ′ ∈ W 1,∞ loc ((0, ∞)) instead of C 2 ((0, ∞)), see [DFTW20, Definition 28].We collected some properties of uniformly convex N-functions in the Appendix, Section A.

Main result.
In this section we introduce our system and present our main result on global regularity.
We denote by ∂ i the partial derivative and by ∇ = (∂ 1 , . . ., ∂ n ) the gradient of either a scalar function, a vector field or a tensor field.If A : R n → R n×n is a function then by div A : R n → R n we mean the row-wise divergence of A, i.e. div A = n k=1 ∂ k A j,k .For matrices M, N ∈ R n×n we define the Kronecker product A : B := n i,j=1 A ij B ij and for vectors u, v we set u ⊗ v := (u i v j ) ij .For a vector field u : R n → R n we define the symmetric gradient εu via εu := 1 2 ∇u + (∇u) T .In the following let φ be a fixed uniformly convex N-function with indices of uniform convexity p − and p + .We define the two important quantities A and V as with A(P ) = V (P ) = 0 for P = 0.The quantities A, V and the N-function φ are closely related by the following formula In (2.2) and (2.3) all implicit constants depend only on p − and p + .These and other useful equivalences can be found in the Appendix, Section A.
Our main result is the following.
Theorem 2.4.Let Ω ⊂ R n be an open and bounded domain with C 2,1 -regular boundary.Let φ be a uniformly convex N-function and A and V defined as in (2.1).
Then for every f ∈ W 1,φ * 0 (Ω) the weak solution u ∈ W 1,φ (Ω) to the system of equations where c > 0 is a constant depending only on n, Ω and the upper and lower indices of uniform convexity of φ.

Basic properties of weak solutions
In this section we begin with the basic properties of weak solutions from existence to Caccioppoli estimates.We also present Korn's inequalities that are needed in the context of the symmetric gradient.
We need the following Korn-type inequality in Orlicz spaces that can be found for example in [DRS10, Theorem 6.13].
As a simple consequence we will obtain the following version for partial zero boundary values.
Proof.We will reduce the claim to an application of Theorem 3.1.For this we extend u by zero to a W 1,1 (B) function on the full ball B. Then with We have Combining this with (3.2), Theorem 3.1 and Jensen's inequality we get This proves the claim.□ In similar fashion we have the following Poincaré type inequality.The Korn and Poincaré inequality enable us to work in spaces depending on the full gradient instead of the symmetric gradient.In case of non-uniformly convex φ optimal embeddings for the spaces with the symmetric gradient can be found in [BC21].
Note that Propositions 3.1, 3.2 and 3.3 also hold for cylinders instead of balls.We will use them in this form later.
3.2.Weak solutions.For every f ∈ L φ * (Ω) there exists a unique weak solution u ∈ W 1,φ 0 (Ω) of system (2.4) in the sense that for all ψ ∈ W 1,φ 0 (Ω) we have Moreover, testing with ψ = u and using the Poincaré inequality we obtain the standard energy estimate where c > 0 depends only on n, Ω, p + and p − .

Caccioppoli estimates.
In this subsection we derive a Caccioppoli-type inequality.Before that we recall some formulas for the symmetric gradient.
The set of smooth functions u with εu = 0 is given by the space of infinitesimal rigid motions Lemma 3.4 (Caccioppoli inequality).Let u be the weak solution of system (2.4).We then have for every ball B with radius r and 2B ⊂ Ω the estimate The constants c 1 , c 2 > 0 depend only on p − , p + , n and |Ω ∩ B|/|B|.
Proof.We begin with the case 2B ⊂ Ω.Let η be a cut-off function with 1 B ≤ η ≤ 1 2B and ∥∇η∥ ∞ ≲ 1/r.For arbitrary π ∈ R we use the test function ψ = η p + (u − π).This gives us Thus, using the properties of η we get The first integral on the right hand side of (3.7) can be estimated using Young's inequality (A.3) and Lemma A.2, giving us for every λ > 0 Choosing λ > 0 small enough, we can absorb the first term into the left hand side.
The second integral on the right hand side of (3.7) can directly be estimated using Young's inequality (A.3).This shows (3.5).The boundary estimate (3.6) follows in the same way using the test function ψ = η p + u.Note that ψ = 0 on ∂Ω, since u = 0 on ∂Ω.□

Global regularity
In this section we prove our main result.Since the proof for interior regularity is similar but easier than the proof of boundary regularity, we focus in our presentation on the regularity at the boundary.First, we explain how we treat the boundary by flattening.Secondly, we show our main result for a subclass of Orlicz functions, which we call N-functions with quadratic growth, see Definition 4.1.Afterwards, we use an approximation argument to prove the general case.4.1.Flattening of the boundary.We flatten the boundary by using local isomorphisms which change Ω only in one direction.Similar to [BR17] we introduce the following notation.First we fix a point P ∈ ∂Ω and define what it means for ∂Ω to be locally C 2,1 -regular at P .For this we fix a local coordinate system in which P = 0 holds.In this coordinate system we use for every x ∈ R n the notation x ′ := (x 1 , ..., x n−1 ).Furthermore, we denote by B ′ r the (n − 1)-dimensional ball with radius r centered at 0.
We now say that ∂Ω is locally C 2,1 -regular at P if there exists 0 < r < 1 and a C 2,1 -function a : B ′ r → (−r, r) with a(0) = 0, such that ∂Ω is locally the graph of a, i.e. (a ).Furthermore, we assume that the local coordinate system is chosen such that ∇a(0) = 0 holds.This also implies that we can make ∥a∥ C 1 (B ′ r ) arbitrarily small by choosing r small enough.Finally, we say ∂Ω is C 2,1 -regular, if it is locally C 2,1 -regular at every point P ∈ ∂Ω.
For 0 < λ ≤ 1 we define the scaled sets λr , a(x ′ ) < x n < a(x ′ ) + λr}.We use the notation B + := 1 2 Ω P and 2B + := Ω P .We call them cylinders at the boundary.We use the convention that we call the first (n − 1) directions the tangential directions and denote them with Greek letters.The n-th direction we call normal direction.Note however, that while the n-th unit vector e n is indeed orthogonal to ∂Ω at P , it is only approximately orthogonal to ∂Ω at points near P .This is due to the curvature of the boundary.
Next we define curved partial derivatives along the boundary.For this we fix a tangential direction e α and define for a function f : Ω P → R and h > 0 the tangential translation via For vector-or matrix-valued functions we apply this definition component wise.If In case of a (locally) flat boundary this coincides with the partial derivative in the tangential direction e α .Otherwise, the second term accounts for the fact that we take the limit of small steps parallel to the boundary and not in the direction of any fixed unit vector.We now collect some formulas for the curved derivatives ∂.
Firstly, the curved derivatives ∂ α do not commute with the standard partial derivatives.A straight-forward calculation reveals that for f ∈ W 1,1 (Ω, R n ) we have in Ω P (4.1) If f, g ∈ W 1,1 (Ω) with supp f ∩ supp g ⊂ Ω P then we have the integration by parts formula Later on we will frequently use cut-off functions to localize either on a ball in the interior of Ω or on a cylinder at the boundary of Ω.
Let B be a ball such that 2B ⊂ Ω.Then we choose where C > 0 depends only on the dimension n.

4.2.
Boundary regularity for systems with quadratic growth.Our main result is proved for general uniformly convex N-functions.While it is in principle possible to work directly in this generality, it is difficult to justify some steps.It is more convenient to work with N-functions that still have some kind of quadratic growth in the sense that δ ≤ φ ′′ (t) ≤ δ −1 for small δ > 0. The corresponding partial differential equation then has linear growth.This allows us to work in the framework of u ∈ W 2,2 (Ω) simplifying the calculations.The estimates that we will derive will be independent of the choice of δ and depend only on p − and p + .This allows us later by some approximation process of the N-function φ to cover the general case of uniformly convex N-functions.
Note that in (2.2) and (2.3) all constants depend on φ only through p − and p + .Using a limiting process we get for every u ∈ W 2,2 (Ω) and almost every x ∈ Ω the equivalences (4.4) and where the implicit constants again only depend on p − , p + and n.This also holds for the curved derivatives in the sense that for every cylinder B + at the boundary and holds pointwise almost everywhere in B + with implicit constants depending only on p − , p + and n.
The goal is to show the main result for N-function with quadratic growth.
We already know that (2.4) admits a unique weak solution.Since φ has quadratic growth our system is of 2-structure in the sense of [BR20, Definition 2.1].Thus, using [BR20, Theorem 1.1] we get even more regularity, namely for almost every x ∈ Ω.
In order to prove Proposition 4.2, we need to control ∂ i V (εu) for every direction 1 ≤ i ≤ n.In Lemma 4.3 we do this for all directions in the interior and for all tangential derivatives ∂ α at the boundary.The main idea is to test the weak equation (3.3) with ψ = ∂ αα u.Since this is not a test function, we have to multiply it with a cut-off function.
For the normal direction i = n we can not use the same approach, since ∂ nn u does not vanish on ∂Ω.Instead, we show in Lemma 4.5 that we can use the pointwise equation (4.7) to bound the normal derivative ∂ n V (εu) by the tangential derivatives.Lemma 4.4 is an auxiliary Lemma that we need for the proof of Lemma 4.5.
Lemma 4.3.Under the assumptions of Proposition 4.2 we have for every ball B with 2B ⊂ Ω and every where c > 0 does only depend on n, p − , p + .Furthermore, if B + is a cylinder at the boundary of Ω, 1 ≤ α ≤ n − 1 and η a cut-off function supported on 2B + (see Section 4.1), then (4.9) where c > 0 does depend on n, p − , p + and the parametrization a of the boundary, while κ > 0 depends only on n, p − , p + and is independent of a .
Proof.Let us first fix a cylinder B + at the boundary of Ω. Multiplying the pointwise equation (4.7) with ψ := ∂ α (η 2 ∂ α u) and integrating over Ω we get (4.10) We estimate the right hand side in (4.10) using the integration by parts formula (4.2), Young's inequality and Korn's inequality: The left hand side in (4.10) can be rewritten using (4.1) as We now give detailed explanation of how to estimate the term T 2 which appears due to the curved boundary.Subsequently, similar terms will appear and we will just refer back to the calculations here.We have for every λ > 0 (4.11) The second integral can be estimated using the pointwise estimate |εu| ≤ |∇u| and Korn's inequality: (4.12) We now turn to T 1 .We have Using (4.5) we have We can estimate T 1,2 exactly like we estimated T 2 in (4.11).Using |∇η| ≲ r −1 and (4.12), the term T 1,3 can be bounded as follows.
Summarizing our findings so far (and choosing λ small enough), we have shown To conclude (4.9) from here we simply note that Thus, (4.9) follows.The interior estimate (4.8) can be shown in the same way but with many simplifications, since we do not need to use the curved derivatives ∂ α .□ Lemma 4.4.Let B + be a cylinder at the boundary and η a cut-off function supported on 2B + .Let 1 ≤ i, j ≤ n with (i, j) ̸ = (n, n).Then (4.13) where c > 0 depends only on p − , p + , n and the parametrization a of the boundary, while κ > 0 depends only on p − , p + , n and is independent of a.
Proof.Assume i ̸ = n.The case j ̸ = n can be treated in exactly the same way.By (4.4) we have for every λ > 0 We can now use Lemma 4.3 to bound the first integral.The choice λ = ∥∇a∥ L ∞ finishes the proof.□ Lemma 4.5.Let 2B + be a cylinder at the boundary and η a cut-off function supported on 2B + .Then (4.14) where c > 0 depends on p − , p + , n and the parametrization a of the boundary, while κ depends only on p − , p + , n and is independent of a.
Proof.We adopt the method used in [SS97].Recall that we use the convention that the greek indices α and β run only over the tangential directions 1, ..., n − 1.We have We start by estimating S 2 .Using the pointwise equation (4.7) we have The first sum can be estimated using Lemma 4.4.To estimate the second sum we use integration by parts.Because η 2 f vanishes on all of ∂(2B + ), we get no boundary terms.Thus Using Young's inequality A.3, we can estimate both summands.Indeed, we have We can bound S 3 using the same strategy.We are thus left with bounding S 1 .For this we use the identity This identity allows us to split S 1 into three terms, the first two of which can again be estimated using Lemma 4.4.This leaves us only with the third term, which we can write as (4.15) where Here we used integration by parts to move the ∂ n away from A(εu) and the ∂ α away from u n , and the fact that where we used (4.4) and (4.12).For A 2 we have Both terms can be estimated with the help of Lemma 4.3.We can estimate A 3 as in (4.12) Finally, for A 4 we have Collecting our estimates for A 1 ,...,A 4 , we get the claimed estimate for S 1 which together with the estimates for S 2 and S 3 completes the proof.□ We are now ready to prove the main Proposition of this section.
Proof of Proposition 4.2.For every cylinder B + at the boundary and every cut-off function η supported on 2B + we know by Lemma 4.3 and Lemma 4.5 that there exists a κ > 0 independent of the local parametrization a of the boundary such that If we choose B + sufficiently small we can assume that (4.16) .
With this assumption and using r < diam(Ω) and Lemma A.2 we get (4.17) The same bound holds for all interior balls B with 2B ⊂ Ω.We now cover Ω by choosing a finite collection of balls (B i ) i∈I resp.cylinders at the boundary (B + j ) j∈J such that (a) The B i together with the B + j cover Ω, (b) 2B i ⊂ Ω resp.2B + j ⊂ Ω holds for all i ∈ I resp.j ∈ J, (c) There exists a natural number C = C(n) such that every x ∈ Ω is contained in at most C of the B i or B + j .Additionally we assume that (4.16) holds for all of the B + j .We denote by η i resp.η j the cut-off function corresponding to B i resp.B + j .We have Using (4.17) and the finite overlap of our covering we get Together with the energy estimate (3.4) this implies the claimed bound.□ 4.3.Passage to the limit.Proposition 4.2 already contains the results of Theorem 2.4 but is restricted to N-functions with quadratic growth.In this section we will show how to approximate every uniformly convex N-function by those with quadratic growth and how to pass to the limit.For this let φ be a uniformly convex N-function.For every 0 ≤ δ − ≤ δ + ≤ ∞ and δ := (δ − , δ + ) we define the truncated N-functions and φ δ (t) := t 0 φ ′ δ (s) ds.These functions were introduced in [DFTW20].See the Appendix for more details.For example by Lemma A.6 we know that φ δ is a uniformly convex N-function with quadratic growth.
Accordingly, we define for We can apply Proposition 4.2 only to systems where f is contained in W 1,2 (Ω).Thus we also need to approximate f .For this we perform a Lipschitz truncation, i.e. we approximate f by functions f δ ∈ W 1,∞ (Ω) which agree with f outside of a small set.
As an auxiliary tool we need the Hardy-Littlewood maximal operator |f | dy.
To approximate f we will use the following approximation result.
The hidden constants only depend on the indices p − and p + of φ and Ω.
Proof.The uniform convexity of φ ensures that 1 < p − ≤ p + < ∞.Now, the proof follows exactly as in [DKS13;BDS16].However, the boundedness of the maximal operator M in L p has to be replaced by the corresponding results on Orlicz spaces, see [KR61]. □ We are now ready to prove Theorem 2.4.
Proof of Theorem 2.4.In order to apply Proposition 4.2 we approximate f by applying the above Lipschitz truncation is chosen in such a way that the truncation of φ at δ + (near ∞) does not see the difference between ∇f δ and ∇f .
Proof of Theorem 1.1.We restrict our proof to the case δ = 0, the case δ > 0 can be treated in the same way.The claimed estimate follows immediately from Theorem 2.4.We are thus left to prove the claimed Sobolev regularity of u.In the case p ≥ 2 we note that The other estimate follows by replacing φ with φ * .□ For N-functions φ one has for t ≥ 0 φ(t) ≤ φ ′ (t)t ≤ φ(2t).
One easily calculates This proves that ψ is a uniformly convex N-function and the indices of ψ and φ are related by the formula The quantities A, V and φ are linked by the following important equivalence.
Lemma A.5.Let φ be a uniformly convex N-function.Then for all P, Q ∈ R n×n Proof.By the definition of φ δ we have

These observations readily imply (a)-(c).
To see (d) we note that This proves (d).