A Newmark space-time formulation in structural dynamics

In this contribution, we present a space-time formulation of the Newmark integration scheme for linear damped structures under both harmonic and transient excitations. The incremental set of equations of motion and the Newmark approximations are transformed into their corresponding space-time equivalents. The dynamic system is then represented by one algebraic space-time equation only. This equation is projected into a coupled pair of space-time equations, which is solved via the fixed point algorithm. The solution is iteratively assembled by enrichments, each of which is decomposed by a dyadic product of spatial and temporal enrichment vectors. The evolution of the spatial enrichment vectors is investigated during convergence and interpreted by comparing them to the set of linear modes of vibration. The new method is demonstrated by means of four numerical examples, presenting not only the excellent convergence behavior and the numerical efficiency but also the limits of the proposed approach.


Introduction
Evaluating the response time history of engineering structures subjected to transient excitation has become a key issue in structural dynamics. The wide range of conventional methods includes analytic strategies for linear systems [1] and step-by-step time integration schemes [2,3]. Especially with emphasis on computational efficiency, several approximation schemes and model order reduction strategies have been proposed in the literature [4][5][6][7].
As a statistical pattern technique in incremental structural dynamics, the proper orthogonal decomposition is used to find the dominant motion patterns for an optimal low-order description of linear and nonlinear systems [8]. Here, the transformation into low-order subspaces depends on the set of proper orthogonal modes that are obtained from parts of B Franz Bamer Bamer@iam.rwth-aachen.de Nima Shirafkan Shirafkan@iam.rwth-aachen.de 1 Institute of General Mechanics, RWTH Aachen University, Eilfschornsteinstraße 18, 52062 Aachen, Germany 2 Univ. Lille, CNRS, Centrale Lille, UMR 9013 -LaMcube -Laboratoire de Mécanique, Multiphysique, Multiéchelle, F-59000 Lille, France the structural response, the so called observation time period [9]. The applicability of the proper orthogonal decomposition to transient excitations mainly depends on the choice of this observation time period [10]. If nonlinearities are not captured in this observation time period, they cannot be reproduced in the low-order model. However, if the nonlinear response behavior is captured, a surprisingly small number of proper orthogonal modes is sufficient to accurately describe the dynamic response behavior of the structure [11,12]. Another approach is to divide the structure of interest into various substructures. The main goal is to provide low-order representations of the inside areas of the substructures while keeping the number of degrees of freedom of the boundaries unaltered. This strategy is pursued both in the Craig-Bampton method [13] and in the dual Craig-Bampton method [14]. Although these two methods are mainly applicable to linear systems, a hierarchical formulation of the Craig-Bampton method has lately been presented for vibroimpact problems with linear contact formulations [15,16] and highly nonlinear contact laws [17]. Recently, also a combination of the proper orthogonal decomposition and substructure techniques has been proposed in the literature [18]. However, the success and efficiency of these strategies is mainly based on the chosen modes or deformation patterns that must be a priori evaluated.

123
Another type of approach is to apply non-incremental structural analysis. In this context, one methodology is to apply finite elements simultaneously in space and time [19]. Using this strategy, different space-time element shapes have been introduced, investigating the numerical stability [20,21]. Furthermore, the application of space-time finite elements in elastodynamics and structural dynamics has been introduced [22][23][24]. It has also been proposed that the space-time finite element approach could be highly effective in countering discontinuous dynamic problems like fracture [25]. Furthermore, space-time finite elements have been applied to contact mechanics [26][27][28], multiscale modeling [29,30], landslide dynamics [31], problems including fluid-structure interaction [32] and modeling of viscoelastic materials [33][34][35]. Lately, an efficient solution algorithm for space-time finite element problems has been suggested that is based on the summation of Kronecker products of temporal and spatial submatrices [36].
The proper generalized decomposition (PGD) is an another strategy that does not depend on a priori response identification [37][38][39], as required for proper orthogonal decomposition based methods [11]. It allows one to consider the whole problem simultaneously in space and time, while using conventional space discretization strategies such as the finite element method. In this way, the method uses an evolving summation of functionals, iteratively evaluated at once within the whole time domain [40]. This solution procedure was firstly introduced by Ladevéze [41]. Subsequently, the PGD was introduced as model order reduction method for multi dimensional problems [42,43]. The strategy has then been applied in many fields such as chemistry [44], thermo-mechanics [45], fluid dynamics [46][47][48], geophysics [49], material science [50], solid and fracture mechanics [51][52][53][54], as well as error estimations and optimization strategies in mechanics [55][56][57]. Especially in structural mechanics, first approaches have been proposed to use the proper generalized decomposition ansatz to solve elastoplastic and viscoplastic structures subjected to quasistatic loading [58,59]. Boucinha et al. gave a first introduction to how to deal with linear undamped dynamic systems using nonincremental formulations [60]. Shirafkan et al. physically interpreted the evolution of spatial and temporal functions of a quasistatic problem, which they refer to as proper generalized decomposition modes, and analyzed these modes during convergence [61]. Additionally, the solution of parametric dynamic structures in the frequency domain has been investigated [62,63]. However, space-time formulations in structural dynamics are far from well understood. In particular, a detailed description of numeric space-time integration approaches including a dynamical interpretation during convergence is, to date, missing from the literature.
In this paper, we present a new space-time solution strategy in structural dynamics. We provide the fundamental basis of the continuous and discrete space-time decomposition, based on which we present the space-time equivalents of the set of equations of motion and the incremental Newmark equations. This leads to a coupled space-time matrix equation pair. In order to solve the coupled space-time problem, we propose an iterative approach based on the fixed point algorithm and investigate the evolution of spatial and temporal modes after each enrichment. This paper aims to provide the basis for future space-time approaches applicable to nonlinear transient structural dynamics.
In Sect. 2, we introduce the problem of the conventional incremental elastodynamics in a nutshell, including the strong form, the weak form, the space discretization, the matrix set of equations of motion and the Newmark equations. In Sect. 3, we present the new Newmark space-time integration strategy, introducing the continuous space-time formulation, space-time discretization, space-time decomposition of the equations of motion and the Newmark equations as well as the implementation of the new strategy. In Sect. 4, the Newmark space-time algorithm is demonstrated on four numerical examples, showing the accuracy, efficiency and limits of the method. Finally, in Sect. 5, the conclusions are drawn.

Problem formulation
We define the initial boundary value problem starting with the linear momentum balance at the material point in a certain body Ω in the time domain of interest [0, T ] [64]: where σ (x, t) is the continuous second order Cauchy stress tensor field, f v denotes the vector field of the volume forces and u(x, t) is the continuous displacement field. The fieldü denotes the second derivative of u(x, t) with respect to time. The scalar field ρ(x) describes the density distribution of the material within the body Ω. Dirichlet boundary conditions are defined on the surface Γ u as: while Neumann boundary conditions are defined on the surface Γ σ as: Furthermore, two necessary conditions hold on the whole boundary ∂Ω of the body Ω: The initial conditions are defined by the initial displacement and initial velocity fields on the whole body at the time instant t = 0: Multiplying Eq. (1) by a weight function that obeys the kinematic boundary conditions and integrating over the body Ω, leads to the weak form of the initial boundary value problem [64]. Interpreting the weight functions as virtual displacements and considering a symmetric stress tensor field leads to the principle of virtual work [64]:

Space discretization
The body Ω is discretized by finite elements in this paper. Using interpolation functions N (e) j (x) for every element e, the nodal element displacement vector u (e) and its second derivative with respect to time as well as the virtual nodal element displacement vector are introduced [65]: Equivalently, the strain field ε and the virtual strain field δε are expressed by the nodal element displacements and the strain-displacement relationship B (e) j [65]: Inserting relations (8)- (12) into equation (7) and summing up all element volumes leads to the discrete version of the principle of virtual work. The assembling procedure over all leads from the local element node coordinates u (e) to the global nodal coordinate u. The spatially discretized version of equation (7) is written in the global nodal coordinate as: The assembly procedure leads from the element mass matrix M (e) to the global mass matrix M ∈ IR n s ×n s , from the internal restoring element force f (e) int to the global internal force vector f int ∈ IR n s , from the element volume force f (e) v to the global volume force vector f v ∈ IR n s and from the external element force vector f (e) ext to the external global force vector f ext ∈ IR n s . In this paper, we present the new strategy applying a linear elastic material law. Thus, the internal restoring element force vector is described by the relationship: whereC is the elasticity matrix following Voigt's notation and K (e) is the element stiffness matrix. Assembly of all element stiffness matrices leads then to the global stiffness matrix:

The undamped set of equations of motion and modal analysis
Considering Eqs. (13) and (15) leads to the matrix set of equations of motion:

123
where we combined the excitation on the right hand side by adding up the external and the volume force vector f(t) = f v + f ext . In order to obtain the linear modes of vibration, we are looking at the homogeneous problem only. Thus, setting the force at the right hand side of Eq. (16) to the zero vector and defining the solution ansatz function u(t) = e iωt , leads to the eigenvalue problem: As the matrix M −1 K is positive definite for the problem in this paper, one obtains n s real and positive eigenvalues (squared eigenfrequencies) ω 2 i (i = 1 . . . n s ). Magnitudewise arrangement of the eigenvalues from low to high values leads to the corresponding sorted set of the modes of vibration

Structural damping
In order to model dissipative effects, linear viscous damping is applied in this paper. A structural damping matrix C is included, using a linear combination of stiffness and mass matrix: where a 0 and a 1 are the Rayleigh damping coefficients to be evaluated as elaborately discussed in [1]. The linear viscously damped matrix set of equations of motion is then written as [1]:

Time discretization using the Newmark scheme
Firstly, we assume that Eq. (19), i.e. dynamic equilibrium, is fulfilled at fixed time instants t i , i = 1 . . . n t . Every time interval Δt = t i − t i−1 is chosen to be constant in this paper. The set of equations of motion at fixed time instants t i reads: where The conventional strategy is an incremental procedure in time, evaluating the solution u i based on the already known solution of all preceding time steps in a step-by-step manner, as shown in Fig. 1.
In between the fixed time instants t i , polynomial approximations are chosen. Within each time interval [t i−1 , t i ], a second order approximation of the solution function is chosen using the implicit Newmark ansatz functions, cf. [66,67]: Fig. 1 Conventional step-by-step integration procedure The parameters β and γ are chosen as 1 4 and 1 2 in this paper. These settings result in constant acceleration within each interval [t i−1 , t i ] and unconditional stability with respect to the size of the chosen time step length Δt [1,68,69].
The displacement increment and its derivatives at time instant t i are introduced as: and the incremental force at time instant t i is evaluated as: The set of discretized equations of motion (20) at time instant t i is rewritten in incremental form as: Equivalently, the Newmark formulations (21) and (22) are rewritten in incremental form as: The three unknowns Δu i , Δu i and Δü i are solved using the three equations (27), (28) and (29). Fig. 2 Illustration of a cantilever beam in the space-time framework Ω × [0, T ]; space-time displacement vector field U(x, t); spatial and temporal functions u s (x) and u t (t) highlighted in red and blue color, respectively t, as indicated by the continuous space-time response surface U(x, t). The goal is to provide this space-time response surface introducing a methodology that does not proceed in a step-by-step manner but that solves the whole space-time problem at once. In doing so, a non-incremental Newmark formulation in elasto-dynamics is proposed.

Continuous space-time decomposition
Let us consider a continuous system of which the evolution is described by a space-time displacement vector field (x, t) → U(x, t). Let x → u s (x) and t → u t (t) be spatial and temporal functions, respectively. Firstly, we define the tensor product u s ⊗ u t such that where the right hand side product is being performed component-wise indicated by the symbol < • >.
Observing Fig. 2, one can easily see that, in a first attempt, the space-time displacement vector field U(x, t) can be approximated by the tensor product of two functions u s (x) and u t (t), indicated by the red and blue lines in Fig. 2. As the displacement vector field U(x, t) can be arbitrary complex, an approximation of only one tensorial product could result in a bad representation of the solution. Thus, we claim that one can define the exact response time history by an infinite sum of displacement enrichments: where each enrichment term is a tensor product The response patterns in elastodynamics strongly depend on the type of the structure as well as the excitation time history. Following this line of thought, it is reasonable to hypothesize that a finite or even small number of enrichments will lead to a satisfying approximation of the displacement vector field U(x, t). This inevitably leads to the idea of model order reduction.
The separation ansatz allows one to easily define the derivatives of the displacement fields with respect to time, i.e. velocity and acceleration fields: Equivalently to the displacement velocity and acceleration fields, we define the decomposition of the force vector field and t → f t (t) be the spatial and temporal force functions. Applying the definition of the tensor product from Eq. (30), the force vector field is also defined as a sum of enrichments: The case of the stiffness is less straightforward because it is a functional linear operator mapping the space of displacement fields onto the space of the force fields. Considering a purely static problem, it provides the resulting force at a position x from a displacement acting at position x . It may be seen as the kernel (x, x , t) → K(x, x , t) (or Green's function) of a functional operator u → f such that: where Ω is the body defined in equation (1), on which the contracted product of the two rank tensor field K and the vector field u occurs. Let (x, x ) → K s (x, x ) and t → k t (t) be spatial and temporal functions, respectively. In general, the tensor product K s ⊗ k t is such that: where the product in the right hand side being performed component-wise. The evolution of the stiffness is then also performed in terms of an infinite sum of enrichments each of which is decomposed by In this context, we note that the decomposition of the mass is equivalent to the decomposition of the stiffness.

Discrete space-time decomposition
In its discrete version, there are n s degrees of freedom known at n t instants of time. The discrete representative of the continuous field U is a n s × n t matrix U: where u i = u(t i ) is the vector containing the displacement data at the consecutive time instants t i of the discrete model. The discrete representatives of the spatial and temporal functions u s and u t are the column vectors: The discrete counterpart of the tensor product U ( j) (x, t) = (u s ⊗ u t ) ( j) of functions is the tensor product of vectors: The object U ( j) is one n s ×n t matrix of which the element u kl is the product of (u s ) k and (u t ) l . The discrete representation of equations (31) and (32) is then the summation of M tensor products: Each summation term can be seen as an enrichment by a discrete space-time solution field. An upper boundary for the number M, necessary for an exact representation of U by the sum of dyadic vector product enrichments (43), is the number of the spatial degrees of freedom times the number of the time steps. Thus, the condition holds: However, it is the objective of this paper to find a significantly smaller number m M that enables for an approximation of the space-time displacement matrix U, while preserving the required level of accuracy for problems in structural dynamics. In Sect. 4, we will show that indeed a surprisingly small number of enrichments is sufficient for a reliable response identification of the dynamical system.
Equivalently, the discrete representations of the continuous velocity and acceleration fields, (33) and (34), are written as: The discrete version of the continuous operator K(x, x , t) is the usual stiffness matrix K, as presented in Eq. (15). Resuming the framework developed for the displacement and taking into account that the stiffness operator depends on two spatial variables and one time variable, we claim that the discrete representatives of the functions (K s ) ( j) and Then, the discrete counterpart of the tensor product Thus, we claim that the stiffness in space-time is a tensor of third order evaluated as a sum of decomposed enrichments: Equivalently, the space-time mass is defined as a tensor of third order: As discussed in Appendix A, the following rule can be applied: This is a general rule and, in this paper, it is applicable to all third order tensors and the respective tensor product of vectors. It will be applied in the upcoming sections of the paper.
In the special case of linear elastodynamics, stiffness and mass are both invariant with respect to time. The third order tensors are defined by one tensor product only: Here, 1 t ∈ R n t demonstrates a vector with every entry equal to one. Taking Rayleigh damping into account, as defined in Eq. (18), also the damping tensor is defined as:

Space-time set of equations of motion
Taking into account the stiffness, mass and damping terms, presented in Eqs. (52)- (54), as well as the discrete representation of the decomposed displacement, velocity and acceleration fields, presented in Equations (43)-(46), the set of equations of motion in tensorial space-time formalism is written as: where . . , f(t n t )] denotes the discretized excitation history matrix. As the mass, the damping and the stiffness terms are invariant with respect to time, they can be moved into the sum of the respective time-acceleration-, time-velocity-and time-displacement decomposition terms. Applying the rule for the product of tensors, defined in equation (51), leads to a simplified version of Eq. (55) that is only represented by tensors of rank two:

Newmark space-time formulations
Firstly, the displacement increment Δu is transformed into its space-time equivalent. In matrix form, it is written as: The transformation from the incremental formulation into the space-time equivalent is then defined as: where I s denotes the space identity matrix n s × n s , A t is a time coefficient matrix with the dimension n t × n t and a 0 is an n t -dimensional coefficient vector considering the initial displacement u 0 : For clarification, the Newmark interpolations (28) and (29) are represented by five terms, t 1 -t 5 , each of which are transformed into their space-time equivalent, adapting the above proposed approach. In doing so, the Newmark formulations, taken from Eqs. (21) and (22), are rewritten as: Equivalent to the transformation (58), the space-time formulations of the defined terms in Eqs. (60) and (61) are defined as: The matrices A 1 , A 2 , A 3 and A 4 are space-time Newmark coefficient matrices, containing β and γ , as well as the time step size Δt. They are all of dimension n t × n t and listed in Appendix B. Furthermore, the space-time equivalents of the initial conditions u 0 ,u 0 ,ü 0 are defined by the vectors a 5a 8 . Regarding all four vectors only the first entry is not equal to zero. They are listed in Appendix C. We note that A 1 -A 4 are sparse matrices as most of their components are equal to zero. This property can be considered when implementing the proposed strategy in order to further increase the level of efficiency.
Inserting the space-time formulations in equations (62-66) into the Newmark equations leads to the space-time equivalent of the Newmark equations: We note again that, in these equations, U,U andÜ are identified as rank two tensors, represented as n s × n t matrices. The discretized space time equation We identify three space time Eqs. (56), (67) and (68) with three unknowns, i.e. space-time displacement, U, space-time velocityU and space-time accelerationÜ. However, we take advantage of the fact that, velocity and acceleration can be written as a function of displacement. Therefore, velocity and acceleration are expressed from the two space-time Newmark interpolations, (67) and (68), with respect to displacement and inserted into the space-time equation of motion (55). In doing so, the linear dynamic problem is expressed by one algebraic space-time equation only: In this equation, matrices W and Y are evaluated as: where The n s × n t matrix L is identified as a space-time excitation term and is evaluated as: with The vectors e, g, n and b are dependent on the coefficient matrices A 1 -A 4 and the coefficient vectors a 5 -a 8 : We note that the right hand side of equation (69) is still not decomposed using the space-time separation ansatz. Therefore, the proper orthogonal decomposition [15,70,71] is applied. In doing so, the space-time excitation term on the right hand side of the space-time set of equations of motion (69) is decomposed into a sum of dyadic products of uncorrelated space and time vectors. This is realized using the singular value decomposition [72]: Here, l s are orthonormal left singular vectors, taken from the matrix L s with dimension n s × n s and l (q) t are the right singular vectors taken from the matrix L t with dimension n t × n t . The matrix S is an n s × n t rectangular diagonal matrix with the singular values s (q) in descending order. In this paper, the number of singular basis vectors n tr is chosen so that 99.9% of the total energy, evaluated by the sum of all singular values n t q=1 s (q) , is captured [10]. In structural dynamics it is, generally, the case that an already surprisingly small number of basis vectors is sufficient to obtain such a high level of accuracy. [11,12,61].

The space-time solution strategy
In this paper, the fixed point algorithm is applied to solve the problem in the space-time framework. The solution at the m th enrichment, U (m) , is divided into already evaluated m − 1 known enrichments and the current m th enrichment solution to be calculated: (83) Figure 3 visualizes the iterative enrichment procedure as a step-wise adaption of the solution in the whole time domain. The solution ansatz (83) as well as the decomposition of the space-time excitation term (82) are inserted into the space-time equation (69): Fig. 3 Iterative space time procedure Equation (84) Analogously, multiplication with u (m) t T leads to a projection into space : fixed. This procedure is indicated by the first horizontal arrow in Fig. 4. These two steps are alternately repeated until convergence is achieved, indicated by the intersection point (u * s , u * t ) of the red and blue curves in Fig. 4.
In Eqs. (85) and (86) the scalar quantities a, b, c, e, f and g are evaluated by: The excitation vectors, projected into space and time for the m th enrichment, are written as: The vectors v s1 -v s3 in equation (85) and v t1 -v t3 in Eq. (86) are written as: Convergence is achieved if the enrichment vectors u (m) t and u (m) s fulfill the Frobenius norm of the space-time equation of motion with sufficient accuracy. This is defined through a scalar parameter ε f p : where ε enr is the global convergence limit to be defined. Additionally, the proposed space-time algorithm is summarized in Algorithm 1.

Numerical examples
For the purpose of demonstration, we introduce a twodimensional frame structure subjected to four different excitation functions. Regarding all the numerical examples, the proposed Newmark space-time approach is presented and compared to the conventional Newmark integration scheme. In doing so, the focus is not only on the development of the spatial and temporal enrichment vectors during convergence but also on the limits of the proposed strategy. Figure 5a presents an illustration of the frame structure subjected to the different chosen excitations. In particular, three different horizontal ground excitationsẍ g and one impact force function f (t) are considered.
Introducing ground excitation, the external force on the right hand side of the set of equations of motion (19) is substituted by the term of inertia induced by the seismic acceleration of the ground [1] :  In this paper, f H is an influence vector, defining all horizontal degrees of freedom that are affected by the ground acceleration time historyẍ g (t). As first ground excitation example, a sinusoidal function x (h) g (t) = a 0 sin νt is chosen with an amplitude of a 0 = 0.1 m/s 2 and an excitation frequency of ν = 1.61 rad/s. As second ground excitation example, also a sinusoidal function is chosen with the same amplitude of a 0 = 0.1 m/s 2 but with a significantly higher excitation frequency compared to the first example. In this context, the excitation frequency for the second example is chosen as 5.0 rad/s. As third ground excitation example, a transient time history is chosen that is taken from an earthquake recorded in Kobe, Japan [73], as shown by the plot in the middle of Fig. 5b. The time duration of the record is 22 seconds and the peak ground acceleration (PGA) is 6.8 m/s 2 . Finally, as fourth example, an impact force function f (t), as shown in the bottom plot of Fig. 5b, is presented. In particular, a sinusoidal function is considered during the first 0.1 seconds followed by 0.9 seconds of a force function equal to zero. All excitations are also summarized in Table 1.
The frame structure is discretized by two-node truss-beam elements, applying the Euler-Bernoulli beam theory [74]. Both columns and beams of each story are discretized by five elements, which leads to a total number of 128 degrees of freedom. The time step Δt is chosen as 10 −2 seconds and the parameters of the numerical example are chosen as: stiffness E = 210 GPa, density ρ = 7850 kg/m 3 . All cross sections of the columns and beams are modeled by quadratic hollow rectangular shapes with the dimension 0.3 m × 0.1 m and a thickness of t = 10 mm.
The output is compared to the conventional step-by-step Newmark solution. Regarding the first three examples, the output quantity considered for comparison is the horizontal displacement history in the left corner of the frame structure, as indicated by the magenta cross and the magenta arrow in Fig. 5a. Regarding the fourth example, the output quantity of interest considered for comparison is the degree of freedom v(t) indicated by the red arrow in Fig. 5a. The convergence limit ε (FN) enr is calculated as the average value over all time steps. Thus, if the quantity ε (FN) enr /n t falls below the value of 1.0e − 4, convergence is achieved.

Example 1: harmonic excitation-resonance
Firstly, the structure is excited by a sinusoidal ground excitation at a frequency of ν = 1.61 rad/s. This frequency is equal to the first eigenfrequency ω (1) that corresponds to the first eigenmode, depicted in Fig. 6. For this numerical example, damping is neglected. The results are depicted in Fig. 7. The Newmark space-time algorithm is presented by the dashed blue line in the bottom plot of Fig. 7, while the Newmark benchmark solution is presented by the black line. The space-time algorithm converges already after the first enrichment m = 1. This is caused by the fact that, in the special case of resonance with respect to the first eigenfrequency, the structure responds in the first mode shape only. In this special situation, the first normalized spatial enrichment vector coincides with the first normalized mode of vibration (1) . The corresponding temporal enrichment vector is depicted in the right top subplot in Fig. 7, clearly following the response characteristics of the resonance effect excited by a sinusoidal excitation. The full solution is then exactly evaluated by the dyadic product of the first spatial and temporal enrichment vectors u (1) s ⊗ u (1) t . 123 s ≡ Φ(1) ) ) ) ) ) ) ) ) Results applying a sinusoidal ground excitationẍ g (t) = a 0 sin νt with ν = 1.61 rad/s, which is equal to the first eigenfrequency of the system; spatial enrichment vectors in red, temporal enrichment vectors highlighted in blue; bottom plot: response function of the output degree of freedom h(t) applying the Newmark space-time approach (blue dashed line) and the conventional Newmark algorithm (black solid line); the Newmark space-time approach already converges after the first enrichment. Thus, the solution is decomposed of the dyadic product u (1) s ⊗ u (1) t . In this special situation it is:

Example 2: harmonic excitation-no resonance
Secondly, the structure is excited by a sinusoidal ground excitation, but, in this case, with a frequency that is higher than the first eigenfrequency. It is chosen as ν = 5.0 rad/s. Also, regarding this example, damping is neglected. The results are depicted in Fig. 8. Convergence is already achieved after two enrichments. The spatial enrichments modes are presented in red, while the corresponding temporal enrichments vectors are presented in blue. Both spatial modes, u (1) s and u (2) s , are linearly independent and they show a high similarity with the first two linear modes of vibration, shown in Fig. 6 . This phenomenon is not surprising as, naturally, the response is mainly composed of the first two eigenmodes if the excitation frequency lies between the first and the second eigenfrequency. However, looking at the corresponding temporal enrichment functions shows that the proposed strategy converges differently than the conventional modal truncation method would, as both the first and the second spatial enrichment function contain more than only one frequency. The horizontal response time history h(t) is shown in the bottom plot. The spacetime solution is depicted by the dashed blue line, while the Newmark benchmark solution is represented by the black solid line. Also regarding this example, excellent agreement is observed when comparing the new strategy with the benchmark solution after the two enrichments. The whole convergence behavior is also presented in the bottom plot of Fig. 8, where the time step average of the Frobenius norm ε (FN) enr /n t versus the number of considered enrichments is depicted.

Example 3: transient excitation
Thirdly, the structure is excited by an earthquake, recorded in Kobe, Japan [73], as shown in the middle plot of Fig. 5b. In this numerical example, Rayleigh damping is introduced using a damping ratio of 4% for the second and sixth mode of vibration [1]. The proposed Newmark space-time strategy already converges after the first two enrichments. However, we also demonstrate the results of the third and fourth enrichment for the purpose of illustration. The results are depicted in Fig. 9.
The horizontal response time history h(t) after two enrichments is shown in the bottom plot of Fig. 9. The blue dashed line shows the Newmark space-time solution, while the black solid line shows the conventional Newmark solution. Also regarding this example, all spatial modes are linearly independent and a high similarity with the first four linear modes of vibration, shown in Fig. 6, is observed.

Example 4: impact force
Finally, the structure is excited by a nodal force time history f (t), as shown in Fig. 5a. This excitation can be interpreted as an impact load with a duration of 0.1 seconds. While the three other examples induce global response patterns, exciting all horizontal degrees of freedom equally introducing a mass dependent term of inertia, the last example is characterized by a rather localized effect. The results regarding the fourth example are presented in Fig. 10.
Also here, the spatial enrichment modes are depicted in red color, while the temporal enrichment modes are depicted in blue color. All spatial modes are linearly independent. While, regarding the other three examples, the first enrich- enr /n t is considerably less distinct compared to the convergence trend of the other examples. However, the converged solution after the 20 th enrichment is also here in excellent agreement with the step-by-step Newmark benchmark response, which indicates the robustness of the proposed algorithm. The classical method of modal truncation would, regarding the present example, not lead to a satisfying approximation of the response function as the local deformation patterns could kinematically not be represented by a spectral basis, evaluated by considering a small number of low-frequency modes from equation (17). Thus, an investigation of the mode shapes would a priori be necessary in order to extract precisely those modes that are responsible for these localized deformation pattern around the impact point. However, using the Newmark space-time approach, the localized deformation patterns are automatically included within the first spatial enrichment vectors as shown in Fig. 10 and no a priori investigation of the dynamic system is required. In this context, this example also shows the clear advantage of the Newmark space-time approach to the classical method of modal truncation.

Computational efficiency
The computational cost of the new strategy depends on the number of enrichments necessary for convergence and, therefore, on the convergence interval ε (FN) enr chosen. As a scalar error quantity U , we introduce the trace of the covariance matrix of the displacement matrix after the m th enrichment, U (m) , which is equivalent to the Frobenius norm of the global error matrix: where U is the response time history obtained using the conventional step-by-step Newmark method. The computational efficiency is investigated, using the numerical example 3, presented in Sect. 4.3, by varying the number of enrichments considered and comparing the corresponding calculation time. The result is depicted in Fig. 11.  The blue line presents the error quantity U versus the number of enrichments considered. A clear strong downward trend is observed within the first few enrichments. We note that the convergence criterion chosen for all examples (ε (FN) enr ≤ 1.0e − 4) is here already satisfied after the second enrichment. However, we increase the number of enrichments for the sake of illustration and do not stop the algorithm after the second enrichment. The magenta line denotes the calculation time versus the number of enrichments. The black dashed line presents the calculation time that is needed when using the conventional Newmark algorithm. For comparison, all other calculation times are normalized with respect to that value. Using four or less enrichments, the proposed method is faster than the classical Newmark method. Increasing the number of enrichments, the calculation time grows linearly with the number of enrichments. Regarding the example discussed in Sect. 4.3, a speed-up by a factor of 2.08 is observed. Additionally, the results regarding computational efficiency are presented in Table 2.

Conclusion
In this paper, a Newmark space-time formulation for linear transient structural analysis is proposed.
The proper generalized decomposition is applied to decompose the incremental formulations into their spacetime equivalent. Instead of proceeding incrementally in time, the system is solved in the whole time domain in terms of enrichment iterations. Within every enrichment iteration a system of two coupled space-time equations is solved using the fixed point algorithm. During convergence, the spatial and temporal enrichment functions are investigated.
The new method is demonstrated, using four different numerical examples. A frame structure is excited by a sinusoidal ground excitation at the first eigenfrequency, a second sinusoidal ground excitation that is larger than the first eigenfrequency, a transient ground excitation and an impact force function that is applied locally.
Exciting the structure by a sinusoidal ground acceleration with the excitation frequency equal to the first eigenfrequency leads to exactly one enrichment step until convergence is achieved. In the special case of resonance, the normalized spatial enrichment vector is equivalent to the normalized first linear mode of vibration and the shape of the temporal enrichment vector follows the response time history of the reference solution.
The structure is then excited by a sinusoidal ground excitation with an excitation frequency that lies between the first and the second eigenfrequency as well as by a transient ground excitation. In both cases the strategy leads to convergence already after two enrichments and the spatial enrichment functions show a high similarity with the corresponding low-frequency eigenmodes.
Finally, the structure is excited by a local impact-type force function. In this case, the first spatial enrichment modes do not reveal a high correlation with the first few modes of vibration. Instead deformation patterns, containing distinct localized shapes, are observed regarding the spatial modes describing both the local and global structural response behavior caused by the impact force function. This indicates that the proposed strategy is superior to classical strategies like the method of modal truncation, as it reveals an automatic adaption of the response patterns in form of an appropriate choice of spatial and temporal enrichment modes. In this regard, the method also outperforms the proper orthogonal decomposition method, as the performance of the latter depends on an a priori system identification for the generation of an appropriate reduced order basis [10,11].
Generally, the computational efficiency of the proposed strategy depends on the number of enrichments necessary to obtain a required level of accuracy. The numerical demonstrations in this paper show that a computational speed-up factor up to a value of four compared to the conventional Newmark method can be achieved with a high level of accuracy regarding problems in structural dynamics. While the conventional Newmark integration scheme proceeds in a step-by-step manner in time and is, therefore, only realiz-able by a considerably high number of calculation steps, it is shown in this paper that the Newmark space-time strategy needs only a surprisingly small number of enrichment iterations until convergence is achieved.
The proposed strategy will serve as a guidance for future space-time approaches in nonlinear structural dynamics. The next step should be an extension of the proposed approach to nonlinear structures. It is essential to extend the framework to geometrically nonlinear problems as, e.g., the consideration of large deformation strain-displacement relationships or contact constraints, as well as problems considering materially nonlinear, history-dependent material behavior. Finally, it would be of high interest to establish the space-time formulations of other integration schemes and answer crucial research questions regarding applicability and efficiency. In particular, the elaboration of a space-time equivalent of the popular central difference integration scheme [2] that reveals only conditional stability in its classical formulation could provide a new perception of explicit integrators in structural dynamics.
Acknowledgements This research was financially supported by Deutsche Forschungsgemeinschaft (DFG) and the Agence Nationale de la Recherche (ANR) with project number MA2233/10-1.
Funding Open Access funding enabled and organized by Projekt DEAL.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/.

A Multiplication rule of decomposed rank-three tensors tensors with decomposed rank-two tensors
In terms of components for a fix enrichment number j we get: In the following, we need to marry the tensor products for the displacement and the stiffness with the contracted product between stiffness and displacement. In index notation we get: or in abstract notations: Combining the relations (42) and (48) allows to define the rule:

B Newmark coefficients in space-time
The four Newmark space-time coefficients are n t × n t matrices. They are defined as:

C Space-time initial conditions
The space-time equivalents of the initial conditions are defined by four vectors a 5 -a 8 that are all of the dimension n t . Only the first component of each of the four vectors is not equal to zero: