---
abstract: 'Two new approaches to solving first-order quasilinear elliptic systems of PDEs in many dimensions are proposed. The first method is based on an analysis of multimode solutions expressible in terms of Riemann invariants, based on links between two techniques, that of the symmetry reduction method and of the generalized method of characteristics. A variant of the conditional symmetry method for constructing this type of solution is proposed. A specific feature of that approach is an algebraic-geometric point of view, which allows the introduction of specific first-order side conditions consistent with the original system of PDEs, leading to a generalization of the Riemann invariant method for solving elliptic homogeneous systems of PDEs. A further generalization of the Riemann invariants method to the case of inhomogeneous systems based on the introduction of specific rotation matrices enabled us to weaken the integrability condition. It allows us to establish the connection between the structure of the set of integral elements and the possibility of the construction of specific classes of simple mode solutions. These theoretical considerations are illustrated by the examples of an ideal plastic flow in its elliptic region and a system describing a nonlinear interaction of waves and particles. Several new classes of solutions have been obtained in explicit form including the general integral for the latter system of equations.'
author:
- |
    A.M. Grundland[^1],\
    Centre de Recherche Mathématiques, Université du Montréal,\
    C.P. 6128, Succc. Centre-ville, Montréal, (QC) H3C 3J7, Canada\
    and Département de mathématiques et informatiques, Université du Québec,\
    Trois-Rivières, (QC) G9A 5H7, Canada\
- |
    V. Lamothe[^2],\
    Département de Mathématiques et Statisque, Université de Montréal,\
    C.P. 6128, Succc. Centre-ville, Montréal, (QC) H3C 3J7, Canada\
title: 'Multimode solutions of first-order quasilinear systems obtained from Riemann invariants. Part I. '
---

[Running Title: Multimode solutions of quasilinear systems. Part I.\
PACS numbers: 02.30.Jr; Secondary 02.70.-c\
Keywords: symmetry reduction method, generalized method of characteristics, Riemann invariants, multimode solutions]{}

Introduction {#intro}
============

Riemann waves represent a very important class of solution of nonlinear first-order systems of partial differential equations (PDEs). They are ubiquitous in the differential equations of mathematical physics, since they appear in all multidimensional hyperbolic systems and constitute their elementary solutions. Their characteristic feature is that, in most cases, they are expressible only in implicit form. For a homogeneous hyperbolic quasilinear system of first-order PDEs,

\[eq:1\] \_\^[i]{}(u)=0,i=1,…,p,=1,…, q,=1,…, m

(where $\mathcal{A}^1,\ldots,\mathcal{A}^p$ are $q\times m$ matrix functions of an unknown function $u$ and we adopt the convention that repeated indices are summed unless one of them is in a bracket), a Riemann wave solution is defined by the equation $u=f(r(x,u))$, where $f:{\mathbb{R}}{\rightarrow}{\mathbb{R}}^q$, and the function $r(x,u)={\lambda}_i(u)x^i$ is called the Riemann invariant associated with the vector $\lambda$ satisfying the equation $\ker{\left( {\lambda}_i\mathcal{A}^i(u) \right)}\neq 0$. These solutions have rank at most equal to one. They are building blocks for constructing more general types of solutions describing nonlinear superpositions of many waves ($k$-waves), which are very interesting from the physical point of view. Until now, the only way to approach this task was through the generalized method of characteristics (GMC) [@Burnat:1972; @CourantHilbert:1962; @Jeffrey:1976; @JohnKlainerman:1984; @Perad:1985; @Rozdestvenski:1983] and more recently through the conditional symmetry method (CSM) [@AblowitzClarkson; @ConteGrundHuard:2009; @Fushchych:1991; @GrundHuard:2006; @GrundHuard:2007; @OlverVorobev1995]. The GMC relies on treating Riemann invariants as new dependent variables (which remain constant along the appropriate characteristic curves of the initial system (\[eq:1\]) and constitute a set of invariants of the Abelian algebra of some vector field $X_a=\xi_a^i(u){\partial}_x^i$ with ${\lambda}_i^a\xi_a^i=0$ for $1\leq
a\leq k<p$. This leads to the reduction of the dimension of the problem. The most important theoretical results obtained with the use of the GMC or CSM [@GrundHuard:2007] include the finding of necessary and sufficient conditions for the existence of Riemann $k$-waves in multidimensional systems. It was shown [@Perad:1985] that these solutions depend on $k$ arbitrary functions of one variable. Some criteria were also found [@Perad:1985] for determining the elastic or nonelastic character of the superposition of Riemann waves described by hyperbolic systems, which is particularly useful in physical applications. In applications to fluid dynamics and nonlinear field theory, many new and interesting results were obtained [@Boillat:1965; @Burnat:1972; @DubrovinNovikov:1983; @FerapontovKhus:2004:1; @FerapontovKhus:2004:2; @JohnKlainerman:1984; @Mises:1958; @Rozdestvenski:1983; @Whitham:1974; @Zakharov:1998]. Both the GMC and CSM methods, like all other techniques for solving PDEs, have their limitations. This fact has motivated the authors to search for the means of constructing larger classes of multiple wave solutions expressible in terms of Riemann invariants by allowing the introduction of complex integral elements in the place of real simple integral elements (with which the solutions of hyperbolic systems are built [@Perad:1985]). This idea originated from the work of S. Sobolev [@Sobolev:1934] in which he solved the wave equation by using the associated complex wave vectors. We are particularly interested in the construction of nonlinear superpositions of elementary simple mode solutions, and the proposed analysis indicates that the language of conditional symmetries is an effective tool for this purpose. This approach is applied to the nonstationary irrotational flow of a ideal plastic material in its elliptic region. A further extension of the proposed method to the case of inhomogeneous systems is proposed in order to be applicable either in elliptic or hyperbolic regions. This allows for a wider range of physical applications. The approach is based on the use of rotation matrices which obey certain algebraic conditions and allow us to write the reduced system in terms of Riemann invariants in the sense that each derivative of dependent variables is equal to an algebraic expression (see equation (\[eq:ms:12\])). We discuss in detail the sufficient conditions for the existence of multimode solutions. This approach is applied to a system describing a propagation of shock waves intensity in the nonlinear interaction of waves and particles. The general integral of this system has been constructed in an explicit form depending on two arbitrary functions of one variable.

#### 

The organization of this paper is as follows. Section \[sec:2\] contains a detailed account of the generalized method of characteristics for first-order quasilinear systems of PDEs in many dimensions based on complex characteristic elements. In Section \[sec:3\] we formulate the problem of multimode solutions expressible in terms of Riemann invariants by means of a group theoretical approach. This allows us to formulate the necessary and sufficient conditions for constructing these types of solutions. In Section \[sec:4\], the usefulness of the method developed in Section \[sec:3\] is illustrated by an example of the ideal plasticity in $(2+1)$ dimensions, in which we find several bounded solutions. Moreover, we have drawn extrusion dies and the flow inside it (limiting ourselves to the region where the gradient catastrophe does not occur). Sections \[sec:5\] and \[sec:6\] comprise a new approach to solving inhomogeneous elliptic systems to obtain simple wave and simple mode solutions. In Section \[sec:7\], we have shown a example of a simple mode solution using the method presented in Section \[sec:6\]. The technique is applied on a system describing the propagation of shock waves for the nonlinear interaction of waves and particles. We obtain its general solution. Section \[sec:8\] summarizes the results obtained and contains some suggestions regarding further developments.

The method of characteristics for complex integral elements. {#sec:2}
============================================================

The methodological approach assumed in this section is based on the generalized method of characteristics which has been extensively developed (*e.g. in* [@Burnat:1972; @DoyleGrundland:1996; @Grundland:1974; @GrundlandTafel:1996; @Perad:1985] and references therein) for multidimensional homogeneous and inhomogeneous systems of first-order PDEs. A specific feature of that approach is an algebraic and geometric point of view. An algebraization of systems of PDEs was made possible by representing the general integral elements as linear combinations of some special elements associated with those vector fields which generate characteristic curves in the spaces of independent variables $X$ and dependent variables $U$, respectively (see [@GrundlandZelazny:1983; @Perad:1985]). The introduction of these elements (called simple integral elements) proved to be very useful for constructing certain classes of rank-$k$ solutions in closed form. These integral elements proved to correspond to Riemann wave solutions in the case of nonelliptic systems and serve to construct multiple waves ($k$-waves) as a superposition of several single Riemann waves.

#### 

The generalized method of characteristics for solving quasilinear hyperbolic first-order systems can be extended to the case of complex characteristic elements (see *e.g.* [@Perad:1985; @Sobolev:1934]). These elements were introduced not only for elliptic systems but also for hyperbolic systems by allowing the wave vectors to be the complex solutions of the dispersion relation associated with the initial system of equations. The starting point is to make an algebraization, according to [@Burnat:1972; @Grundland:1974; @Perad:1985], of a first-order system of PDEs (\[eq:1\]) in $p$ independent and $q$ dependent variables written in its matrix form

\[eq:2.1\] &\^i(u)u\_i=0,i=1,…,p,\
&x=(x\^1,…,x\^p)X\^p, u=(u\^1,…,u\^q)U\^q,

where $\mathcal{A}^i(u)={\left( \mathcal{A}^{\mu i}_\alpha(u) \right)}$ are $m\times q$ matrix functions of $u$ and we denote the partial derivatives by $u^\alpha_i={\partial}u^\alpha/{\partial}x^i$. The matrix $L_i^\alpha$ satisfying the conditions[@Burnat:1972]

\[eq:intelem\] u\_i\^[{ L\_i\^:\_\^[i]{}L\_i\^=0, =1,…,q }]{}

at some open given point $u_0\in U$ is called an integral element of the system (\[eq:2.1\]). This matrix $L={\left( {{\partial}u^\alpha}/{{\partial}x^i} \right)}$ is a matrix of the tangent mapping $du:X{\rightarrow}T_uU$ given by the formula $$X\ni (\delta x^i){\rightarrow}(\delta u^\alpha)\in T_uU,\quad \text{where }\delta u^\alpha=u_i^\alpha \delta x^i.$$The tangent mapping $du(x)$ determines an element of linear space $L(X,T_uU)$, which can be identified with the tensor product $T_uU\otimes X^\ast$, (where $X^\ast$ is the dual space of $X$, [*i.e.* ]{}the space of linear forms). It is well known [@Burnat:1972; @Grundland:1974; @Perad:1985] that each element of this tensor product can be represented as a finite sum of simple tensors of the form $$L=\gamma\otimes \lambda,$$where $\lambda\in X^\ast$ is a covector and $\gamma=\gamma^\alpha\frac{{\partial}}{{\partial}u^\alpha}\in
T_uU$ is a tangent vector at the point $u\in U$. Hence, the integral element $L_i^\alpha$ is called a simple element if ${\operatorname{rank}}(L_i^\alpha)=1.$ To determine a simple integral element $L_i^\alpha$ we have to find a vector field $\gamma\in T_uU$ and a covector $\lambda\in X^\ast$ satisfying the so-called wave relation

\[eq:2.2\] [( \_i\^[i]{}\_(u) )]{}\^=0,=1,…,m.

The necessary and sufficient condition for the existence of a nonzero solution $\gamma$ for the equation (\[eq:2.2\]) is

\[eq:2.3\] [( \_i\^i(u) )]{}&lt;(m,q).

This relation is known as the dispersion relation. If the covector $\lambda=\lambda_i(u)dx^i$ satisfies the dispersion relation (\[eq:2.3\]) then there exists a polarization vector $\gamma\in
T_uU$ satisfying the wave relation (\[eq:2.2\]). The algebraic approach which has been used in [@Burnat:1972; @Grundland:1974; @Perad:1985] for hyperbolic systems of equations (\[eq:2.1\]) allows the construction of certain classes of $k$-wave solutions admitting $k$ arbitrary functions of one variable. The replacement of the matrix of derivatives $u_i^\alpha$ in the system of equations (\[eq:2.1\]) by the simple real element $L_i^\alpha$ allows us to construct more general classes of solutions by replacing the real elements with complex ones. A specific form of solution $u(x)$ of an elliptic system (\[eq:2.1\]) is postulated for which the tangent mapping $du(x)$ is a sum of a complex element and its complex conjugate

\[eq:2.4\] du\^(x)=(x)\^(u)\_i(u)dx\^i+|(x)|\^(u)|\_i(u)dx\^i,=1,…,q,

where $\gamma=(\gamma^1,\ldots, \gamma^q)\in{\mathbb{C}}^q$ and $\lambda=(\lambda_1,\ldots,\lambda_p)\in {\mathbb{C}}^p$ satisfy

\[eq:star\] \_i\_\^[i]{}(u)\^=0,|\_i\_\^[i]{}(u)|\^=0.

Here the quantity $\xi(x)\neq 0$ is treated as a complex function of the real variables $x$. In what follows we assume that the vectors $\gamma$ and $\bar{\gamma}$ are linearly independent. The proposed form of solution (\[eq:2.4\]) is more general than the one proposed in [@Jeffrey:1976; @Rozdestvenski:1983] for which the derivatives $u^\alpha_i$ are represented by a real simple element leading to a simple Riemann wave solution. To distinguish this situation from the one proposed in (\[eq:2.4\]), we call the real-valued solution associated with a complex element and its complex conjugate a simple mode solution in accordance with [@Perad:1985]. This means that all first-order derivatives of $u^\alpha$ with respect to $x^i$ are decomposable in the following way

\[eq:star2\] =(x)\^(u)\_i(u)+|(x)|\^(u)|\_i(u),

or alternatively $$\frac{{\partial}u^\alpha}{{\partial}x^i}=i{\left( \xi(x) \gamma^\alpha(u)\lambda_i(u)-\bar{\xi}(x)\bar{\gamma}(u)^\alpha\bar{\lambda}_i(u) \right)},$$where a set of functions $(\lambda,\gamma)$ and their conjugates $(\bar{\lambda},\bar{\gamma})$ on $U$ satisfy the wave relation (\[eq:star\]). Similarly, as in the case of $k$-waves, for hyperbolic systems [@Perad:1985], we have to find the necessary and sufficient conditions for the existence of solutions of type (\[eq:2.4\]). These conditions are called involutivity conditions.

#### 

First we derive a number of necessary conditions on the vector fields $(\gamma,\lambda)$ and their complex conjugate $(\bar{\gamma},\bar{\lambda})$ as a requirement for the existence of rank-$2$ (simple mode) solutions of the homogeneous system (\[eq:2.1\]). Namely, closing (\[eq:2.4\]) by exterior differentiation, we obtain the following 2-forms

\[eq:2.5\] (d+d)+d+|(d||+|d|)+|d||=0,

which have to satisfy (\[eq:2.4\]). Using (\[eq:2.4\]) we get

\[eq:2.6\] d=||\_[,|]{}+|\_[,]{},d=|\_[,|]{}|+\_[,]{}d,

where we have used the following notation $$\lambda_{,\gamma}=\gamma^\alpha\frac{{\partial}}{{\partial}u^\alpha}\lambda,\quad \gamma_{,\bar{\gamma}}=\bar{\gamma}^\alpha\frac{{\partial}}{{\partial}u^\alpha}\gamma.$$Substituting (\[eq:2.6\]) into the prolonged system (\[eq:2.5\]) we obtain

\[eq:2.8\] +|+||=0,

whenever the differential (\[eq:2.4\]) holds. The commutator of vector fields $\gamma$ and $\bar{\gamma}$, is denoted by $${\left[ \gamma,\bar{\gamma} \right]}={\left( \gamma,\bar{\gamma} \right)}_u+\gamma_{,\lambda_i}\lambda_{i,\bar{\gamma}}-\bar{\gamma}_{,\lambda_i}\lambda_{i,\gamma},$$while by ${\left( \gamma,\bar{\gamma} \right)}_u$ we denote a part of the commutator which contains the differentiation with respect to the variables $u^\alpha$, [*i.e.* ]{}$${\left( \gamma,\bar{\gamma} \right)}_u=\bar{\gamma}^\alpha\frac{{\partial}}{{\partial}u^\alpha}\left.\gamma(u,\lambda)\right|_{\lambda=const.}
-\gamma^\alpha\frac{{\partial}}{{\partial}u^\alpha}\left.\bar{\gamma}(u,\bar{\lambda})\right|_{\bar{\lambda}=const.}.$$Let $\Phi$ be an annihilator of the vectors $\gamma$ and $\bar{\gamma}$, [*i.e.* ]{}$$<\omega \lrcorner\gamma>=0, \qquad <\omega \lrcorner \bar{\gamma}>=0, \quad \omega\in\Phi=\operatorname{An}{\left\{ \gamma,\bar{\gamma} \right\}}.$$Here, by the parenthesis $<\omega\lrcorner \gamma>$, we denote the contraction of the 1-form $\omega\in T_u^\ast U$ with the vector field $\gamma\in T_uU$ and similarly for the vector field $\bar{\gamma}\in T_uU$. Multiplying the equations (\[eq:2.8\]) by the 1-form $\omega\in\Phi$ we get

\[eq:2.9\] |&lt;&gt;|=0.

We look for the compatibility condition for which the system (\[eq:2.9\]) does not provide any algebraic constraints on the coefficients $\xi$ and $\bar{\xi}$. This postulate means that the profile of the simple modes associated with $\gamma\otimes \lambda$ and $\bar{\gamma}\otimes
\bar{\lambda}$ can be chosen in an arbitrary way for the initial (or boundary) conditions. It follows from (\[eq:2.9\]) that the commutator of the vector fields $\gamma$ and $\bar{\gamma}$ is a linear combination of $\gamma$ and $\bar{\gamma}$, where the coefficients are not necessarily constant, [*i.e.* ]{}$${\left[ \gamma,\bar{\gamma} \right]}\in\operatorname{span}{\left\{ \gamma,\bar{\gamma} \right\}}.$$So, the Frobenius theorem is satisfied. That is, at every point $u_0$ of the space of dependent variables $U$, there exists a tangent surface $S$ spanned by the vector fields $\gamma$ and $\bar{\gamma}$ passing through the point $u_0\in U$. Moreover the above condition implies that there exists a complex-valued function $\alpha$ such that

\[eq:2.10\] =-||,

since $\overline{{\left[ \gamma,\bar{\gamma} \right]}}={\left[ \bar{\gamma},\gamma \right]}=-{\left[ \gamma,\bar{\gamma} \right]}$ holds. Next, from the vector fields $\gamma$ and $\bar{\gamma}$ defined on $U$-space, we can construct the coframe $\Psi$, that is the set of 1-forms $\sigma_1$ and $\sigma_2$ defined on $U$ satisfying the conditions

\[eq:2.11\] &&lt;\_1&gt;=1,&&&lt;\_1|&gt;=0,\_1,\_2,\
&&lt;\_2&gt;=0, &&&lt;\_2|&gt;=1.

Substituting (\[eq:2.6\]) and (\[eq:2.10\]) into the prolonged system (\[eq:2.8\]) and multiplying by $\sigma_1$ and $\sigma_2$ respectively, we get

\[eq:2.11b\] &(i)&& d+[( ||\_[,|]{}+|\_[,]{} )]{}+||=0,\
&(ii) &&d||+|[( |\_[,]{}+||\_[,|]{} )]{}+|||=0.

If $d\xi=d\bar{\xi}$ then equation (\[eq:2.11b\].i) is a complex conjugate of (\[eq:2.11b\].ii). So one can consider one of them, say (\[eq:2.11b\].i). We look for the condition of integrability such that the system (\[eq:2.11b\].i) does not impose any restriction on the form of the coefficients $\xi$ and $\bar{\xi}$. By exterior multiplication of (\[eq:2.11b\].i) by $\lambda$, we obtain

\[eq:2.12i\] \_[,|]{}|=0,

and its respective complex conjugate equation

\[eq:2.12ii\] |\_[,]{}|=0.

Note that if conditions (\[eq:2.10\]), (\[eq:2.12i\]) and (\[eq:2.12ii\]) are satisfied, then from the Cartan lemma[@Cartan:1953] the system (\[eq:2.11b\]) and consequently the set of 1-forms (\[eq:2.8\]) has nontrivial solutions for $d\xi$ and $d\bar{\xi}$. So under these circumstances the following result holds [@Perad:1985].

\[prop:1\] The necessary condition for the Pfaffian system (\[eq:2.4\]) to possess a simple mode solution is that the following two constraints on the complex-valued vector fields $\gamma$ and $\lambda$ be satisfied $$(i)\ \quad {\left[ \gamma,\bar{\gamma} \right]}=\alpha \gamma -\bar{\alpha} \bar{\gamma},$$for any complex-valued function $\alpha\in{\mathbb{C}}$, and $$(ii)\quad \lambda_{,\bar{\gamma}}, \bar{\gamma}_{,\gamma}\in\operatorname{span}{\left\{ \lambda,\bar{\lambda} \right\}}.$$

#### 

Note that the conditions (\[eq:2.10\]), (\[eq:2.12i\]) and (\[eq:2.12ii\]) are strong requirements on the functions $\gamma={\left( \gamma^1,\ldots,\gamma^q \right)}$ and $\lambda=(\lambda_1,\ldots,\lambda_p)$ and their complex conjugates which satisfy the wave relation (\[eq:2.2\]). Consequently, the postulated form of a one-mode solution $u(x)$ of the elliptic system (\[eq:2.1\]), required by the generalized method of characteristics (GMC), is such that all first-order derivatives of $u(x)$ with respect to $x^i$ are decomposable in the form (\[eq:2.4\]). This restriction is a strong limitation on the admissible class of solutions of (\[eq:2.1\]). So far, there have been few known examples of such solutions. Therefore, it is worth developing this idea by weakening the integrability conditions (\[eq:2.10\]), (\[eq:2.12i\]) and (\[eq:2.12ii\]) in order to construct multi-mode solutions (considered to be nonlinear superpositions of simple modes). In the next section we propose an alternative way of constructing multi-mode solutions (expressed in terms of Riemann invariants) which are obtained from a version of the conditional symmetry method (CSM) [@GrundHuard:2007] by adapting it to the elliptic systems. In Section \[sec:5\] and \[sec:6\] a further refinement of the above technique will be presented, which consist of generalising the idea of mode solutions for quasilinear systems (\[eq:2.1\]). This idea is based on the specific factorization of integral elements (\[eq:intelem\]) through the introduction of some rotation matrices by weakening the restrictions impose by the wave relations (\[eq:2.2\]). This approach allows us to go deeper into the geometrical aspects of solving systems (\[eq:2.1\]) and enables us to obtain new results. This constitutes the objective of this paper.

Conditional symmetry method and multimode solutions in terms of Riemann Invariants. {#sec:3}
===================================================================================

\[sec:condsym\] In this section, we examine certain aspects of the conditional symmetry method in the context of Riemann invariants for which the wave vectors $\lambda$ and $\bar{\lambda}$ are complex solutions of the dispersion relation (\[eq:2.3\]) associated with the original system (\[eq:2.1\]). Let us consider a nondegenerate first-order quasilinear elliptic system of PDEs (\[eq:2.1\]) in its matrix form

\[eq:3.1\] \^1(u)u\_1+…+\^p(u)u\_p=0,

where $\mathcal{A}^1,\ldots,\mathcal{A}^p$ are $m\times q$ real-valued matrix functions of $u$. Let us suppose that there exist $k$ linearly independent non-conjugate complex-valued wave vectors

\[eq:wv\] \^(u)=[( \_1\^A(u),…,\_p\^A(u) )]{}\^p,\^A|\^BA,B[{ 1,…,k&lt;p }]{},

which satisfy the dispersion relation (\[eq:2.3\]). One should note that in (\[eq:wv\]) we do not require indices $A\neq B$, which means that real wave vectors are excluded from our consideration (we do not consider here the mixed case for wave vectors involving real and complex wave vectors). Under the above hypotheses, the $k$ wave vectors (\[eq:wv\]) and their complex conjugates $$\bar{{\lambda}}^A(u)={\left( \bar{{\lambda}}_1^A(u),\ldots,\bar{{\lambda}}_p^A(u) \right)}\in{\mathbb{C}}^p,\quad A=1,\ldots,k,$$satisfy the dispersion relation (\[eq:2.3\]). In what follows it is useful to introduce the notation $c.c.$ which means the complex conjugate of the previous term or equation. This notation is convenient for computational purposes allowing the presentation of some expressions in abbreviated form.

#### 

Let us suppose that there exists a unique solution $u(x)$ of the system (\[eq:3.1\]) of the form

\[eq:3.5\] u=f(r\^1(x,u),…, r\^k(x,u),|[r]{}\^1(x,u),…,|[r]{}\^k(x,u))+c.c.,

where the complex-valued functions $r^A,\bar{r}^A:{\mathbb{R}}^p\times {\mathbb{R}}^q{\rightarrow}{\mathbb{C}}$ are called the Riemann invariants associated respectively to wave vectors $\lambda^A$, $\bar{\lambda}^A$ and are defined by

\[eq:3.2\] r\^A(x,u)=\_i\^A(u)x\^i,|[r]{}\^A(x,u)=|\_i\^A(u)x\^i,A=1,…,k,

where $c.c.$ means complex conjugated previous term. Note that the functions $u(x)$ are defined implicitly in terms of $u^\alpha, x^i$, $r^A$ and $\bar{r}^A$. For any function $f:{\mathbb{C}}^k{\rightarrow}{\mathbb{C}}^q$ and its complex conjugate, the equation (\[eq:3.5\]) determines a unique real-valued function $u(x)$ on a neighborhood of the origin $x=0$. Note also that an analogue analysis as presented below can be performed if one replace the postulated form of the solution (\[eq:3.5\]) written in the Riemann invariants by the expression $$u=i{\left( f(r^1(x,u),\ldots,r^k(x,u),\bar{r}^1(x,u),\ldots,\bar{r}^k(x,u))-c.c. \right)}.$$Therefore we omit this case.

#### 

The Jacobi matrix of derivatives of $u(x)$ is given by

\[eq:3.6\] u=(u\^\_i)=[( I\_q- )]{}\^[-1]{} [( [( + )]{}+c.c. )]{},

or equivalently as

\[eq:3.7\] u=[( + )]{}[( M\^1+M\^2| )]{}+c.c.,

where c.c. means the complex conjugate of the previous term. The $k\times k$ matrices $M^1$ and $M^2$ are defined by

\[eq:3.7b\] M\^1=&\^[-1]{},\
M\^2&=-M\^1 [( + )]{}[( I\_k- )]{}\^[-1]{},

\[eq:3.8\] =[( )]{}\^[qk]{},=[( )]{}\^[qk]{},=(\_i\^A)\^[kp]{},

\[eq:3.9\] =[( )]{}=[( x\^i )]{}\^[kq]{},r=(r\^1,…,r\^k),

and their respective conjugate equations. The matrices $I_q$ and $I_k$ are the $q\times q$ and $k\times k$ identity matrices respectively. We use the implicit function theorem to obtain the following conditions ensuring that functions $r^A,\bar{r}^A$ and $u^\alpha$ are expressible as graphs over some open subset $\mathcal{D}\subset {\mathbb{R}}^p$,

\[eq:3.10\] [( I\_q- )]{}0

or

\[eq:3.11\] &[( I\_k- )]{}0,0.

So the inverse matrix in (\[eq:3.6\]) or in (\[eq:3.7\]) is well-defined in the vicinity of $x=0$, since

\[eq:3.12\] =0,=0 x=0.

Note that on the hypersurface defined by equation (\[eq:3.5\]) and one of the expressions (\[eq:3.10\]) or (\[eq:3.11\]) equal to zero, the gradient of the function $u(x)$ becomes infinite for some value of $x^i$. So the multimode solution, given by expression (\[eq:3.5\]) loses its sense on this hypersurface. Consequently, some types of discontinuities, [*i.e.* ]{}shock waves can occur. In what follows, we search for solutions defined on a neighborhood of $x=0$ under the assumption that the conditions (\[eq:3.10\]) or (\[eq:3.11\]) are different from zero. This means that if the initial data is sufficiently small, then there exists a time interval ${\left[ t_0,T \right]}$, $T>t_0$ (where we denote the independent variables by $t=x^0, \tilde{x}=x^1,\dots,
x^n$ where $p=n+1$) in which the gradient catastrophe for the solution $u(t,\tilde{x})$ of the system (\[eq:2.1\]) does not take place [@GrundlandVassiliou:1991; @Rozdestvenski:1983].

#### 

Note that the Jacobian matrix of $u(x)$ has at most rank equal to $2k$. It follows that the proposed solution (\[eq:3.5\]) is also at most of rank-$2k$. Its image is a $2k$-dimensional submanifold ${\mathbb{S}}_{2k}$ in the first jet space $J^1=J^1(X\times U)$.

#### 

Let us introduce a set of $p-2k$ linearly independent vectors $\xi_a:{\mathbb{R}}^q{\rightarrow}{\mathbb{C}}^p$ defined by

\[eq:3.13\] \_a(u)=[( \_a\^1(u),…,\_a\^p(u) )]{},a=1,…,p-2k,

satisfying the orthogonality conditions

\[eq:3.14\] \_i\^A\_a\^i=0,|\_i\^A\_a\^i=0,A=1,…,k,a=1,…,p-2k,

for a set of $2k$ linearly independent wave vectors ${\left\{ \lambda^1,\ldots,\lambda^k,
\bar{\lambda}^1,\ldots, \bar{\lambda}^k \right\}}$. It should be noted that the set ${\left\{ \lambda^1,\ldots,\lambda^k,\bar{\lambda^1},\ldots,\bar{\lambda}^k,\xi^1,\ldots,\xi^{p-2k} \right\}}$ forms a basis for the space of independent variables $X$. Note also, that the vectors $\xi_a$ are not uniquely defined since they obey the homogeneous conditions (\[eq:3.14\]). As a consequence of equation (\[eq:3.6\]) or (\[eq:3.7\]), the graph $\Gamma={\left\{ x,u(x) \right\}}$ is invariant under the family of first-order differential operators

\[eq:3.15i\] X\_a=\_a\^i(u),a=1,…,p-2 k,

defined on $X\times U$ space. Since the vector fields $X_a$ do not include vectors tangent to the direction of $u$, they form an Abelian distribution on $X\times U$ space, [*i.e.* ]{}

\[eq:3.16\] =0,a, b=1,…,p-2k.

The set ${\left\{ r^1,\ldots, r^k,\bar{r}^1,\ldots,\bar{r}^k,u^1,\ldots
u^q \right\}}$ constitutes a complete set of invariants of the Abelian algebra $\mathcal{L}$ generated by the vector fields (\[eq:3.15i\]). So geometrically, the characterization of the proposed solution (\[eq:3.5\]) of equations (\[eq:3.1\]) can be interpreted in the following way. If $u(x)$ is a $q$-component function defined on a neighborhood of the origin $x=0$ such that the graph of the solution $\Gamma={\left\{ (x,u(x)) \right\}}$ is invariant under a set of $p-2k$ vector fields $X_a$ with the orthogonality property (\[eq:3.14\]), then for some function $f$ the expression $u(x)$ is a solution of equation (\[eq:3.5\]). Hence the group-invariant solutions of the system (\[eq:3.1\]) consist of those functions $u=f(x)$ which satisfy the overdetermined system composed of the initial system (\[eq:3.1\]) together with the invariance conditions

\[eq:3.17i\] \_a\^iu\_i\^=0,i=1,…,p,a=1,…,p-2k,

ensuring that the characteristics of the vector fields $X_a$ are equal to zero.

#### 

It should be noted that, in general, the conditions (\[eq:3.17i\]) are weaker than the differential constraints (\[eq:2.4\]) required by the generalized method of characteristics, since the latter method is subjected to the algebraic conditions (\[eq:2.2\]). In fact, equations (\[eq:3.17i\]) imply that there exist complex-valued matrix functions $\Phi_A^\alpha(x,u)$ and $\overline{\Phi}_A^\alpha(x,u)$ defined on the first jet space $J=J(X\times U)$ such that all first derivatives of $u$ with respect to $x^i$ are decomposable in the following way

\[eq:3.18\] u\_i\^=\_A\^(x,u)\_i\^A+\_A\^(x,u)|\_i\^A,

where

\[eq:3.19\] \_A\^&=[( I\_q- )]{}\^[-1]{} [( + )]{},\
\_A\^&=[( I\_q- )]{}\^[-1]{} [( + )]{},

or

\[eq:3.20\] \_A\^&=[( + )]{}M\^1+[( + )]{}\^2,\
\_A\^&=[( + )]{}\^1+[( + )]{}M\^2.

The matrices $\Phi_A^\alpha\lambda_i^A$ and $\overline{\Phi}_A^\alpha\bar{\lambda}_i^A$ appearing in equation (\[eq:3.18\]) do not necessarily satisfy the wave relation (\[eq:2.2\]). As a result of this fact, the restrictions on the initial data at $t=0$ are eased, so we are able to consider more diverse types of modes in the superpositions than in the case of the generalized method of characteristics described in Section \[sec:2\].

#### 

Let us now proceed to solve the overdetermined system composed of equations (\[eq:3.1\]) and differential constraints (\[eq:3.17i\])

\[eq:3.21\] \^[i]{}\_(u)u\_i\^=0,\_a\^i(u)u\_i\^=0.

Substituting (\[eq:3.6\]) or (\[eq:3.7\]) into (\[eq:3.1\]) yields the trace condition

\[eq:3.22\] [( \^\^[-1]{} [( [( + )]{}+c.c. )]{} )]{}=0,

or

\[eq:3.23\] [( \^ )]{}=0,

on the wave vectors $\lambda$ and $\bar{\lambda}$ and on the functions $f$ and $\bar{f}$, where $\mathcal{A}^1,\ldots,\mathcal{A}^q$ are $p\times q$ matrix functions of $u$ ([*i.e.* ]{}$\mathcal{A}^\mu=(\mathcal{A}^{\mu i}_\alpha(u))\in
{\mathbb{R}}^{p\times q}, \mu=1,\ldots, m$). For the given initial system of equations (\[eq:3.1\]), the matrices $\mathcal{A}^\mu$ are known functions of $u$ and the trace conditions (\[eq:3.22\]) or (\[eq:3.23\]) are conditions on the functions $f$, $\bar{f}$, $\lambda$, $\bar{\lambda}$ (or on $\xi$ due to the orthogonality conditions (\[eq:3.14\])). From the computational point of view, it is useful to split $x^i$ into $x^{i_A}$ and $x^{i_a}$ and to choose a basis for the wave vector $\lambda^A$ and $\bar{\lambda}^A$ such that

\[eq:3.24\] \^A=dx\^[i\_A]{}+\_[i\_a]{}\^Adx\^[i\_a]{},|\^A=dx\^[i\_A]{}+|\_[i\_a]{}\^Adx\^[i\_a]{}, A=1,…, k,

where $(i_A,i_a)$ is a permutation of $(1,\ldots, p)$. So, the expressions (\[eq:3.9\]) become

\[eq:3.25\] =x\^[i\_a]{},=x\^[i\_a]{}.

Substituting (\[eq:3.25\]) into (\[eq:3.22\]) (or (\[eq:3.23\])) yields

\[eq:3.26\] [( \^\^[-1]{}[( + )]{})]{}=0,=1…,m,

or

\[eq:3.27\] [( \^\^[-1]{})]{}=0,=1…,m,

where $R=(r^1,\ldots, r^k,\bar{r}^1,\ldots,\bar{r}^k)^T$ and

\[eq:3.28\] Q\_a=[( + )]{}+[( + )]{}=[( + )]{} \_a\^[qq]{},\
K\_a=

(

[cc]{} &\
&

)

=\_a[( + )]{}\^[2k2k]{},

and for simplicity of notation, we note

\[eq:3.28b\] =

(

&\
&|

)

\^[2kp]{},\_a=

(

&\
&

)

\^[2kq]{},=[( , )]{}\^[q2k]{},

for $i_A$ fixed and $i_a=1,\ldots, p-1$. In (\[eq:3.28\]) the $2k\times 2k$ matrix $K_a$ is defined in terms of the $k\times k$ subblocks of the form ${{\partial}\lambda_{i_a}}/{{\partial}u}{\left( {{\partial}f}/{{\partial}r}+{{\partial}\bar{f}}/{{\partial}r} \right)}$, where $\eta_a$ is a matrix form of the block ${\partial}\lambda_{a}/{\partial}u$ over the block ${\partial}\bar{\lambda}_a/{\partial}u$. The notation ${\left( {{\partial}f}/{{\partial}r}, {{\partial}f}/{{\partial}\bar{r}} \right)}$ represents the matrix formed of the left block ${\partial}f/{\partial}r$ and the right block ${\partial}f/{\partial}\bar{r}$. Note that the functions $r^A$, $\bar{r}^A$ and $x^{i_a}$ are functionally independent in a neighborhood of $x=0$ and the matrix functions $\mathcal{A}^\mu$, ${\partial}f/{\partial}r$, ${\partial}f/{\partial}\bar{r}$, ${\partial}\bar{f}/{\partial}r$, ${\partial}\bar{f}/{\partial}\bar{r}$, $Q_a$ and $K_a$ depend on $r$ and $\bar{r}$ only. So, equations (\[eq:3.26\]) (or (\[eq:3.27\])) have to be satisfied for any value of coordinates $x^{i_a}$. As a consequence, we have some constraints on these matrix functions. From the Cayley-Hamilton theorem, we know that for any $n\times n$ invertible matrix $M$, the expression $(M^{-1}\det M)$ is a polynomial in $M$ of order $(n-1)$. Thus, using the tracelessness of the expression $\mathcal{A}^\mu{\left( I_q-Q_ax^{i_a} \right)}^{-1}({\partial}f/{\partial}R)\Lambda$, we can replace equations (\[eq:3.26\]) by the following condition

\[eq:3.29\] [( \^Q[( + )]{})]{}=0,Q=(I\_q-Q\_ax\^[i\_a]{})\^[qq]{},

where $\operatorname{adj}M$ denotes the adjoint of the matrix $M$. As a consequence the matrix $Q$ is a polynomial of order $(q-1)$ in $x^{i_a}$. Taking (\[eq:3.29\]) and all its partial derivatives with respect to $x^{i_a}$ (with $r$, $\bar{r}$ fixed at $x=0$), we obtain the following conditions for the matrix functions $f(r,\bar{r})$ and $\lambda(f(r,\bar{r}))$

\[eq:3.30\] [( \^)]{}=0,=1,…,m,

\[eq:3.31\] [( \^Q\_[(a\_1.]{}…Q\_[a\_s)]{}[( + )]{})]{}=0,

where $s=1,\ldots,q-1$ and $(a_1,\ldots, a_s)$ denotes the symmetrization over all indices in the bracket. A similar procedure can be applied to system (\[eq:3.27\]) to yield (\[eq:3.30\]) and

\[eq:3.32\] [( \^K\_[(a\_1]{}…K\_[a\_s)]{})]{}=0,

where now $s=1,\ldots,2k-1$. Equations (\[eq:3.30\]) represent the initial value conditions on a surface in the space of independent variables $X$, given at $x^{i_a}=0$. Note that equations (\[eq:3.31\]) (or (\[eq:3.32\])) form the conditions required for the preservation of the property (\[eq:3.30\]) along the flows represented by the vector fields (\[eq:3.15i\]). Equation (\[eq:3.24\]) allows us to express $X_a$ in the form

\[eq:3.33\] X\_a=\_[i\_a]{}-\_[i\_a]{}\^A\_[i\_A]{}-|\_[i\_a]{}\^A\_[i\_A]{},A=1,…,k.

Substituting expressions (\[eq:3.28\]) into (\[eq:3.31\]) or (\[eq:3.32\]) and simplifying gives the unified form

\[eq:3.34\] [( \^[( + )]{}\_[(a\_1.]{}[( + )]{}…\_[[a\_s]{})]{}[( + )]{})]{}=0.

The index $s$ is either $\max(s)=q-1$ or $\max(s)=2k-1$, we choose the one which is more convenient from the computational point of view. In this case, for $k\geq 1$ the two approaches, CSM and GMC, become essentially different and, as we demonstrate in the following example, the CSM can provide rank-$2k$ solutions which are not Riemann $2k$-waves as defined by the GMC since we weaken the integrability conditions (\[eq:2.1\]) for the wave vectors $\lambda$ and $\bar{\lambda}$.

#### 

A change of variable on $X\times U$ allows us to rectify the vector fields $X_a$ and considerably simplify the structure of the overdetermined system (\[eq:3.21\]) which classifies the preceding construction of multimodes solutions. For this system, in the new coordinates, we derive the necessary and sufficient conditions for the existence of rank-$2k$ solutions of the form (\[eq:3.5\]). Suppose that there exists an invertible $2k\times 2k$ subblock matrix

\[eq:3.35\] H=(\^s\_t), 1s,t2k,

of the larger matrix $\Lambda\in {\mathbb{C}}^{2k\times p}$, then the independent vector fields $X_a$ can be written as

\[eq:3.36\] X\_a=\_[x\^[a+2k]{}]{}-(H\^[-1]{})\^s\_t\^t\_[a+2k]{}\_[x\^A]{},

which have the required form (\[eq:3.15i\]) and for which the orthogonality conditions (\[eq:3.14\]) are fulfilled. We introduce new coordinate functions

\[eq:3.37\] z\^1=r\^1(x,u),…, z\^[k]{}=r\^[k]{}(x,u),z\^[k+1]{}=|[r]{}\^1(x,u),…, z\^[2k]{}=|[r]{}\^[k]{}(x,u),\
z\^[2k+1]{}=x\^[2k+1]{},…,z\^p=x\^p,v\^1=u\^1,…, v\^q=u\^q,

on $X\times U$ space which allow us to rectify the vector fields (\[eq:3.36\]). As a result, we get

\[eq:3.38\] X\_1=,…, X\_[p-2k]{}=.

The $p$-dimensional submanifold invariant under $X_{1},\ldots,X_{p-2k}$, is defined by equations of the form

\[eq:3.39\] v=f(z\^1,…, z\^k, |[z]{}\^1,…, |[z]{}\^k)+c.c.

for an arbitrary function $f:X{\rightarrow}U$. The expression (\[eq:3.39\]) is the general solution of the invariance conditions

\[eq:3.40\] v\_[z\^[2k+1]{}]{},…, v\_[z\^p]{}=0.

In general, the initial system (\[eq:3.1\]) described in the new coordinates $(x,v)\in X\times U$ is a nonlinear system of first-order PDEs, of the form

\[eq:3.42\]

&\^[l]{}\_(v)=0,&&i=1,…, k,\
&\^[l]{}\_(v)=0,&&i=k+1,…, 2k,\
&\^[i]{}\_(v)=0, &&i=2k+1,…,p.

We obtain the following Jacobi matrix in the coordinates $(z,\bar{z},v)$

\[eq:3.43\] =[( \^[-1]{} )]{}\^j\_l\^l\_i\^[pp]{},=\_s\^i-x\^l,

whenever the invariance conditions (\[eq:3.40\]) are satisfied. Appending to the system (\[eq:3.42\]) the invariance condition (\[eq:3.40\]), we obtain the quasilinear reduced system of PDEs

\[eq:3.44\] [( \^(v)[( I\_[q]{}- )]{}\^[-1]{})]{}=0,,…,=0,=1,…,m,

or

\[eq:3.45\] [( \^(v)\^[-1]{})]{}=0,,…,=0,=1,…,m.

#### 

We now provide some basic definitions that we require for the use of the conditional symmetry method in order to encompass the use of Riemann invariants.

#### 

A vector field $X_a$ is called a conditional symmetry of the original system (\[eq:3.1\]) if $X_a$ is tangent to the manifold ${\mathbb{S}}={\mathbb{S}}_\Delta\cap {\mathbb{S}}_Q$, [*i.e.* ]{}

\[eq:3.46\] .\^[(1)]{}X\_a|\_T\_[(x,u\^[(1)]{})]{},

where the first prolongation of $X_a$ is given by

\[eq:3.47\] \^[(1)]{}X\_a=X\_a-\^i\_[a,u\^]{}u\_j\^u\_i\^,a=1,…,p-2k

and the submanifolds of the solution spaces are given by

\[eq:3.48\] \_=[{ (x,u\^[(1)]{}):\_\^[i]{}u\_i\^=0,=1,…,m }]{},

and

\[eq:3.49\] \_Q=[{ (x,u\^[(1)]{}):\_a\^i(u)u\_i\^=0,=1,…, q, a=1,…,p-2k }]{}.

Consequently, an Abelian Lie algebra $\mathcal{L}$ generated by the vector fields $X_1,\ldots, X_{p-2k}$ is called a conditional symmetry algebra of the original system (\[eq:3.1\]) if the conditions

\[eq:3.50\] .\^[(1)]{}X\_a[( \^iu\_i )]{}|\_=0,a=1,…, p-2k,

are satisfied.

#### 

Supposing that $\mathcal{L}$, spanned by the vector fields $X_1,\ldots,
X_{p-2k}$, is a conditional symmetry algebra of the system (\[eq:3.1\]), a solution $u=f(x)$ is said to be a conditionally invariant solution of the system (\[eq:3.1\]) if the graph $\Gamma={\left\{ (x,f(x)) \right\}}$ is invariant under the vector fields $X_1,\ldots, X_{p-2k}$.

#### **Proposition** 2.

*A nondegenerate quasilinear hyperbolic system of first-order PDEs (\[eq:3.1\]) in $p$ independent variables and $q$ dependent variables admits a $p-2k$-dimensional symmetry algebra $\mathcal{L}$ if and only if ($p-2k$) linearly independent vector fields $X_1,\ldots, X_{p-2k}$ satisfy the conditions (\[eq:3.30\]) and (\[eq:3.34\]) on some neighborhood of $(x_0,u_0)$ of $S$. The solutions of (\[eq:3.1\]) which are invariant under the Lie algebra $\mathcal{L}$ are precisely $k$-mode solutions of the form (\[eq:3.5\]).*

#### Proof.

The proof of this proposition is essentially similar to that of the proposition in [@GrundHuard:2007]. For the sake of completeness we now give it. Let us express the vector fields $X_a$ in the new coordinates $(z,v)$ on $X\times U$. Equations (\[eq:3.38\]) and (\[eq:3.47\]) imply that

\[eq:3.60\] \^[(1)]{}X\_a=X\_a,a=1,…, p-2k.

The symmetry criterion that has to be satisfied for $G$ to be the symmetry group of the overdetermined system (\[eq:3.44\]) (or (\[eq:3.45\])) requires that the vector fields $X_a$ of $G$ satisfy $X_a(\Delta)=0$, whenever equation (\[eq:3.44\]) (or (\[eq:3.45\])) is satisfied. Thus the symmetry criterion applied to the invariance conditions (\[eq:3.40\]) vanishes identically. Applying this criterion to the system (\[eq:3.42\]) in the new coordinates, carrying out the differentiation and taking into account the conditions (\[eq:3.30\]) and (\[eq:3.34\]), we obtain the equations which are identically satisfied.

#### 

The converse is also true. The assumption that the system (\[eq:3.1\]) is nondegenerate means, according to [@Olver:Application_of_Lie], that it is locally solvable and is of maximal rank at every point $(x_0,u_0)\in S$. Therefore, the infinitesimal symmetry condition is a necessary and sufficient condition for the existence of the symmetry group $G$ of the overdetermined system (\[eq:3.21\]). Since the vector fields $X_a$ form an Abelian distribution on $X\times U$, it follows that conditions (\[eq:3.30\]) and (\[eq:3.34\]) are satisfied. The solutions of the overdetermined system (\[eq:3.21\]) are invariant under the algebra $\mathcal{L}$ generated by the $p-2k$ vector fields $X_1,\ldots,X_{p-2k}$. The invariants of the group $G$ of such vector fields are provided by the functions ${\left\{ r^1,\ldots, r^k,\bar{r}^1,\ldots, \bar{r}^k,u^1,\ldots,u^q \right\}}$. So the general multimode solution of (\[eq:3.1\]) takes the required form (\[eq:3.5\]).  $\Box$

The ideal plastic flow. {#sec:4}
=======================

In this section we would like to illustrate the proposed approach for constructing multimode solutions with the example of the ideal nonstationary irrotational planar ideal plastic flow subjected to an external force due to a work function $V$ (potential if $V_t=0$). Under the above assumptions the examined model is governed by a quasilinear elliptic homogenous system of five equations in (2+1) dimensions of the form [@Chak:2006; @Hill:1998; @Kat:1]

\[eq:4.1\] &(a) &&\_x- [( \_x 2+\_y 2)]{}+[( V\_x-u\_t - u u\_x - v u\_y )]{}=0,\
&(b) &&\_y- [( \_x2- \_y 2)]{}+[( V\_y-v\_t - u v\_x - v v\_y )]{}=0,\
&(c) &&(u\_y+v\_x)2+ (u\_x-v\_y)2=0,\
&(d) && u\_x+v\_y=0,\
&(e) && u\_y-v\_x=0.

Equations (\[eq:4.1\].a) and (\[eq:4.1\].b) involve the independent variables $\sigma, \theta,
u, v$ and the potential $V$ is a given function of $(t,x,y)\in{\mathbb{R}}^3$. The stress tensor is defined by the mean pressure $\sigma$ and the angle $\theta$ relative to the $x$-axis minus $\pi/4$. Equation (\[eq:4.1\].c) represents the Saint-Venant-Von Mises plasticity equation. Equations (\[eq:4.1\].d) and (\[eq:4.1\].e) for the velocities $u$ and $v$ (along the $x$-axis and $y$-axis respectively) correspond to the incompressibility and the irrotationality of the flow of the plastic material. The independent variables are denoted by $(x^i)=(t,x,y)\in X\subset {\mathbb{R}}^3$ and the unknown functions by $(u^\alpha)=(\sigma,\theta, u,v)\in U\subset {\mathbb{R}}^4$. The first-order partial derivatives of $\sigma$ with respect to the independent variables $x$ and $y$ are denoted $\sigma_x$ and $\sigma_y$ respectively. The first-order partial derivatives of $u$, $v$ and $\theta$ are denoted similarly.

#### 

Before using the approach described in Section \[sec:condsym\] to obtain solutions of the system (\[eq:4.1\]), we begin by expressing the system (\[eq:4.1\]) in a form which is more convenient for the purpose of computation. In the first place, we note that the equations (\[eq:4.1\].a) and (\[eq:4.1\].b) governing the average pressure $\sigma$ can be solved by quadrature when the angle $\theta$ satisfies the compatibility conditions for the mixed derivatives of $\sigma$ with respect to $x$ and $y$. This equation takes the form

\[eq:4.2\] 2[( \_[x]{}\^2-\_[xy]{}-\_y\^2 )]{}(2)+[( \_[xx]{}-4\_x\_y-\_[yy]{} )]{}(2)=0.

If $\theta$ satisfies (\[eq:4.2\]), then the pressure $\sigma$ can be expressed in terms of the velocities $u$, $v$, of the angle $\theta$ and of the potential $V$ with the form

\[eq:4.3\] (t,x,y)&=-V(t,x,y)+[12]{}(2(t,x,y))++u\_t(t,x,y)dx\
& +\_y(t,x,y)(2(t,x,y))dx+c\_0(t).

It should be noted that in order to obtain the compatibility condition (\[eq:4.2\]), one has to make use of equation (\[eq:4.1\].e), that is, we must suppose that the flow is irrotational. In order to find a solution of the original system (\[eq:4.1\]), it is still necessary to find a solution of the overdetermined system consisting of equations (\[eq:4.2\]), (\[eq:4.1\].c)-(\[eq:4.1\].e) for the three unknown functions $u,v$ and $\theta$. In order to reduce the order of the PDE (\[eq:4.2\]), we introduce the dependent variables $\phi$ and $\psi$, which are defined by the equations

\[eq:4.4\] =\_x,=\_y.

Therefore, we have to solve the following well-determined first-order system of 4 equations in the 4 unknowns $\phi$, $\psi$, $u$, $v$

\[eq:4.5\] (a)&&&[( 2\^2-(\_y-\_x)-2\^2 )]{}u\_y+[( \_x+4-\_y )]{}u\_x=0,\
(b)&&&\_y-\_x=0,\
(c)&&&u\_x+v\_y=0,\
(d)&&&u\_y-v\_x=0.

Time appears in system (\[eq:4.5\]) only as a parameter since the system does not involve any derivatives with respect to time. Therefore, when solving the system (\[eq:4.5\]) we treat the unknowns $\phi, \psi, u, v$ in (\[eq:4.5\]) as functions of $x$ and $y$ only, and subsequently we replace the integration constant by arbitrary function of $t$. Since equation (\[eq:4.5\].a) is not quasilinear, the CSM as presented in Section \[sec:condsym\] cannot be applied to the entire system (\[eq:4.5\]). Nevertheless, the CSM can be used for the subsystem $\Delta$ of 3 equations in the 4 unknowns consisting of (\[eq:4.5\].b)-(\[eq:4.5\].d), which is then underdetermined. Next, we use the degrees of freedom of the obtained solution to satisfy equation (\[eq:4.5\].a). The subsystem $\Delta$ can be written in the following matrix form:

\[eq:4.6\]

(

[cccc]{} 0 & -1 & 0 & 0\
0 & 0 & 1 & 0\
0 & 0 & 0 & -1

)

(

[c]{}\
\
u\
v

)

+

(

[cccc]{} 1 & 0 & 0 & 0\
0 & 0 & 0 & 1\
0 & 0 & 1 & 0

)

(

[c]{}\
\
u\
v

)

=

(

[c]{} 0\
0\
0\
0

)

.

The dispersion relation (\[eq:2.3\]) takes the form

\[eq:4.7\]

(

[cccc]{} \_2 & -\_1 & 0 & 0\
0 & 0 & \_1 & \_2\
0 & 0 & \_2 & -\_1

)

&lt;3.

By solving the dispersion relation (\[eq:4.7\]), we obtain the wave vector $\lambda=(\lambda_1,\lambda_2)=(1,i)$ and its complex conjugate $\bar{\lambda}=(1,-i)$. Hence, the Riemann invariants associated with the wave vectors $\lambda$ and $\bar{\lambda}$ are

\[eq:4.9\] r=x+iy,|[r]{}=x-iy.

We look for a nontrivial real solution of (\[eq:4.5\].b)-(\[eq:4.5\].d) of the form

\[eq:4.10\] =&f\_1(r,|[r]{})+\_1(r,|[r]{}),&&=f\_2(r,|[r]{})+\_2(r,|[r]{}),\
u=&f\_3(r,|[r]{})+\_3(r,|[r]{}), &&v=f\_4(r,|[r]{})+\_4(r,|[r]{}).

Introducing the notation $\phi=u^1$, $\psi=u^2$, $u=u^3$, $v=u^4$ and $$\Lambda={\left( \begin{array}{cc}
                1 & i \\
                1 & -i
              \end{array}
 \right)},$$we observe that the matrices $\eta_{a_j}$ defined by (\[eq:3.28b\]) all vanish. Consequently, all trace conditions (\[eq:3.34\]) used to obtain solutions of the form (\[eq:4.10\]) are identically satisfied. We still have to consider the trace conditions (\[eq:3.30\]) where the matrices $\mathcal{A}^\mu$, $\mu=1,2,3$, are given by $$\mathcal{A}^1={\left( \begin{array}{cccc}
           0 & -1 & 0 & 0 \\
           1 & 0 &  0 & 0
         \end{array}
 \right)},\qquad \mathcal{A}^2={\left( \begin{array}{cccc}
           0 & 0 & 1 & 0 \\
           0 & 0 & 0 & 1
         \end{array} \right)},\qquad \mathcal{A}^3={\left( \begin{array}{cccc}
           0 & 0 & 0  & -1 \\
           0 & 0 & -1 & 0
         \end{array} \right)}.$$Conditions (\[eq:3.30\]) take the form

\[eq:4.11\] &(a)+--+i[( +++ )]{},\
&(b)++++i[( +-- )]{},\
&(c)+--+i[( +++ )]{}.

Condition (\[eq:4.11\].a) implies that

\[eq:4.12\] e[( )]{}=-m[( )]{},

while conditions (\[eq:4.11\].b) and (\[eq:4.11\].c) imply that the solutions for velocities $u$ and $v$ take the particular form

\[eq:4.13\] u=h(r)+|[h]{}(|[r]{}),v=i[( h(r)-|[h]{}(|[r]{}) )]{},

where $r$ and $\bar{r}$ are given by (\[eq:4.9\]). Note that the most general solution for $u$ and $v$ of the trace conditions (\[eq:3.30\]-\[eq:3.32\]) is exactely the one given by the expression (\[eq:4.13\]). The function $h$ and its complex conjugate are also determined using condition (\[eq:4.12\]) and the PDE (\[eq:4.5\].a). This allows us to also determine the functions $\phi$ and $\psi$. For computational purposes, it is convenient to introduce velocities of the form (\[eq:4.13\]) in equation (\[eq:4.1\].c) in order to determine the angle $\theta$ in the form

\[eq:4.14\] =[( i )]{},

where we denote $h^{(n)}=d^nh/dr^n$. Next, we substitute (\[eq:4.14\]) into equation (\[eq:4.2\]) in order to determine the function $h$ and its complex conjugate. Proceeding in this way, we find that $h$ and $\bar{h}$ satisfy the third-order ODE which is separable in $r$ and $\bar{r}$. It is therefore equivalent to the system

\[eq:4.15\] (a)2-3=,(b)2 -3=,

where $\Omega$ is a real separation constant. Equation (\[eq:4.15\].b) can be solved in a way similar way to (\[eq:4.15\].a). Defining

\[eq:4.16\] g(r)=h\^[(1)]{}(r),

equation (\[eq:4.15\].a) can be written as

\[eq:4.17\] gg”-3/2 (g’)\^2+(/2) g\^3=0.

According to [@Kamke:1979] equation (\[eq:4.17\]) admits a first integral which can be rewritten in the equivalent form

\[eq:4.19\] =d.

By the change of variable

\[eq:4.20\] g=[( )]{},

equation (\[eq:4.19\]) can be transformed to $2\exp{\left( (2\Omega)^{-1}(\omega^2-c_1) \right)}d\omega=\pm
\Omega d\xi,$ for which the solution, in term of the inverse error function $\operatorname{erf}^{-1}$, is

\[eq:4.21\] =,

where $c_2$ is a constant of integration. We substitute (\[eq:4.20\]) into (\[eq:4.16\]) with $\omega$ given by (\[eq:4.21\]) and integrate the result. We obtain the function $h(r)$ and its complex conjugate in terms of the error functions $\operatorname{erfi}$ and $\operatorname{erf}^{-1}$ as

\[eq:4.22\] h(r)&=-+c\_3,\
|[h]{}(|[r]{})&=-+|[c]{}\_3,

where the $c_i\in{\mathbb{C}}$ are integration constants and the real separation constant $\Omega$ is non zero. Since equations (\[eq:4.22\]) solve the system (\[eq:4.15\]), we have that $\theta$ given by (\[eq:4.14\]) satisfies the compatibility condition (\[eq:4.2\]) for $\sigma$. One should note that no derivatives with respect to time appear in the PDE (\[eq:4.2\]). Consequently, the equation (\[eq:4.2\]) is still satisfied if, in equations (\[eq:4.13\]) and (\[eq:4.14\]), we replace the function $h(r)$ and its complex conjugate respectively by

\[eq:4.22bis\] h(t,r)&=-+c\_3(t),\
|[h]{}(t,|[r]{})&=-+|[c]{}\_3(t),

which are obtained by substituting the arbitrary complex functions $c_i(t)$ in place of the integration constant $c_i$ and the arbitrary real function $\Omega(t)$ in place of the separation constant $\Omega$. Hence, the general solution of the system (\[eq:4.1\]) takes the form

\[eq:4.23\] u(t,x,y)&=h(t,r)+|[h]{}(t,|[r]{}),\
v(t,x,y)&=i[( h(t,r)-|[h]{}(t,|[r]{}) )]{},\
(t,x,y)&=[( i )]{},\
(t,x,y)&=-V(t,x,y)+[12]{}(2(t,x,y))++u\_t(t,x,y)dx\
& +\_y(t,x,y)(2(t,x,y))dx+\_0(t),

where the functions $h$ and $\bar{h}$ are defined by (\[eq:4.22bis\]). Note that the solution (\[eq:4.23\]) is a real-valued solution.

#### 

Let us now consider the situation when $\Omega$ is identically equal to zero in the system of ODEs (\[eq:4.15\]). This leads to two possible independent solutions for $h(r)$ depending on whether its second order derivative $h^{(2)}(r)$ vanishes or not. These solutions are respectively

\[eq:4.24\] (i) &h(r)=c\_1r+c\_2,&&h\^[(2)]{}(r)=0,\
(ii)&h(r)=+c\_3,&&h\^[(2)]{}(r)0,

and their complex conjugates, where the $c_i$ are integration constant. Let us consider separately two cases: when the function $h$ is given by (\[eq:4.24\].i) and when $h$ is of the form (\[eq:4.24\].ii).

#### Case i.

In this case, the real solution of the original system (\[eq:4.1\]) takes the form

\[eq:4.25\] u(t,x,y)=&4[( e(c\_1(t))x+m(c\_1(t))y+e(c\_2(t)) )]{},\
v(t,x,y)=&4[( m(c\_1(t))x-e(c\_1(t))y+m(c\_2(t)) )]{},\
(t,x,y)=& {

&-[( )]{},&&m(c\_1(t))0,\
& && m(c\_1(t))=0,

.\
(t,x,y)=&-V(t,x,y)+[( 2|c\_1(t)|\^2+e(\_1(t)) )]{}x\^2 -2m(\_1(t))xy\
&+[( 2|c\_1(t)|\^2-e(\_1(t)) )]{}y\^2+2[( 2e[( c\_1(t)|[c]{}\_2(t) )]{}+e(\_2(t)) )]{}x\
&-2[( 2m[( c\_1(t)|[c]{}\_2(t) )]{}+m(\_2(t)) )]{}y+\_0(t),

where the real function $\sigma_0(t)$ and the complex functions $c_1(t)$, $c_2(t)$ are arbitrary functions of time and $\dot{c}_1(t)$, $\dot{c}_2(t)$ their respective derivatives with respect to $t$. It should be noted that if the potential $V$, the real function $\sigma_0(t)$ and the complex functions $c_1(t)$, $c_2(t)$, are all bounded and are their derivatives $\dot{c}_1(t)$, $\dot{c}_2(t)$ are also bounded, then the solution (\[eq:4.19\]) is bounded. Moreover, if we take these functions to be of the form

\[eq:4.25\] \_0(t)=a\_0e\^[-s\_0t]{},c\_i(t)=a\_je\^[-s\_jt]{}+ib\_je\^[-q\_jt]{},j=1,2,

where $0<q_j,s_j\in {\mathbb{R}}$, $a_j,b_j\in {\mathbb{R}}$, $j=0,1,2$, then we obtain a damped solution. Another physically interesting situation occurs when the arbitrary functions $c_i(t)$ are constant in solution (\[eq:4.19\]). In this case, we obtain a stationary solution for velocities $u$ and $v$. This allows us to draw the shape of those extrusion dies which are admissible by requiring that the tool walls coincide with the lines of flow generated by the velocities $u$ and $v$. This is necessary since we suppose that the flow is incompressible. From the practical point of view, it is convenient to press the plastic material through the die rectilinearly with constant speed. If we feed the die in this way, the plasticity region limits are curves defined by the ODE [@Lamothe:1; @Lamothe:2]

\[eq:4.26\] =,

where $U_0$ and $V_0$ are constants representing the feeding velocity (or extraction velocity) of the die along the $x$-axis and $y$-axis respectively. Equation (\[eq:4.26\]) is a consequence of the hypothesis that the flow is incompressible and the mass is conserved. This condition reduces to the boundary conditions described in [@Czyz:1974] when we require that the limits of the plasticity region correspond to the slip lines (which correspond to the characteristic curves of the original system (\[eq:4.1\])), that is, when we require that $dy/dx=\tan\theta$ or $dy/dx=-\cot\theta$. Here, we use the weakened condition (\[eq:4.26\]) because no constraints are imposed on the flow lines which are in contact with the tool walls. This is so because, for a given solution and given parameters, we choose the walls of the extrusion die to lie along the flow lines. For the purpose of illustrating the applicability of the method, we have drawn in figure \[fig:1\] the shape of a tool and the flow of matter inside the extrusion die for the following parameters : $\mathcal{R}e(c_1(t))=1$, $\mathcal{I}m(c_1(t))=0$, $\mathcal{R}e(c_2(t))=0$, $\mathcal{I}m(c_2(t))=0$. The feeding velocity of the tool is $U_0=5.95$, $V_0=0$ and the tool expels the matter at a velocity of $U_1=24.05$, $V_1=0$. The plasticity region at the opening is bounded by curve $C_1$ and at the exit by the curve $C_2$. This extrusion die can thin a plate or rod of ideal plastic material.

![Extrusion die corresponding to the solution (\[eq:4.19\]).[]{data-label="fig:1"}](outil_cas_i.jpg){width="3.5in"}

#### Case ii.

If $h$ is defined by (\[eq:4.17\].ii), then the corresponding nontrival solution of system (\[eq:4.1\]) takes the form

\[eq:4.27\] u(t,x,y)=&e(c\_3(t))+2,\
v(t,x,y)=&e(c\_3(t))+2,\
(t,x,y)=&-(\^[-1]{}\
&),

where the mean pressure $\sigma$ is given by (\[eq:4.3\]), where we substitute the values of functions $u,v,\theta$ given by (\[eq:4.21\]). The complex functions $c_i(t)$, $i=1,2,3$, and the real function $\sigma_0(t)$ which appear in this solution are arbitrary. For any time $t_0$, solution (\[eq:4.3\]), (\[eq:4.27\]) has a singularity at point $(x_0,y_0)$ which satisfies the equation $$x_0^2+y_0^2+2{\left( \mathcal{R}e(c_2(t_0))x_0+\mathcal{I}m(c_2(t_0))y_0 \right)}+|c_2(t_0)|^2=0.$$This singularity is stationary if the function $c_2(t)$ is constant. Otherwise, its position varies with time. If the functions $\sigma_0$ and $c_i$, $i=1,2,3$, are of the form (\[eq:4.25\]) defined in a region of the $xy$-plane on a time interval $[T_0,T)$ where the gradient catastrophe does not occur, then the solution is bounded and damped. In figure \[fig:2\], we have drawn the shape of an extrusion die corresponding to the solution (\[eq:4.27\]) for the following choice of parameters: $\mathcal{R}e{\left( c_1(t) \right)}=0$, $\mathcal{I}m{\left( c_1(t) \right)}=0$, $\mathcal{R}e{\left( c_2(t) \right)}=0$, $\mathcal{I}m{\left( c_2(t) \right)}=-0.5$, $\mathcal{R}e{\left( c_3(t) \right)}=-0.5$ and $\mathcal{I}m{\left( c_3(t) \right)}=0$. The feeding velocity has component $U_0=0.2$, $V_0=0.2$, and the extraction of material is performed at the velocity $U_1=0.2$, $V_1=-0.2$. This type of tool can be used to bend a rod by extrusion without having to fold it. Finally, we should emphasize that the flow changes considerably when the parameters are varied and we have the freedom to choose the walls of the tool among the flow lines for certain fixed parameters $c_i(t)$. Moreover, the velocity and the orientation of the feeding (extraction) can vary somewhat for a given shape of the tool. Consequently, many types of extrusion dies can be drawn.

![Extrusion die corresponding to the solution (\[eq:4.21\]).[]{data-label="fig:2"}](outil_cas_ii.jpg){width="3.5in"}

Simple wave solutions of inhomogeneous quasilinear system. {#sec:ms}
==========================================================

\[sec:5\] Consider an inhogeneous first-order system of $q$ quasilinear PDEs in $p$ independent variables and $q$ unknowns of the form

\[eq:SW:1\] \^[i]{}\_(u) u\^\_i=b\^(u),,=1,…,q,i=1,…,p.

Let us underline that the system can be either hyperbolic or elliptic. We are looking for real solutions describing the propagation of a simple wave which can be realized by the system (\[eq:SW:1\]). We postulate a form of the solution $u$ in terms of a Riemann invariant $r$, [*i.e.* ]{}

\[eq:SW:2\] u=f(r),r=\_i(u) x\^i, i=1,…,p,

where $\lambda(u)={\left( \lambda_1(u),\ldots,\lambda_p(u) \right)}$ is a real-valued wave vector. We evaluate the Jacobian matrix $u^\beta_i$ by applying the chain rule: $$u^\beta_i=\frac{{\partial}f^\beta}{{\partial}r}{\left( r_{x^i}+r_{u^\alpha}u^\alpha_i \right)}=\frac{{\partial}f^\beta}{{\partial}r}{\left( \lambda_i+\lambda_{i,u^\alpha}x^iu^\alpha_i \right)}.$$We assume that the matrix

\[eq:SW:4\] =[( I\_q- )]{}\^[qq]{},

is invertible, where we have used the following notation ${\partial}f/{\partial}r={\left( {\partial}f^1/{\partial}r,\ldots, {\partial}f^q/{\partial}r \right)}^T$, ${\partial}r/{\partial}u ={\left( {\partial}r/{\partial}u^1,\ldots,{\partial}r/{\partial}u^q \right)}$. The Jacobian matrix ${\partial}u$ takes the form

\[eq:SW:5\] u=\^[-1]{}\^[qq]{}.

Replacing the Jacobian matrix (\[eq:SW:5\]) into the original system (\[eq:SW:1\]), we get

\[eq:SW:7\] \^[i]{}\_[( \^[-1]{} )]{}\^\_\_i=b\^.

Note that the expression $\Phi^{-1}\frac{{\partial}f}{{\partial}r}\in {\mathbb{R}}^q$ is a contravariant vector as well as $b\in{\mathbb{R}}^q$. Hence there exists a nonzero scalar function $\Omega=\Omega(x,u)$, a rotation matrix $L=L(x,u)\in SO(q)$ and a vector $\tau=\tau(x,u)\in{\mathbb{R}}^q$ such that

\[eq:SW:8\] \^[-1]{}=L b+,\^i\_i=0.

It should be noted that we assume appropriate levels of differentiability of the functions $\Omega$, $\lambda$, $L$, $\tau$ and $b$, as necessary in order to justify all the following steps. Using relation (\[eq:SW:8\]), we eliminate the vector $\Phi^{-1}\frac{{\partial}f}{{\partial}r}\in {\mathbb{R}}^q$ from equation (\[eq:SW:7\]), which allows us to factor out the vector $b$ on the right after regrouping all terms on the left of the resulting equation. Therefore, we obtain the condition

\[eq:SW:13\] [( \^[i]{}\_i L-I\_q )]{}b=0,

on the scalar function $\Omega$, the wave vector $\lambda$ and the rotation matrix $L$. This implies that we have the following dispersion relation

\[eq:SW:14\] [( \^[i]{}\_i L-I\_q )]{}=0.

Once a scalar function $\Omega$, a wave vector $\lambda$ and a matrix $L$ satisfying (\[eq:SW:13\]) have been obtained, equation (\[eq:SW:8\]) must be used in order to determine the function $f$. Replacing the expression (\[eq:SW:4\]) for the matrix $\Phi$ into equation (\[eq:SW:8\]) and solving for the vector $\frac{{\partial}f}{{\partial}r}$ and taking into account the relation ${\partial}r/{\partial}u=({\partial}\lambda_i/{\partial}u) x^i$, we find that

\[eq:SW:15\] =,

which cannot admit the gradient catastrophe. Indeed, if we suppose that $1+\Omega({\partial}\lambda_i/{\partial}u) x^i L b=0$ when we proposed from equation (\[eq:SW:8\]) to equation (\[eq:SW:15\]), we easily conclude that $\Omega L b=0$. Consequently, since the matrix $\Phi$ is invertible, we conclude from (\[eq:SW:8\]) that the solution for $f(r)$ is constant, so it cannot admit the gradient catastrophe if $\det \Phi\neq 0$.

#### 

Up till now, we cannot be sure that system (\[eq:SW:15\]) for $f(r)$ is well-defined in the sense that it represents a system for $f(r)$ express in terms of $r$ only. To ensure this, we begin by introducing vector fields orthogonal to the wave vector $\lambda$, that is, vector fields of the form

\[eq:SW:16\] X\_a=\^i\_a(u)\_[x\^i]{},a=1,…,p-1,i=1,…,p,

where

\[eq:SW:17\] \_a\^i\_i=0,a=1,…,p-1,(\_a\^i)=p-1.

Consequently, the wave vector $\lambda$ together with the vectors $\xi_a$, $a=1,\ldots,p-1$ form a basis for ${\mathbb{R}}^p$. Moreover, let us note that the vector fields $X_a$ form an Abelian Lie algebra of dimension $p-1$. Application of the vector field (\[eq:SW:16\]) to equation (\[eq:SW:15\]) cancels out the left side since we suppose that $f=f(r)$, so it must be the same on the right side. Therefore, the conditions for the system of ODEs (\[eq:SW:15\]) to be well-defined in terms of $r$ are

\[eq:SW:18\] X\_a=0,a=1,…,p-1.

#### 

In summary, system (\[eq:SW:1\]) admits a simple wave solution if the following conditions are satisfied:

-   there exist a scalar function $\Omega(x,u)$, a wave vector $\lambda(u)$ and a rotation matrix $L(x,u)$ satisfying the wave equation (\[eq:SW:13\]);

-   there exist $p-1$ vector fields (\[eq:SW:16\]) which satisfy the orthogonality relation (\[eq:SW:17\]);

-   the right-hand side of equation (\[eq:SW:15\]) is annihilated by the vector fields (\[eq:SW:16\]), [*i.e.* ]{}conditions (\[eq:SW:18\]) are satisfied;

-   $\det \Phi\neq 0$, where $\Phi$ is given by (\[eq:SW:4\]);

-   $1+\Omega \frac{{\partial}\lambda_i}{{\partial}u}Lbx^i\neq 0$, $\lambda\in C^1$.

These conditions are sufficient, but not necessary, since the vanishing of $\det\Phi$ does not imply that solutions of form (\[eq:SW:2\]) do not exist. We finish the present analysis of the simple wave solutions of system (\[eq:SW:1\]) by several remarks:

-   In general, there are more parameters (arbitrary functions) defining $\Omega$ and the rotation matrix $L$ than the minimum number required to satisfy condition (\[eq:SW:13\]). The remaining arbitrary quantities are arbitrary functions of the $x$’s and $u$’s, which have to be used to satisfy the conditions (\[eq:SW:18\]). However, in the particular case when the system (\[eq:SW:1\]) has two equations in two dependent variables, the two-dimensional matrix $L$ is defined by a single parameter and $\Omega$ is the only other parameter available to satisfy the condition (\[eq:SW:13\]). Thus, since the system is autonomous (it can be expressed solely in terms of the $u$’s and their derivatives), it is clear that, in that case, $\Omega$ and $L$ depend only on the dependent variables $u$’s.

-   Supposing that $b$ is continous and assuming that conditions (i-v) are satisfied, the right side of equation (\[eq:SW:15\]) is continuous, which ensures the existence and uniqueness of the solution for $f(r)$ (see for example [@Ince]).

-   A sufficient condition for the system (\[eq:SW:15\]) to be expressible in terms of $r$ only ([*i.e.* ]{}well-defined) is: $$\frac{{\partial}\lambda_i}{{\partial}u}L b=0,\qquad i=1,\ldots, p.$$This condition is trivially satisfied when the vector $\lambda$ is constant. More generally, it is sufficient to require that $\frac{{\partial}\lambda_i}{{\partial}u}Lbx^i$ be proportional to $r=\lambda_i x^i$, [*i.e.* ]{}that there exist $\beta(x,u)$ such that $$\frac{{\partial}\lambda_i}{{\partial}u}L b=\beta(x,u)\lambda_i,\qquad i=1,\ldots,p.$$

-   If the matrix $A^i\lambda_i$ is invertible, then $$\frac{{\partial}\lambda_i}{{\partial}u}Lbx^i=\frac{{\partial}\lambda_i}{{\partial}u}(\mathrm{A}^j\lambda_j)^{-1}(\mathrm{A}^j\lambda_j) Lb x^i=\frac{1}{\Omega}\frac{{\partial}\lambda_i}{{\partial}u}(\mathcal{A}^j\lambda_j)^{-1}bx^i$$is satisfied due to condition (\[eq:SW:13\]). From the previous remark, we deduce $$\frac{1}{\Omega}\frac{{\partial}\lambda_i}{{\partial}u}(\mathcal{A}^ j\lambda_j)^{-1}b=\beta(u) \lambda_i,\quad i=1,\ldots,p,$$which is a sufficient condition where the rotation matrix $L$ is not involved when the matrix $A^j\lambda_j$ is invertible.

We note that the approach presented above generalizes the results obtained in [@GrundlandZelazny:1983] where the simple states have been constructed with wave vectors $\lambda$ of constant direction for hyperbolic inhomogeneous systems (\[eq:SW:1\]).

Simple mode solutions for inhomogeneous quasilinear system. {#sec:6}
===========================================================

We now generalize the concept of simple wave solutions for inhomogeneous, quasilinear, systems of form (\[eq:SW:1\]). As in the case of the simple wave, the system can be either hyperbolic or elliptic. We look for a real solution, in terms of a Riemann invariant $r$ and its complex conjugate $\bar{r}$, of the form

\[eq:ms:1\] u=f(r,|[r]{}),r(u,x)=\_i(u)x\^i,|[r]{}(u,x)=|\_i(u)x\^i,i=1,…, p,

where $\lambda(u)={\left( \lambda_1(u),\ldots,\lambda_p(u) \right)}$ is a complex wave vector and $\bar{\lambda}(u)$ its complex conjugate. The Jacobian matrix takes the form

\[eq:ms:2\] u=\^[-1]{}[( + | )]{}\^[qp]{},

where we assume that the matrix

\[eq:ms:3\] =[( I\_q-- )]{}\^[qq]{}.

is invertible. Replacing the Jacobian matrix (\[eq:ms:2\]) into the system (\[eq:SW:1\]), we obtain

\[eq:ms:4\] \^[i]{}\^[-1]{}[( + | )]{}=b,\^i=[( A\^[i]{}\_)]{}\^[qq]{}.

We introduce a rotation matrix $L=L(x,u)\in SO(q,{\mathbb{C}})$ and its complex conjuguate $\bar{L}=\bar{L}(x,u)\in SO(q,{\mathbb{C}})$ such that the relations

\[eq:ms:5\] \^[-1]{}=L b + ,\^[-1]{}=| |[L]{} b+ |,

hold, where $\Omega(x,u)$ and its conjugate $\bar{\Omega}$ are scalar complex functions, while $\tau(x,u)$ and its complex conjugate vector $\bar{\tau}(x,u)$ satisfy

\[eq:ms:7\] \^i\_i +\^i|\_i |=0.

The vectors $\tau$ and $\bar{\tau}$ can be seen as characteristic vectors of the homogeneous part of equation (\[eq:ms:4\]). We eliminate the vectors $\Phi^{-1}({\partial}f/{\partial}r)$ and $\Phi^{-1}({\partial}f/{\partial}\bar{r})$ from equation (\[eq:ms:4\]) using equations (\[eq:ms:5\]). Considering the equation (\[eq:ms:7\]), we obtain, as a condition on functions $\Omega$, $\lambda$, $L$ and their complex conjugates, the relation

\[eq:ms:8\] [( \^i[( \_iL+|\_i||[L]{} )]{}-I\_q )]{}b=0.

The vector $b$ can be identified with an eigenvector of equation (\[eq:ms:8\]). Since $b$ is known from the initial system (\[eq:SW:1\]), the equation has to be satisfied by an appropriate choice of functions $\Omega$, $\bar{\Omega}$, of wave vectors $\lambda$, $\bar{\lambda}$ and of rotation matrices $L$ and $\bar{L}$. In particular, this requires that the dispersion relation

\[eq:ms:10\] [( \^i[( \_i L+|\_i|[L]{} )]{}-I\_q )]{}=0

holds. Equation (\[eq:ms:10\]) is an additional condition on $\lambda$ and $\bar{\lambda}$, $L$ and $\bar{L}$ for the equation (\[eq:ms:8\]) to have a solution. Multiplying equations (\[eq:ms:5\]) on the left by the matrix $\Phi$, writing the matrix $\Phi$ explicitly using the notations (\[eq:ms:3\]), and then solving for ${\partial}f/{\partial}r$ and ${\partial}f/{\partial}\bar{r}$, we obtain the system

\[eq:ms:12\] =,=,

where the scalar functions $\sigma_1$ and $\sigma_2$ and their complex conjugates are defined by the equations

\[eq:ms:13\] \_1=(L b+),\_2=(L b+).

In order to ensure that system (\[eq:ms:12\]) is well-defined in the sense that it can be expressed as a system for $f$ in terms of $r$ and $\bar{r}$, we introduce the vector fields

\[eq:ms:14\] X\_a=\^i\_a(u)\_[x\^i]{},a=1,…, p-2,

where the complex coefficients $\xi^i_a(u)$ satisfy the orthogonality relations

\[eq:ms:15\] \^i\_a\_i=0,\^i\_a|\_i=0.

Next, we apply them to equations (\[eq:ms:12\]), which gives us the following conditions

&X\_a=0,\
&X\_a=0,

$a=1,\ldots, p-2$, since the vector fields $X_a$ annihilate all the functions $f(r,\bar{r})$.

#### Remark 1:

The system (\[eq:ms:12\]) for $f(r,\bar{r})$ is not necessarily integrable. However, it is possible to use the arbitrary functions defining the function $\Omega(x,u)$, the wave vectors $\lambda(u)$, $\bar{\lambda}(u)$, the vectors $\tau(x,u)$, $\bar{\tau}(x,u)$ and the matrices of rotation $L(x,u)$, $\bar{L}(x,u)$, in order to satisfy the compatibility conditions of the system (\[eq:ms:12\]).

Example. Nonlinear interaction of waves and particles. {#sec:7}
======================================================

Let us consider the inhomogeneous system describing the propagation of shock waves intensity in nonlinear interaction of waves and particles [@Luneburg:1964]

\[eq:e1:1\] u\_x+\_y=2\^[1/2]{}a (u/2)(/2),u\_y-\_x=-2\^[1/2]{}a(u/2)(/2).

It should be noted that the compatibility condition of the mixed derivatives of $\phi$ corresponds to the Liouville equation

\[eq:e1:2\] u\_[xx]{}+u\_[yy]{}=a\^2.

Therefore, each solution of the system (\[eq:e1:1\]) also gives us a solution of the Liouville equation (\[eq:e1:2\]). We use the methods presented in Section \[sec:6\] to obtain the general solution of the system (\[eq:e1:1\]). First, write the system (\[eq:e1:1\]) in the matrix form

\[eq:e1:3\]

(

[c c]{} 1 & 0\
0 & -1

)

(

[c ]{} u\_x\
\_x

)

+

(

[c c]{} 0 & 1\
1 & 0

)

(

[c]{} u\_y\
\_y

)

=

(

[c]{} b\_1\
b\_2

)

,

where $$b_1=2^{1/2}a \exp(u/2)\sin(\phi/2),\qquad b_2=2^{1/2}a\exp{u/2}\cos(\phi/2).$$The matrices $\mathcal{A}^i$ introduced in equation (\[eq:ms:4\]) are given by $$\mathcal{A}^1={\left( \begin{array}{c c}
1 & 0\\
0 & -1
\end{array} \right)},\qquad \mathcal{A}^2={\left( \begin{array}{c c}
0 & 1\\
1 & 0
\end{array} \right)}.$$Condition (\[eq:ms:8\]) is then satisfied by the scalar function $\Omega$, the wave vector $\lambda$ and the rotation matrix $L$, defined by $$\Omega=12^{1/4}(1-\epsilon^i),\quad \lambda=(1,i),\quad L={\left( \begin{array}{c c} l_{11} & l_{12}\\
l_{21} & l_{22}\end{array} \right)},\quad \epsilon=\pm 1,$$together with their complex conjugates, where $$\begin{aligned}
l_{11}&=l_{22}=-108^{1/4}i{\left( \frac{3^{1/2}\epsilon (b_1+ib_2)^2+i(b_1-ib_2)^2}{6(1-\epsilon i) (b_1^2+b_2^2)} \right)}\\
l_{12}&=-l_{21}=108^{1/4}i{\left( \frac{3^{1/2}\epsilon i (b_1+ib_2)^2+(b_1-ib_2)^2}{6(1-\epsilon i) (b_1^2+b_2^2)} \right)}
\end{aligned}$$Since the wave vectors $\lambda$ and $\bar{\lambda}$ are constant, the quantities $\sigma_1$, $\sigma_2$ (and their conjugates) are zero. So, we obtain the system (\[eq:ms:12\]) written in terms of the Riemann invariants $r=x+ i y$ and $\bar{r}=x- i y$

\[eq:e1:5\] u\_r&=24\^[-1/2]{}a[( -i(3\^[1/2]{}+3 i )(/2)+(3\^[1/2]{}-3 i )(/2) )]{}(u/2)+\
\_r&=24\^[-1/2]{}a[( (3\^[1/2]{}-3 i )(/2)-i(3\^[1/2]{}+3 i )(/2) )]{}(u/2)+i\
u\_[|[r]{}]{}&=24\^[-1/2]{}a[( i(3\^[1/2]{}-3 i )(/2)+(3\^[1/2]{}+3 i )(/2) )]{}(u/2)+|\
\_[|[r]{}]{}&=24\^[-1/2]{}a[( i(3\^[1/2]{}-3 i )(/2)+(3\^[1/2]{}+3 i )(/2) )]{}(u/2)-i|\

where $\sigma(r,\bar{r})$ is an arbitrary scalar function defining the vector $\tau=(\sigma,i \sigma)^T$, which satisfies condition (\[eq:ms:7\]). The system (\[eq:e1:5\]) can be written in the more compact form

\[eq:e1:6\] &a)&&[( u+i )]{}=2\^[-1/2]{}a i [( )]{},\
&b)&&[( u-i )]{}=-2\^[-1/2]{}a i [( )]{},\
&c)&&[( u-i )]{}=-(3/2)\^[1/2]{}a [( )]{}+2,\
&d)&&[( u+i )]{}=-(3/2)\^[1/2]{}a [( )]{}+2|.\

Since, the equations (\[eq:e1:6\].c) and (\[eq:e1:6\].d) are complex conjugates they can be satisfied by defining the quantity $\sigma$ annihilating (\[eq:e1:6\].c). Consequently, there are only the two equations (\[eq:e1:6\].a) and (\[eq:e1:6\].b) to solve. After the change of variable

\[eq:e1:7\] f=u+i ,|[f]{}=u-i ,

equations (\[eq:e1:6\].a) and (\[eq:e1:6\].b) take the form

\[eq:e1:8\] f\_r=2\^[-1/2]{} i a|f|\^2,f\_[|[r]{}]{}=-2\^[-1/2]{} i a|f|\^2,

which have general solution

\[eq:e1:8\] f(r,|[r]{})=,

where $\psi$ and $\bar{\psi}$ are arbitrary functions of $r$ and $\bar{r}$, respectively. Replacing the solution (\[eq:e1:8\]) into equations (\[eq:e1:7\]), solving for $u$ and $\phi$, and using the Riemann invariants $r$ and $\bar{r}$, we obtain the general solution of the system (\[eq:ms:1\]). This solution depends on one arbitrary complex function $\psi$ of one complex variable $r$ and its complex conjugate function.

u&=2[( )]{},\
&=-i[( - )]{}+2 n,

where $n$ is odd and the term $2 n\pi$ is associated with the admissible branches of the logarithm. It is easily verified that the solution for $u$ satisfies the Liouville equation (\[eq:ms:2\]).

Final remarks. {#sec:8}
==============

The generalized method of characteristics was originally devised for solving first-order quasilinear hyperbolic systems. The proposed techniques described in Sections \[sec:3\], \[sec:5\] and \[sec:6\] allow us to extend the applicability of this approach not only to hyperbolic systems but also to encompass elliptic and mixed (parabolic) type systems, both homogeneous and inhomogeneous. A variant of the conditional symmetry method for obtaining multimode solutions has been proposed for these types of systems. We have demonstrated the usefulness of this approach through the examples of nonlinear interaction of waves and particles and of the ideal plasticity in $(2+1)$ dimensions in its elliptic region. New classes of real solutions have been constructed in closed form, some of them bounded. Some of the obtained solutions described a stationary flow for an appropriate choice of parameters. For these solutions, we have drawn extrusion dies and the vector fields which define the flow inside a region where the gradient catastrophe does not occur.

#### 

The proposed approach for constructing multimode solutions can be used in several potential applications arising from systems describing nonlinear phenomena in physics. It should be noted that in the multidimensional case, for many physical models, there are few known examples of multimode solutions written in terms of Riemann invariants for elliptic systems. This is a motivating factor for the elaboration of the generalization of the methods presented in Sections \[sec:5\] and \[sec:6\] through the introduction of rotation matrices in the factorization of the Jacobian matrices (\[eq:ms:2\]). This fact weakens the integrability condition required in the expression (\[eq:ms:12\]). The approach proposed in this paper offers a new and promising way to construct and investigate such types of solutions. This makes our approach attractive since it can widen the potential range of applications leading to more diverse types of solutions.

Acknowledgements {#acknowledgements .unnumbered}
================

This work was supported by a research grant from the Natural Sciences and Engineering Council of Canada (NSERC).

[10]{} , , London Math.Soc. Lecture note ser 149, Cambrige Univ. Press. (1992).

G. Boillat, , Gauthier-Villars, Paris, (1965).

, The method of Riemann invariants and its applications to the theory of plasticity, ,part 1, 26, 6, 1974, 817-838 and part 2, 24, 1, (1972), 3-26.

, , Chapt. 1: Les systèmes différentielles en involution, Gauthier-Villars, Paris, (1953).

J. Chakrabarty, , Elsevier, (2006).

R. Conte, A.M. Grundland, B. Huard, Elliptic solutions of isentropic ideal compressible fluid flow in $(3+1)$ dimensions, ,Math. Theor. 42, 135203, (2009) (14pp)

R. Courant, D. Hilbert, , Interscience, New-York, (1962).

J. Czyz, Construction of a flow of an ideal plastic material in a die, on the basis of the method of [R]{}iemann invariants, , 26(4), (1974), 589-616.

, Simple waves and invariant solutions of quasiliear systems, ,37, 6, (1996), 2969-2979.

, Hamiltonnian formalism of one-dimensional systems of hydrodynamic type, ,27, (1983), 665-669.

, On the integrability of $(2+1)$ dimensional quasilinear systems, ,248, (2004), 187-206.

, The characterization of two-component (2+1)-dimensional integrable systems of hydrodynamic type, ,37, (2004), 2949-2963.

, Conditional symmetry of equations of mathematical Physics, ,43, (1991), 1456-1470.

A.M. Grundland, B. Huard, Riemann invariants and rank-$k$ solutions of hyperbolic systems, , 13,3 , (2006) 393-419.

A.M. Grundland, B. Huard, Conditional symmetries and Riemann invariants for hyperbolic systems of PDEs, , Math. Theor. 40, (2007), 4093-4123.

, Symmetry reduction and Riemann wave solutions, ,198, (1996), 879-892.

A.M. Grundland, P. Vassiliou, On the solvability of the Cauchy problem for the Riemann double waves by the Monge-Darboux method. , (1991) 221-278.

A.M. Grundland, R. Zelazny, Simple waves in quasilinear systems, Part I and Part II, , 24, 9, (1983), 2305-2329.

, Riemann invariants for nonhomgeneous systems of quasilinear partial differential equations, , Ser Sci Techn 22, (1974), 273-282.

R. Hill, , Oxford University press, 1998.

E.L. Ince, , Dover Publ. New-York, 1956.

A. Jeffrey, , Pitman Publ., 1976.

, Almost global existence of nonlinear wave equations in three space dimensions, ,37, (1984), 443-455.

E. Kamke, , , I B.G. Teubner Stuttgart, 1983.

L. Katchanov, , Éditions Mir, Moscou, 1975.

, . Publ. Brasil Math. Soc., Sao Paulo, 1967.

V. Lamothe, Symmetry group analysis of an ideal plastic flow. , 53, 3, 033704, (2012).

V. Lamothe, Group analysis of an ideal plasticity model, , 45, 285203, (2012).

, , University of California Press, 1964.

R. Mises, , Acad. Press, New-York, 1958.

, , Springer-Verlag, New-York, 1986.

, , Hankbook of Lie group analysis, Ed. N. Ibragimov, vol.3 chapt., 11 CRC Press, London, 1995

Z. Peradzynski, Geometry of interactions of Riemann waves, , Ed. Debnath, Research notes in mathematics no 111, Pitman Avanced Publ. Boston, 1985.

B.L. Rozdestvenski, N.N. Janenko, Systems of Quasilinear Equations and their Applications to Gas Dynamics, , vol. 55, AMS Providence 1983.

, Functionally invariant solutions of wave equations, ,5, (1934), 259-264

, , John-willey Publ, New-York, 1974.

, , in serie Advances of Modern Mathematics, 1998.

[^1]: grundlan@crm.umontreal.ca

[^2]: lamothe@crm.umontreal.ca
