{'fraction_non_alphanumeric': 0.039691955039913136, 'fraction_numerical': 0.011038947982225081, 'mean_word_length': 4.630476621759311, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 17, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Analytical solution method for rheological problems of solids

In classical continuum theory, Volterra's principle [1, 2] is a long-known method to solve linear rheological (viscoelastic) problems derived from the corresponding elastic ones. Here, we introduce and present another approach that is simpler to apply (no operator inverse is required to compute but only linear ordinary differential equations to solve). Our method starts with the known elastic solution, replaces the elasticity coefficients with time dependent functions, derives differential equations on them, and determines the solution corresponding to the initial conditions. We present several examples solved via this new method, like tunnels and spherical hollows opened in various initial stress states, and pressurizing of thick-walled tubes and spherical tanks. These examples are useful for applications and, in parallel, are suitable for testing and validating numerical methods of various kinds.


Motivation and introduction
A large variety of solid materials-like plastics, rocks, asphalt, biomaterials etc. -possess viscoelastic/rheological characteristics. Correspondingly, one can observe some kind of delayed and damped elastic behaviour.
Rheological behaviour of solid media is well-known, e.g., in civil engineering and in mine industry. A hollow opened in an underground stone block often takes its eventual shape only years after the drilling. The diagrams in Figure 1 demonstrate this. Accordingly, rheological behaviour must be taken into account when designing technical devices and facilities. Nevertheless, this means not only disadvantages and problems but also benefits. For example, one can rely on its effect of damping and absorbing vibrations. Many biology-originated objects and protheses used in medical technology also show rheological properties. For instance, rheological behaviour of knee ligaments is apparent. Similarly, a freshly opened underground tunnel needs initially only a temporary-relatively weak-support, and the eventual support is enough to be established only later, when most of displacements have already been occured, allowing thus a much cheaper eventual support.
However, all these need reliable calculations. Nowadays, the most often applied methods are the discretization based numerical methods, which face at problems. Solutions obtained by such methods may considerably depend on the resolution of the applied discretization. Furthermore, for complex three-dimensional problems calculation times are large. Therefore, analytical solution of a simplified version of the problem may provide a reasonable first approximation and give useful insight. Analytical solutions can also be utilized for validating numerical methods.
Here, we introduce and present an exact analytical method for solving linear rheological problems of solids. The approach is based on the corresponding elastic solutions assumed to be already known: the elasticity coefficients are replaced with time dependent functions, which are determined from the rheological equations. The method has been born in a conceptually simple form, with limited range of applicability, and has been enhanced and generalized subsequently in two further steps. Here, we present these three stages in order of increasing generality. For each stage, we show several examples, which not only illustrate the method but also demonstrate its power and limitations.

Elasticity and rheology
We are going to treat purely mechanical problems of homogeneous and isotropic continuous media, and our aim is to determine the displacement field u, the strain field ε and the stress field σ (where both ε and σ are symmetric tensors). We wish to work in the force equilibrial approximation, i.e., when acceleration is neglected: 1 σ · ← ∇ = − g, (2.1) where is mass density and g is volumetric force density (assumed to be time independent), and ← ∇ and → ∇ are the nabla operators acting to the left and to the right, respectively (reflecting proper tensorial order, see also below).
Concerning ε, we stay in the small-strain approximation, which then imposes the geometric compatibility equation in the form → ∇ × ε × ← ∇ = 0. (2.2) According to mathematics, to a symmetric tensor field ε with property (2.2), there exists a vector field u Cauchycalled hereafter Cauchy vector potential [4]-from which ε can be obtained as where Sym denotes the symmetric part of a tensor. The Cauchy vector potential is not unique for a given ε, and all Cauchy vector potentials can be derived from the strain field according to Cesàro's formula [5], u Cauchy (t, r) = u 0 (t) + Ω(t) (r − r 0 ) + r r0 ε(t,r) + 2 ε(t,r) ⊗ ← ∇ A1,3 (r −r) dr (2.4) with A1,3 denoting antisymmetrization in the first and third indices, where the position vector r 0 , the path of integration, the vector function u 0 (t) and the antisymmetric tensor function Ω(t) are each arbitrary. The displacement field u is one of these Cauchy vector potentials so when we wish to reconstruct u from the strain field then we need to fix these uncertainties using symmetry arguments and other physically plausible considerations.
In case of linear elasticity (for a homogeneous and isotropic medium), connection between stress and strain is provided by Hooke's law, (2.5) in the deviatoric-spherical separation, where σ s = 1 3 (tr σ) 1 denotes the spherical part-which is proportional to the identity tensor 1-, while σ d = σ − σ s is the deviatoric (traceless) part; furthermore, E d = 2G is the deviatoric elasticity coefficient and E s = 3K is the spherical one 2 .
For linear rheological models of solids, one can generalize Hooke's law by replacing the elasticity coefficients with polynomials of the time derivative operator. Namely, S s σ s = E s ε s , (2.6) where the stress related operators S d , S s and the strain related ones E d , E s are
In the case of elasticity, equations (2.1), (2.2) and (2.5) form the system of equations to be solved. Together with appropriate boundary conditions imposed at the boundaries of the spatial domain considered, the solution exists and is unique, however, to obtain this solution is not necessarily simple, since (2.5) poses separate conditions for the deviatoric and spherical parts, while (2.1), (2.2) and the boundary conditions prescribe requirements for the sum of the deviatoric and spherical parts.
When we deal with the above-described rheological generalization of the problem then, in addition to (2.1), (2.2), (2.6) and the boundary conditions (which may be time dependent in general), initial conditions are also required to ensure uniqueness of the solution, since the constitutive equations contain time derivatives [see (2.7) and (2.8)]. In the rheological case, all fields are functions of both time and space, and all equations and boundary conditions have to be satisfied for all time instants, which raises an even more complicated task than for the elastic counterpart. It would considerably simplify the situation if one could utilize the known space dependence of the corresponding elastic problem, leaving only time dependence to address.
Volterra's principle [1,2] provides such an opportunity, according to which principle the constants E d , E s in the elastic solution are to be replaced with the rheological operators S d , S s , E d , E s , and solving the resulting temporal equations leads to the solution of the rheological problem. However, in some cases the application of Volterra's principle is difficult, e.g., in cases of time dependent stress boundary conditions, in addition to the fact that, in Volterra's approach, typically one also has to invert operators.
Below, we present another route, which is also motivated by Volterra's idea to treat the space dependence aspect of the problem via utilizing the known solution of the corresponding elastic problem, but is technically easier to follow since only a set of ordinary differential equations is to be solved for the time dependence aspect.

The analytical solution method for the rheological problem
As the first step, let us separate the effect of the force density by subtracting some such time independent fields σ andε -henceforth: primary fields-that satisfy the equations [note that, for time independent stress and strain fields, (2.6) gets simplified to Hooke's law with Then the difference fields-henceforth: complementary fields- 2 G is the shear modulus and K is the bulk modulus. 3 The small-strain assumption allows to approximate the substantial time derivative with the partial time derivative.
satisfy the homogeneous equationsσ Naturally, in general, this transformation modifies the boundary conditions, which is to be taken into account during the calculations.
If the spatial domain filled with the medium has more than one boundary (more than one connected boundary surface) then the problem can be divided into subproblems in which only one boundary condition is nonzero (the boundary condition is nonzero only on one connected boundary surface). Henceforth, we always analyze one such subproblem. Thanks to linearity of all equations involved, the sum of such subsolutions provides solution for a whole problem.
In this work, we consider problems with time dependent boundary conditions that are prescribed for stress, not for displacement. The applications discussed here will all be related to such boundary conditions. Time dependence of boundary conditions will be allowed with the limitation that time dependence must mean a space independent rescaling of the boundary condition-like gradual loading of a surface where loading may be space dependent but the ratio of normal stress values at two different boundary points is time independent. 4 In notation, time dependence of the boundary condition is of the form of a time dependent multiplier λ(t). This λ(t) can be quite arbitrary, the only restriction being that it be sufficiently many times differentiable. For gradual switching on, like when modelling drilling, this factor can be chosen as Corresponding to such a time dependent homogeneous rescaling of the boundary condition, the solution of the elastic problem also gets rescaled-space independently rescaled-by the factor λ(t) .
Another remarkable property of the elastic solution (at any fixed t) is that, based on dimensional reasoning, dependence of the stress solution on E d and E s must be such that stress depends only on the dimensionless ratio In other words, it depends only on the Poisson's ratio This property is, naturally, apparently visible in the examples considered below. It is to be emphasized that, although usually one focuses only on the space dependence of an elastic solution, for the solution methods described here, dependence on the elasticity coefficients E d , E s will also be of central importance.

First approach: Method of elasticity constants made time dependent
In the first-the simplest-version of our approach, we adapt Volterra's principle in the form that the rheological solution is obtained from the elastic one by replacing the elasticity coefficients E d , E s with time dependent functions (rather than with operators, as in Volterra's methodology). At any instant t, for any E d (t), E s (t) the solution is a valid elastic solution (for the current boundary condition) so (3.5) and (3.6)-as well as the boundary condition-are satisfied.
Apparently, the spatial equations-(3.5), (3.6) and the stress boundary condition-will be satisfied at any time instant, with the actual values E d (t), E s (t). The remaining equations-namely, the rheological ones (3.7)-will yield ordinary differential equations on E d (t) and E s (t).
One can observe that one part of the spatial conditions, (3.5) and the stress boundary condition, refers only to the stress field while the remaining spatial condition, (3.6), refers only to the strain field. Accordingly, one is allowed to detune the elastic solution for stress from the elastic solution for strain: the elastic solution for stress can be utilized with some E d σ (t), E s σ (t) while the elastic solution for strain can contain some separate E d ε (t), E s ε (t). Notably, the only condition forbidding this detuning would be (2.5) but for rheology it is replaced with (3.7) so it is not excluded that some consistent solution can be found.
As can be seen in the examples below, consistent solutions are indeed possible, either with this detuning or even without it.

Cylindrical bore (tunnel) opened in infinite, homogeneous and isotropic stress field
In an infinite, homogeneous and isotropic stress fieldσ, we open an infinite cylindrical bore with radius R (see Figure 2). In cylindrical coordinates, the boundary conditions specifying the elastic solution for the completely open bore are σ rr (R, ϕ, z) = 0, lim r→∞ σ(r, ϕ, z) =σ, (3.11) which are rewritten for the complementary field aŝ The solution of the elastic problem, for this completely opened bore, iŝ Notice that the spherical part of this tensor is zero,σ s = 0, therefore,σ =σ d , from which the strain tensor iŝ Now let us consider the problem with time dependent boundary condition, e.g., with a λ(t) of the form (3.8) via which we can model the drilling process. The corresponding elastic solution is nothing but the previous one rescaled by λ(t):σ For a rheological problem (3.7), the above-described method of time dependent elastic coefficients says to substitute the only elasticity coefficient E d present in the previous equations for an unknown time dependent function E d (t):σ Obviously, the spatial conditions are satisfied at any time t. The only thing left to do is to impose (3.7). The spherical equation is trivially fulfilled, while is generated for the deviatoric part. We have to solve this equation with such initial conditions that the complementary stress and strain fields before the drilling (before the time dependent change in the boundary condition) are zero for an extended time interval-implying that all time derivatives of the fields are also zero. Then the solution of (3.19) is unique so the method has reached the goal. Comparing the outcome found here with the one known in the literature [11], obtained via another approach, we find that the two results are in complete agreement.

The rheological process of a spherical hollow in an infinite, homogeneous and isotropic stress field
The boundary conditions for a completely opened spherical hollow of radius R opened in infinite, homogeneous and isotropic stress fieldσ (see Figure 3) are, in spherical coordinates: The elastic stress solution isσ The spherical part of this tensor is zero again soσ s = 0,σ d =σ, hence, the strain tensor iŝ Also analogously to the previous, cylindrical, case, the solution of the elastic problem for the gradually opened hollow isσ el,λ (t, r) =σ d el,λ (t, r) = λ(t) ·σ d el (r), (3.24) (3.25) 5 The same holds for all the subsequent problems, too: For plots obtained for concrete rheological models, see Section 4.
For rheology, we substitute the only elasticity coefficient E d for an unknown time dependent function: When the rheological operators act on these functions then the spherical equation is trivially fulfilled, while follows for the deviatoric part. With initial conditions as in the previous example-zero initial history-the solution exists and is unique. Comparing (3.28) with (3.19) shows that, from the point of view of our method, these two problems lead to the same rheological equation (there are only spatial differences between the two problems). The same similarity is the reason why we treat the two next problems, pressurizing of a thick-walled tube and of a spherical tank, together.

Pressurizing of a thick-walled tube and of a spherical tank
Our next examples are the pressurizing of a thick-walled tube and of a spherical tank from zero overpressure to overpressure p 0 . Effect of of possible initial pressure can be subtracted by the primary fields so we formulate the problem directly for the complementary fields.
The boundary conditions for a thick-walled tube with the inner and outer radia R 1 , R 2 at overpressure p 0 areσ The elastic stress solution in the wall iŝ which can separated into a deviatoric and a spherical part aŝ Hence, the strain tensor isε Similarly, the boundary conditions for a spherical tank with inner and outer radia R 1 , R 2 arê The elastic stress solution isσ the deviatoric and spherical parts of which arê Comparing the equations (3.32) and (3.36) shows that our method of time dependent elastic coefficients will lead to the same temporal ordinary differential equations. When we model the gradual pressurizing of the tube/tank, we multiply the inner boundary condition by λ(t). The corresponding elastic solution iŝ For the rheological solution, the elasticity coefficients E d and E s are changed to time dependent functions: (3.40) From this ansatz, the rheological operators generate the equations which can be solved for the two unknown functions E d (t), E s (t) (or, more conveniently, for λ(t) E d (t) and λ(t) E s (t) ). Again, initial conditions are taken from that, for t < t 1 , λ(t) E d (t) = 0 and λ(t) E s (t) = 0 .

Cylindrical bore (tunnel) opened in an infinite homogeneous but anisotropic stress field
Now let us consider an infinite and homogeneous, but now anisotropic, stress field, and let us analyize the rheological process caused by drilling a cylindrical bore. This problem is a generalization of our first example (cylindrical bore/tunnel opened in infinite, homogeneous and isotropic stress field, Subsection 3.1.1 with Figure 2). The solution of the elastic problem for fully opened bore can be taken from [12], 6 and is the sum of two terms containing linearly independent space dependent functions, -so to say, 'spatial patterns'-, the first term depending on η and the other being independent of it: where we are using the following notations related to the primary fieldσ: [cf. (3.9)]. Drilling is modelled again via multiplying the final boundary condition by a factor λ(t); then the solution of the elastic problem gets multiplied by the same λ(t): Following the already known recipe, let us replace the elasticity coefficients in the stress and strain solution with unknown time dependent functions. At this point we can realize that it is possible (and will indeed be necessary) to substitute two separate function pairs E d [recall the argument in Section 3.1]. Accordingly, let us replace η in the stress solution with , (3.47) while in the strain solution we use some separate time dependent function , (3.48) in addition to changing in (3.46) the explicite E s to E s ε (t). Since stress depends only on the ratio of E d σ (t) and E s σ (t) [on η σ (t) solely], there is only a one-function freedom in stress. To fix the arbitrariness, we take the simplest choice E s We then have three functions to be determined. The rheological operators generate three conditions on them, corresponding to the fact that the strain solution (3.46) contains three linearly independent tensor fields: For the time dependent coefficient of each of these independent tensor fields, one equation is generated; two equations follow from the deviatoric rheological equation and one from the spherical one [cf. (3.7)]. These three equations read (3.49) At first sight, solving this system of equations for η σ , η ε , E s ε seems difficult. Nevertheless, introducing the auxiliary functions and realizing the relationship one arrives at a system of linear differential equations, This system is solvable, and has a unique solution with, for example, the assumption of 'zero past' (zero complementary fields before the opening, for a whole time interval). 7 At this point, we can see that this method-the method of elasticity constants made time dependent-is limited to three unknown time dependent functions.

Second approach: Method of elastic spatial patterns
Here, we establish a method that enables more than three unknown time dependent functions to be determined. We assume that the solution of the elastic problem is of the form where J is some integer, the coefficient functions c j (η) are linearly independent of each other, and the spatial patterns s j (r) are also linearly independent. Notably, one could allow the sum to be infinite; however, expansion with respect to an infinite function series raises convergence questions, and here we wish to avoid such mathematical complications. Fortunately, the finite sum form already allows us to treat a good number of special cases, as shown below.
In deviatoric-spherical separation, (3.53) readŝ Dimensional reasoning suggests to use, instead of strain, a multiple of it that has the dimension of stress. This can be simply achieved by 8 ζ := E s ε, (3.55) called hereafter stress-dimensioned strain. Correspondingly, The compatibility equation (3.6) remains in the same form forζ: Withζ, Hooke's law is simplified toσ Then it is apparent that, in a solution of an elastic problem with stress boundary condition,ζ depends on the elasticity coefficients only through η, too [see (3.9)]. For the elastic solution, Hooke's law (2.5) leads tô When the stress boundary condition is multiplied by a time dependent factor λ(t), the corresponding elastic stress and stress-dimensioned-strain solutions are multiplied accordingly: The rheological problem characterized by (3.7) imposes (recall that, as seen in Sections 2 and 3, E s 0 plays the role of E s in rheology). As in the previous method, we look for the rheological solution as a time parametrized succession of elastic solutions, in order to satisfy the spatial conditions (3.5), (3.57) and the stress boundary condition. However, while previously η was replaced with time dependent functions, now let us fix J constants η 1 , η 2 , . . . , η k , . . . , η J , consider the corresponding elastic solutions, and take a time dependent linear combination of them as our new ansatz. We allow for separate time dependent coefficients forσ andζ: The explanation for why to use exactly J elastic solutions is that we have J linearly independent spatial patterns, the vector space of spatial patterns is J dimensional, so using more than J pieces of combinations would be redundant, would not increase our freedom. Thanks to working from elastic solutions, the spatial equations (3.5) and (3.57) are satisfied, while satisfying the boundary condition is ensured via the condition What is left is only to fulfil the rheological relationships (3.63). From (3.64)-(3.65), one finds An important advantage of this second method is that, apparently, the system of differential equations to be solved for the functions λ k (t), κ k (t) is linear. The question arises whether the solution depends on the choice of the values η 1 , . . . , η J . Now, the J elastic solutions form a basis 9 in the space of the linear combinations of the J spatial patterns. Therefore, another set of values η 1 , . . . , η J represents another basis, corresponding to which the linear expansion coefficients of the solution are some other functions λ k (t), κ k (t). The solution itself is the same (assuming the appropriate amount of initial conditions, naturally).
To simplify the notations, let us introduce the matrix of constant elements. This matrix is nondegenerate-except probably for certain special choices of η 1 , . . . , η J -, since the functions c j (η) are also linearly independent. Formulae (3.67)-(3.68) can be rewritten as The spatial patterns s 1 (r), . . . , s J (r) are linearly independent functions. If their deviatoric parts are also linearly independent and their spherical parts are also linearly independent, (3.70)-(3.71) require equality of the 9 Apart probably from certain special choices of η 1 , . . . , η J .
coefficients of each component: These are 2J equations, and taking (3.66) also into account, altogether we have 2J + 1 equations for the 2J unknown functions λ 1 (t), . . . , λ J (t), κ 1 (t), . . . , κ J (t). Hence, in general, this second method cannot provide a solution. Hovewer, we can observe in each of the problems analyzed subsequently that the J spherical patterns s s 1 (r), . . . , s s J (r) are not linearly independent but one of them can be expressed as a linear comination of the others. Then (3.71) means only J − 1 independent equations and the method can lead to a solution. If the initial conditions can also be written in the form of (3.64)-(3.65)e.g., the 'zero past history' assumptionis of this form, corresponding to λ = 0 -, then the found solution is the solution of the problem.
Let us now see some examples how this method works and performs in practice. The first three problems have already been solved with the method of time dependent elastic coefficients, too, so comparison can also be made.

Cylindrical bore (tunnel) and spherical hollow opened in infinite, homogeneous and isotropic stress field
As the first example, let us consider a cylindrical bore and a spherical hollow opened in an infinite, homogeneous and isotropic stress field. As happened in Subsections 3.1.1-3.1.2, these two problems are expected to lead to the same rheological ordinary differential equations to be solved.
The elastic stress solutions are in this case traceless-pure deviatoric-tensors [see Subsections 3.1.1-3.1.2], which can be written, in the notation of (3.53), aŝ where c(η) = 1 and s(r) =σ(r), while the stress-dimensioned strain iŝ The elastic solution of both problems, with time dependent boundary conditions, iŝ According to the second method, the rheological solution is looked for in the form of Comparing this with (3.19) and (3.28) shows that we have reached the same equations.

Cylindrical bore (tunnel) opened in infinite, homogeneous but anisotropic stress field
The elastic solution of this problem, (3.42), can be written in the form of (3.53) aŝ where c 1 (η) = 1, and the notations forσ can be found in (3.43).
Since we have two independent spatial patterns, we are looking for two functions κ k (t) and two functions λ k (t), and the introduced matrix C jk in (3.69) reads here Consequently, the rheological solution is looked for, based on (3.64)-(3.65), in the form For the four unknowns, besides these three equations, the boundary condition related (3.66) gives the fourth equation: Solutions for concrete rheological models can be found in Section 4.

Cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight
As our most complicated example, now we discuss the rheological process of a semi-infinite domain loaded by its own weight and weakened by a bore. To describe the rheological process of an underground tunnel, this model can be considered as more accurate than the previous bore models. The elastic solution of the problem can be found in [13,14,15]. The primary stress field-the one before the drilling-, in an appropriate Cartesian coordinate system (see Figure 4), can be written as where γ = g describes the homogeneous force density, d is the depth of the center of the bore from the surface and the parameter k is called lateral pressure factor. [13] gives the solution of the problem for three different values of k: • When we assume hydrostatic pressure distribution for the primary field then k = 1. This is a good approximation for tunnels opened at large depths.
• When one can assume that the dilatation of the medium is laterally inhibited then the strain components ε xx andε zz are zeros so one can derive for the lateral pressure factor k = ν 1−ν = 1−η 1+2η . • When one can assume that the dilatation of the medium is laterally free then the stress componentsσ xx andσ zz are zeros so k = 0.
Transforming (3.87) to cylindrical coordinate system yields the stress components The boundary conditions are prescribed for the contour of the cylinder and for the plane surface-the horizontal boundary-after the drilling; on these boundaries the normal component of stress is zero.
Mindlin gives the solution in form of an infinite series in bipolar coordinate system -which suits to both the cylinder and the plane [13]. If the ratio of the depth d of the center of the bore from the surface and the radius R of the cylinder satisfies d R > 1.5-large-depth approximation-, then it suffices to take the leading order term from the infinite series.
This first term is transformed to cylindrical coordinate system in [14,15]: The entries of the complementary field (having a plane strain situation) arê Now we analyze the rheological process caused by the drilling for the above-mentioned three values of the lateral pressure factor.

Free lateral deformations (k = 0)
In this case, the complementary stress field can be represented again as the sum of three independent spatial patterns:σ where and The matrix C jk is (3.120) Now noticing that s s 3 = 2 (s s 1 − s s 2 ), there is again one spherical equation less than deviatoric: (3.122) Finally, from the boundary condition, we have (3.123)

Pressurizing of a thick-walled tube and a spherical tank?
The The rheological operators generate the set of equations which is an overdetermined system of equations [two equations for the only unknown κ(t)] so in this case the second method cannot provide the solution. (Fortunately, the third method coming soon will cover this case as well.)

Conclusions about the first and second approaches
We have seen that the first method is rather limited since only those problems can be solved via the method where the number of linearly independent spatial pattern deviatoric parts, plus the number of linearly independent spatial pattern spherical parts, is not more than three. In the second approach there is no upper bound on the number of independent spatial patterns. However, a limitation is that only the freedom in stress elastic patterns is utilized, (stress-dimensioned) strain is not considered, and the independent deviatoric and spherical equations are not controlled, either. Therefore, such simple problems as the pressurizing of thick-walled tubes or of spherical tanks cannot be solved via the second approach.
These experiences have inspired us to establish a common generalization of the first two methods.

Common generalization: The third approach: method of four elastic spatial pattern sets
Since the second method prescribes separate equations related to the deviatoric and spherical parts of the spatial patterns s j (r), now the idea is to explore the linear independence content within the deviatoric sector and within the spherical one separately. Moreover, we reveal the linear independence content not only for stress but also for stress-dimensioned strain. We assume again that the elastic stress solution can be written in a finite sum form: where I is some integer, a i (η)'s are linearly independent coefficient functions and α i (r)'s are linearly independent spatial patterns (as before). However, now we also express the elastic stress-dimensioned strain solution in an analogous finite sum:ζ When the boundary condition is multiplied by the time dependent factor then the elastic stress solution iŝ σ el,λ (t, η, r) = λ(t)σ el (η, r). (3.135) As in the first two methods, we wish to obtain the solution of the rheological problem as (time dependent) combination of elastic solutions, hence fulfilling the spatial condition (3.5), (3.57) and the stress boundary condition automatically. Let us fix I values η 1 , η 2 , . . . , η m , . . . , η I for the stress solution and J values η 1 , η 2 , . . . , η n , . . . , η J for the stress-dimensioned strain solution, and consider time dependent linear combinations of the corresponding elastic solutions, separate ones for stress and for the stress-dimensioned strain (as before).
In light of (3.131) and (3.132) , our ansatz is noŵ The rheological conditions generate the system of equations Since the spatial patterns δ(r)'s and γ(r)'s are linearly independent sets, the equality of the corresponding coefficients follows: The boundary condition requires Altogether, we have K + L + 1 equations for the I + J unknowns. This is the most general version of our approach. Compared to the first two methods, this third covers both, and it extends the idea probably to the most general level reachable. Namely, while the first two methods utilize the available elastic spatial patterns in two restricted forms, the third formulation explores all freedom enabled by the stress spatial pattern set, the stress-dimensioned strain one, the deviatoric pattern set and the spherical one each. Again, let us next see via examples how to use it.

Cylindrical bore (tunnel) and spherical hollow opened in infinite, homogeneous and isotropic stress field
As we have seen above, these two problems lead to the same rheological equations so we treat them together again. The elastic stress solution of the cylindrical bore case can be written in the form of (3.131) aŝ with I = 1, a 1 (η) = 1 and α 1 (r) =σ el (r). In parallel, the elastic stress solution of the spherical hollow case iŝ (3.150) with I = 1, a 1 (η) = 1 and α 1 (r) =σ el (r) again. The deviatoric part of the elastic stress in the form of (3.133) isσ d el (r) = K k=1 ηc k (η) γ k (r) (3.151) with K = 1, c 1 (η) = 1 η and γ 1 (r) =σ el (r). Since the spherical part of the elastic stress is zero,σ s el (r) = 0, which means that L = 0. Similarly, with these, the stress-dimensioned strain in the form of (3.132) iŝ with J = 1, b 1 (η) = 1 η and β 1 (r) =σ el (r). The matrices A ik and B jk can be read off easily, are 1 × 1 matrices, namely, A ik = B jk = 1. Therefore, subsituting into (3.146), (3.147) and (3.148) provides the system of equations two equations for the two unknown functions ϕ(t), ψ(t).

Solutions of the previous examples with concrete material models
All results of each of the three methods have so far been presented for an arbitrary linear rheological model. In this section, we show the solutions for certain concrete material models. The function that scales the boundary condition is chosen as [cf. (3.8)]: (4.1) According to experience, rheological material behaviour is usually seen in the deviatoric part, hence, in our examples we always assume Hooke's elasticity model for the spherical part. In parallel, we take the Kelvin model and the Kluitenberg-Verhás model for the deviatoric part.
Three cases are analyzed in all of the considered examples: • when the switch-on λ(t) is very slow compared to the rheological time scales, • when the switch-on time scale t 2 − t 1 is comparable to the rheological time scales, • when the switch-on λ(t) is very fast compared to the rheological time scales.
In the figures below, the functions indicate the time dependent factors of the deviatoric spatial patterns of stress and of stress-dimensioned strain; and the analogous time dependent coefficients of the spherical patterns are Λ s j (t) and K s j (t); corresponding colour codes are defined in Figure 5 (in cases when less functions are enough then the rest are omitted). 10 Κ 2 s (t) Figure 5: Left: colour codes of the coefficient functions belonging to the spatial patterns of stress. Right: colour codes of the coefficient functions belonging to the spatial patterns of stress-dimensioned strain.
Furthermore, in plots obtained by the third method, the notations are used.
The following figures show the time evolution of the patterns belonging to stress on the left the time evolution of the patterns belonging to the stress-dimensioned strain is shown on the right. The cases of slow, medium and fast changes in the boundary condition are displayed in the order up-to-down.

Kelvin -Hooke model
First consider the simplest rheological material model, Kelvin model in the deviatoric part and Hooke model in the spherical part: The rheological time scale of the model is The simplest examples are the cylindrical bore (tunnel) and spherical hollow opened in homogeneous and isotropic stress field, which lead to the same deviatoric equation (3.80), whose solution is plotted in Figure 6. Due to the geometrical simplicity of the setting, we cannot observe any remarkable phenomenon but it is to be noted that the only coefficient function that belongs to stress is equal to the coefficient function λ(t) of the boundary condition, and strain has an observably slower increase so it can be easily recognized that the medium gets deformed even after the end of the drilling.
As this problem is solvable via the first method, too, in Figure 7 we plot the corresponding time dependent function that replaces the static deviatoric elasticity coefficient. As expected, the ratio E d (t)/E s converges to the static-Hookean-value, E d /E s = η = 0.4.
For the example of pressurizing of a thick-walled tube and of a spherical tank, we find one deviatoric coefficient function and two spherical ones ( Figure 8). Correspondingly, in the first method, the deviatoric and spherical elasticity constants are replaced with separate time dependent functions ( Figure 9). Here, the deviatoric ratio tends again to the Hookean value η = 0.4, as it should.
In the case of the cylindrical bore (tunnel) opened in homogeneous but anisotropic stress field, the situation is the opposite: One has two deviatoric coefficient functions and a single spherical one, shown in Figure 10. In this case one can observe an interesting phenomenon: as long as the boundary condition changes, there are strong transients in the solution, and some of them are visible even after the end of the change in the boundary condition. Also, this is a case where, according to the first method, altogether three time dependent functions replace the elasticity constants in stress and in strain ( Figure 11). The two deviatoric ratios, E d σ (t)/E s and E d ε (t)/E s , tend to η = 0.4 and the spherical one, E s ε (t)/E s to 1, following the expectations. Next, we present the solutions of the cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight field. The solutions for the hydrostatic initial stress state (k = 1) can be seen in Figure 12. Here it is also observable that one of the coefficient functions starts with anomalous sign-with opposite sign compared to its eventual sign in the large-time asymptotics-. Consequently, rheology means not only damping and delay but a more complicated mechanical behaviour.
When the primary field does not allow lateral deformations (k = ν 1−ν = 1−η 1+2η ) then there are three independent deviatoric and two independent spherical patterns so we seek for five-five functions; these can be seen in Figure 13.
Finally considered is the case of free lateral deformations (k = 0). The solutions are plotted in Figure 14. Time dependence is more complicated here, three of the five functions starting with anomalous sign.
The question arises what it may cause when the coefficient function of a pattern starts with anomalous sign. We give the answer for this question later, during the analysis of the displacement field.     Figure 9: Time dependence at the place of the elasticity coefficients (orange: deviatoric, red: spherical) according to the first method, for the example of pressurizing of a thick-walled tube and of a spherical tank. In black, the gradual opening λ(t) is also displayed. Up to down: slow, medium and fast changes in the boundary condition.   Figure 11: Time dependence at the place of the elasticity coefficients according to the first method, for the example of a cylindrical bore (tunnel) opened in homogeneous but anisotropic stress field. Black: the gradual opening λ, green: Up to down: slow, medium and fast changes in the boundary condition.   We consider the example of the cylindrical bore (tunnel) opened in infinite, homogeneous but anisotropic stress field. In case of weak oscillation the time evolution of the coefficient functions are plotted in Figure 15, while in case of strong oscillation in Figure 16. In the previous case the fast change in the boundary condition causes stronger transients than it was by the Kelvin model, while in the latter case the solutions show strong oscillations, as we have expected.
We also present the solution of a cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight when no lateral deformations are allowed in the primary field (k = ν 1−ν = 1−η 1+2η ). In case of stongly damped oscillations the time evolution of the coefficient functions are shown in Figure 17, and the case of weakly damped oscillations are plotted in Figure 18. In the previous case one can see weak oscillation in the initial transients, while in the latter case the oscillation can be observed on each coefficient functions.
Since materials with positive index of inertia are currently not known, the corresponding displacement fields are not presented.

Displacement fields
The displacement field 11 can be given via the Cesàro formula (2.4). One needs to substitute into (2.4)ε = 1 E sζ based on (3.56), and then (3.65) or (3.137). Since stress-dimensioned strain is given in finite sum form and, 11 Understood, naturally, with respect to the primary initial state of the continuum. moreover, spatial and time dependences are separated, one can obtain u Cauchy (t, r) = u 0 (t) + Ω(t) (r − r 0 ) + (r −r) dr + (4.7) (r −r) dr + (4.8) where Ψ d k , Ψ s l are defined in (4.4). One should not forget that these are only Cauchy vector potentials so the rigid body like displacement and rotation have to be fixed.
In case of the cylindrical bore (tunnel) opened in infinite, homogeneous and isotropic or anisotropic stress field, these uncertainties can be completely fixed by integrating from the centre and by subtracting the rotation in the infinity.
To fix these uncertainties for a cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight is not so easy. The elastic solution, which is used to derive the rheological solution, is only the first-term approximation of an infinite series so our rheological solution is just an approximation, valid only in a certain neighbourhood of the bore. The Cauchy vector potentials that originate from these solutions diverge at infinity so fixing the uncertainties at infinity cannot be used. 12 Instead of this, at each time instant, we determine the point of the surface where the tangent is zero, and subtract the vertical displacement of this point from the displacement field, which procedure turns out to provide reasonably realistic result (see the figures).
In the example of a cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight, we have chosen the dimensionless ratios d/R = 2, γ/E s = 0.02 (4.9) for the calculated results shown here. For the cylindrical bore (tunnel) opened in homogeneous and isotropic stress field, the used stress componentσ rr has been adjusted correspondingly:σ rr = −γd. Similarly, for the problem of the cylindrical bore (tunnel) opened in homogeneous but anisotropic stress field we took the stress componentσ yy = −γd, and the stress componentσ xx was be calculated fromσ yy via the lateral pressure factor k, applying the value k = ν 1−ν = 1−η 1+2η , which describes the case when lateral deformations are not allowed in the primary field. We choseσ xy = −γd, and zeros for the other stress components because we analyze the problem perpendicular to the axis of the bore. Figures 19-23 show the time evolution of the displacement field for the cases of slow, medium and fast changes in the boundary condition. Time goes by from the blue line towards the red line; the contour of the bore and the surface are plotted in the nondimensional time instants 0, 0.5, 1, 2, 4 and 8. In the figures, the displacement of the contours has been artificially enlarged since the found small-strain displacements themselves would be too small for the eye, and our primary goal here is the visualisation of the tendencies.
Generally, it can be stated that the shape of the bore changes weakly at the beginning of the process for slow changes in the boundary condition, however, when the boundary condition changes fast then the bore rapidly deforms at the beginning of the process. According to our expectations, the elastic final (asymptotic) state is independent of the changes of the boundary condition. Figure 19 shows the time evolution of the displacement field of a cylindrical bore (tunnel) opened in homogeneous and isotropic stress field. The bore is shrunk as time goes by and, at the same time, the surface is displaced in the direction of the center of the bore. Figure 20 shows the time evolution of the displacement field of a cylindrical bore (tunnel) opened in homogeneous but anisotropic stress field. Although this model does not take into account the load of the self weight of the medium, it is possible to calculate the rotation of the cross section using this approximation. Other interesting observations can also be made: we can notice that the geometry is stretched in the direction of the major axis of the rotated ellipse and shrinks only later. Comparing these with Figure 10, we can see that one of the coefficient functions starts with anomalous sign at the beginning of the process so the pattern that belongs to this function initially expands and shrinks only later. Figures 21-23 show the time evolution of the displacement field belonging to the cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight. This model approximates the processes around tunnels the best. One can observe in each figure that the center of the bore moves towards the direction of the surface. In case of hydrostatic initial stress state (Figure 12), the contour of the bore deforms approximately isotropically, and the sinking of the surface is negligible. When lateral deformations in the primary field are not allowed (Figure 22), and in the opposite case ( Figure 23) the sinking of the surface and a swelling-like phenomenon are also observable. Figure 19: Time evolution tendency of the displacement field of a cylindrical bore (tunnel) opened in homogeneous and isotropic stress field. Up to down: slow, medium and fast changes in the boundary condition. Figure 20: Time evolution tendency of the displacement field of a cylindrical bore (tunnel) opened in homogeneous but anisotropic stress field. Up to down: slow, medium and fast changes in the boundary condition. Figure 21: Time evolution tendency of the displacement field of a cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight, hydrostatic initial stress state (k = 1). Up to down: slow, medium and fast changes in the boundary condition. Figure 22: Time evolution tendency of the displacement field of a cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight, no lateral deformations allowed in the primary field (k = ν 1−ν = 1−η 1+2η ). Up to down: slow, medium and fast changes in the boundary condition. Figure 23: Time evolution tendency of the displacement field of a cylindrical bore (tunnel) opened in homogeneous medium loaded by its self weight, free lateral deformations in the primary field (k = 0). Up to down: slow, medium and fast changes in the boundary condition.

Conclusions and outlook
The method presented here (in its most general, third, form) is applicable to many rheological problems. It can be helpful e.g., during the process of designing of tunnels and underground facilities as an insightful approximation, but can also be useful for validating numerical solvers.
As a heuristic summation of the range of problems that can be covered by this method is that problems that are simple enough to be treatable via reasonable analytic means can probably be solved this way.
Mathematically, a limitation of the method is that initial conditions must be expressible in terms of the utilized elastic spatial patterns. Fortunately, the vast majority of practical problems-those when one assumes a stationary, undisturbed, equilibrial initial state for a time interval-is of this kind.
The development of the method can continue in the future in several directions. First, it would be instructive and useful to compare obtained results with ones stemming from other methods (e.g., finite element method).
Second, in each of the considered examples, the process was a plane-strain process. In this case, the form of the strain tensor isε This seems to be the reason why in the considered examples we always happen to find a relationship among the spherical parts of the independent spatial patterns, which reduces the number of independent rheological equations so one arrives at a well-determined system of linear ordinary differential equations. Our conjecture is that the method can be applied for all plane problems. This statement has not proven yet but in all the considered examples this property can be seen so it is plausible to search for a proof of it. When there is no relationship among the spherical patterns then the method results in an overdetermined differential equation system. In this case, a generalized, approximate, version of the method could be used which would provide an optimized approximate solution. In this context, variational formulation of the rheological problem may be fruitful.
Finally, the method could be generalized to cases when the elastic stress solution can be given only in infinite sum form (see e.g., [13]), or the recent remarkable, tunnel related, analytical solution by Cai and coworkers [16]). The approximation of the infinitely many terms with a finite sum can lead to a novel numerical method that may be capable to treat situations with complex geometries. For this purpose, a variational formulation can also come helpful.