# How to impose boundary conditions when solving a nonlinear dynamical system given by the FEM solver

I am solving a nonlinear dynamical system given by a nonlinear elastic problem which takes the following form:

$$\boldsymbol{M} \ddot{u} + \boldsymbol{K}_{\textrm{NL}}u = 0 ,$$

here $$u \in \mathbb{R}^N$$ is the displacement and $$\ddot{u}$$ is the second time derivative of the displacement or the acceleration, $$\boldsymbol{M} \in \mathbb{R}^{N\times N}$$ and $$\boldsymbol{K}_{\textrm{NL}} \in \mathbb{R}^{N\times N}$$ are the mass and stiffness matrices, respectively. In this problem, the stiffness matrix is a function of the current displacement solution at time step $$t_k$$ during the simulation which means that $$\boldsymbol{K}_{\textrm{NL}}$$ varies with time. The problem above is solved together with boundary conditions (for imposing excitation in the system) for the displacement at two points at the boundary $$\Gamma$$ which read as follows:

$$u|_{p_1,p_2} = u_L(t) = V_{\textrm{amp}} \sin(2\pi f_{\textrm{ex}} t) \sin^2(\frac{\pi f_{\textrm{ex}} t}{n}), \quad \forall t,$$ where $$p_1, p_2 \in \Gamma$$ are the excitation points, $$V_{\textrm{amp}}$$ is the amplitude of the excitation, $$f_{\textrm{ex}}$$ is the frequency of the excitation and $$n=5$$

I am assuming the following scenario: given the mass matrix $$\boldsymbol{M}$$ and the stiffness matrices at each time step during the simulation from $$t_i$$ till $$t_f$$, I would like to solve the system for the solution $$u_k$$ at each time step $$t_k \in [t_i,t_f]$$.

The process is done in the literature usually using the Newmark method which is a time-stepping scheme for general 2nd order systems like the one above. The non-linearity of the problem requires an additional treatment using the Newton method for example. I am having trouble ensuring that the solution of the system holds the boundary conditions. Essentially I need to build the Newton-Newmark method in a way that does not vary the value of $$u$$ at the two points $$p_1$$ and $$p_2$$, of course, the value of $$\dot{u}$$ and $$\ddot{u}$$ is known as well. How should one try to do that? This question is similar to a previous one here

Boundary condtions on nonlinear FEM time integration

The author there suggests to remove the rows and columns, which correspond to the degrees of freedom involved in the boundary conditions, from the matrices of the system. However, this option would make it impossible to march a meaningful value of $$u$$ (in my case) since it is basically the boundary condition which is propagating in the domain and causing the evolution of the displacement.

• So you are imposing these nonhomogeneous Dirichlet conditions at exactly 2 points in the boundary? How is that implemented even without the nonlinearity. What are the BCs everywhere else? Feb 13 at 16:16
• I am puzzled by the equation for your dynamic system. $K_{NL}u$ should represent the internal forces in the FEM model. For this to be so, $K_{NL}$ would have to be a secant stiffness matrix which is a concept rarely seen in nonlinear solid mechanics problems. Feb 14 at 1:00
• @whpowell96 I should have mentioned that there is a third point at which the boundary condition is a homogeneous Dirichlet BC, everywhere else at the boundary the solution is assumed to be varying freely. I don't understand the first part of your question but I would like to say that I am using a commercial software so I don't have full access on how the problem is solved or how the BCs are imposed. Feb 14 at 13:21
• @BillGreene Here K is the so-called tangent stiffness matrix, and it is true that the nonlinear part represents the internal forces in the model, for more info on this you can see the book of Chopra 'DYNAMICS OF STRUCTURES', section 5.7 for example which is related to my question. Feb 14 at 13:26
• If $K_{NL}$ is the tangent stiffness matrix and also a function a function of $u$, then $K_{NL}u$ CANNOT be the internal forces in the structure, and your first equation is incorrect. I've seen others make the same mistake which is why I asked. Feb 14 at 20:06

The general approach in a Newton method is to solve for increments $$\delta u$$. Then, if your current solution already has the right values at the boundary, $$\delta u$$ needs to have zero boundary values. You can see this kind of scheme implemented in step-15 of the deal.II library see https://dealii.org/developer/doxygen/deal.II/step_15.html . In your case, the nonlinear to be solved is the system you get at each time step, after time discretizations.