## Theoretical Modeling
### End-Point Algorithm  
The discretized dynamic equilibrium equation at the end of \(n+1\) time step is given by:  
\[\mathbf{q}_{i,n+1}+\mathbf{q}_{mas,n+1}-\mathbf{q}_{e,n+1}=\mathbf{0} \tag{1}\]  
where \(\mathbf{q}_{i,n+1}\), \(\mathbf{q}_{mas,n+1}\) and \(\mathbf{q}_{e,n+1}\) are, respectively, the nodal internal, inertial and external force vectors at time \(n+1\).  
According to [PERSON] [37], \(\mathbf{q}_{mas}\) can be expressed as:  
\[\mathbf{q}_{mas}=\mathbf{U}_{e}^{*}\mathbf{f}_{inter}=\mathbf{U}_{e}^{*}\left( \mathbf{M}\ddot{\mathbf{p}}+\frac{l}{12}\left[\begin{array}{c}\mathbf{0}_{3 \times 1}\\ \mathbf{3}\mathbf{S}(\mathbf{w}_{1})\mathbf{I}_{p}\mathbf{w}_{1}+\mathbf{S}( \mathbf{w}_{1})\mathbf{I}_{p}\mathbf{w}_{2}+\mathbf{S}(\mathbf{w}_{2})\mathbf{ I}_{p}\mathbf{w}_{1}+\mathbf{S}(\mathbf{w}_{2})\mathbf{I}_{p}\mathbf{w}_{2}\\ \mathbf{0}_{3\times 1}\\ \mathbf{S}(\mathbf{w}_{1})\mathbf{I}_{p}\mathbf{w}_{1}+\mathbf{S}(\mathbf{w}_{ 1})\mathbf{I}_{p}\mathbf{w}_{2}+\mathbf{S}(\mathbf{w}_{2})\mathbf{I}_{p} \mathbf{w}_{1}+3\mathbf{S}(\mathbf{w}_{2})\mathbf{I}_{p}\mathbf{w}_{2}\\ \end{array}\right]\right) \tag{2}\]  
where \(l\) is the length of element; \(\mathbf{M}\) is the consistent mass matrix; \(\mathbf{I}_{p}\) is the tensor of mass moment of inertia; \(\mathbf{w}_{i}\) is the vector of angular velocity at node \(i\) in the body-attached elemental frame \(\mathbf{U}_{e}\); and \(\mathbf{S}(\cdot)\) is the \(3\times 3\) skew-symmetric matrix. \(\dot{\mathbf{p}}\) is the vector of nodal accelerations, whose translational components are defined in global reference system \(\mathbf{U}_{g}\) and rotational ones are given in the body-attached frame \(\mathbf{U}_{e}\). Furthermore,  
\[\mathbf{U}_{e}^{*}=diag(\mathbf{I}_{3\times 3},\mathbf{U}_{e},\mathbf{I}_{3 \times 3},\mathbf{U}_{e}). \tag{3}\]  
is the diagonal matrix composed of \(\mathbf{U}_{e}\) and \(\mathbf{I}_{3\times 3}\), where \(\mathbf{I}_{3\times 3}\) is the \(3\times 3\) identity matrix. The superscript (\({}^{*}\)) is used only to distinguish the \(3\times 3\) elemental frame \(\mathbf{U}_{e}\) and the \(12\times 12\) diagonal matrix.  
Equation (2) implies that the initial elemental reference frame \(\mathbf{U}_{0}\) should be in the same direction as the global one. It is no longer applicable when the two frames point in different directions. Thus, to make it more generic and adaptive, the inertial force vector \(\mathbf{q}_{mas}\) in the present framework is modified as:  
\[\mathbf{q}_{mas}=\mathbf{U}_{e}^{*}\mathbf{f}_{iner}=\mathbf{U}_{e}^{*}( \mathbf{M}\mathbf{U}_{0}^{*}\mathbf{\ddot{p}}+\frac{l}{12}\mathbf{\hat{f}}) \tag{4}\]  
and \(\mathbf{U}_{e}^{*}\) is modified as:  
\[\mathbf{U}_{e}^{*}=diag(\mathbf{U}_{0},\mathbf{U}_{e},\mathbf{U}_{0},\mathbf{ U}_{e}) \tag{5}\]  
where  
\[\mathbf{U}_{0}^{*}=diag(\mathbf{U}_{0}^{T},\mathbf{I}_{3\times 3},\mathbf{U}_{0 }^{T},\mathbf{I}_{3\times 3}) \tag{6}\]  
and \(\mathbf{U}_{0}\) is the initial elemental reference frame expressed in the global reference system.  
To solve such a system of algebraic non-linear equations with the Newton-Raphson method, it is necessary to compute the variations of vectors \(\mathbf{q}_{i,n+1}\) and \(\mathbf{q}_{mas,n+1}\) to take into account the contributions of internal and inertial forces to the generalized tangent stiffness matrix (details on the derivation of \(\delta\mathbf{q}_{i,n+1}\), which leads to the static tangent stiffness matrix \(\mathbf{K}_{t,stat}\), are given in Refs. [38; 39]).  
The variation of \(\mathbf{q}_{mas,n+1}\) can be expressed as:  
\[\delta\mathbf{q}_{mas,n+1}=\mathbf{K}_{t,mas}\delta\mathbf{p}=(\mathbf{K}_{mas 1}+\mathbf{K