import theme
theme.load_style()
This lecture by Tim Fuller is licensed under the Creative Commons Attribution 4.0 International License. All code examples are also licensed under the MIT license.
$n$ is the normal to the boundary, and it allows us to track the direction of the traction force applied to the domain boundary.
|
Let us now re-consider the wire problem, except this time the right side of wire is no longer fixed but is fixed to a spring with spring constant $k$, as shown below
The governing equation for this case is
$$ T\frac{d^2u}{d{x}^2}+q=0 $$$$ u(0)=0, \qquad \left( T\frac{du}{dx}\right)\Bigg|_{x=L}=ku_L $$For a mesh of two linear elements, the global system of equations is:
$$ \frac{2T}{L}\begin{pmatrix} 1 & -1 & 0 \newline -1 & 2 & -1 \newline 0 & -1 & 1 \end{pmatrix} \begin{pmatrix} u_1 \newline u_2 \newline u_3 \end{pmatrix} = \frac{qL}{4}\begin{pmatrix}1 \newline 2 \newline 1 \end{pmatrix} + \begin{pmatrix} Q_1^1 \newline Q_2^1+Q_1^2 \newline Q_2^2 \end{pmatrix} $$Applying boundary conditions
$$ \frac{2T}{L}\begin{pmatrix} 1 & -1 & 0 \newline -1 & 2 & -1 \newline 0 & -1 & 1 \end{pmatrix} \begin{pmatrix} 0 \newline u_2 \newline u_3 \end{pmatrix} = \frac{qL}{4}\begin{pmatrix} 1 \newline 2 \newline 1 \end{pmatrix} + \begin{pmatrix} Q_1^1 \newline 0 \newline ku_3 \end{pmatrix} $$We re-write so that all unknowns are on the LHS: $$ \frac{2T}{L}\begin{pmatrix} 1 & -1 & 0 \newline -1 & 2 & -1 \newline 0 & -1 & 1 - \frac{kL}{2T} \end{pmatrix} \begin{pmatrix} 0 \newline u_2 \newline u_3 \end{pmatrix} = \frac{qL}{4}\begin{pmatrix} 1 \newline 2 \newline 1 \end{pmatrix}
\frac{2T}{L}\begin{pmatrix} 2 & -1 \newline -1 & 1-\frac{kL}{2T} \end{pmatrix} \begin{pmatrix} u_2 \newline u_3 \end{pmatrix} = \frac{qL}{4}\begin{pmatrix} 2 \newline 1 \end{pmatrix} $$
which can then be solved for $u_2$ and $u_3$. We then find the reaction $R_1$
$$ R_1=-\frac{2T}{L}u_2-\frac{qL}{2} $$Let us revisit the problem of the tensioned wire with elastic support, shown below
The boundary condition at $x=L$ becomes $F=ku_L$, where, by balancing torque
$$ F = \frac{1}{L} \int_0^L x q(x) dx $$In this case, neither $u$ nor $u'$ is known at $x=L$. If we substitute $q(x) = -T u''(x)$, integration by parts gives
$$ F = -\frac{T}{L}\left(L u'(L) - u(L)\right) = T\frac{u(L)}{L} - Tu'(L) $$The boundary condition is of the form
$$ \alpha u + \beta u' = \gamma $$where $\alpha$, $\beta$, and $\gamma$ are known. The text (Ch 3.7) calls this type of boundary condition "generalized boundary conditions". This type of boundary condition arises in applications, such as convection, or elastic support of elastic bars.
Note: prescribed displacements corresponds to $\alpha=1$, $\beta=0$ and prescribed force corresponds to $\alpha=0$, $\beta=1$.
For the mixed case in which $\beta$ is nonzero at both ends of the wire, the boundary conditions may be written
$$ u'(0) = \frac{\gamma_0}{\beta_0} - \frac{\alpha_0}{\beta_0} u(0) $$$$ u'(L) = \frac{\gamma_L}{\beta_L} - \frac{\alpha_L}{\beta_L} u(L) $$Substituting $u(x) = u_i \phi_i(x)$ and $w(x) = w_i \phi_i(x)$ eventually leads to the following system of equations for the unknown nodal displacements
$$ K^*_{ij} u_j = f^*_i $$where $K^*_{ij}$ and $f^*_i$ are the same as $K_{ij}$ and $f_i$ except
$$ \boxed{ \begin{eqnarray} K^*_{11} &= K_{11} - \frac{\alpha_0}{\beta_0} \\ K^*_{NN} &= K_{NN} + \frac{\alpha_L}{\beta_L} \\ f^*_{1} &= f_{1} - \frac{\gamma_0}{\beta_0} \\ f^*_{N} &= f_{N} + \frac{\gamma_L}{\beta_L} \end{eqnarray} } $$This holds for nonzero $\beta_0$ and $\beta_L$. In the case where $\beta_0=\beta_L=0$, replace $\beta$ witha tiny number $(10^{-9})$ and the result will be the penalty method!
This if boundary conditions are of general form, no special case coding is needed! Actually, supporting more complicated boundary conditions, in this case, actually simplifies the structure of the code.
Sidenote
A code should include a solvability check to ensure that
1D equations for heat flow, diffusion and elasticity are all two‐point boundary problems of the form
$$ \frac{d}{dx}\left(a \frac{du}{dx}\right) + bu + c = 0, \quad \text{on }\Omega $$We can introduce a Generalized Boundary Condition statement
$$ \left(\alpha n \frac{du}{dx}-\overline{\Phi}\right) + \beta\left(u-\overline{u}\right)=0, \quad \text{on } \Gamma_{\Phi} $$With the penalty method, we can apply this condition to the entire domain boundary, selecting parameters to reduce this expression to Dirichlet, Neumann, or Robin boundary conditions.
For very large $\beta$, $u=\overline{u}$. |
For $\beta=0$, $\frac{du}{dx} = \frac{\Phi}{\alpha}$. |
Options for modeling a self‐weighted bar, supported by a spring
The generalized boundary condition approach can be extended to less trivial problems with elastic support and provides a more efficient method to defining the boundary interaction.
|
For the simple uniform bar with constant density shown, we can compute the total spring force, and use the spring elongation to define the prescribed displacement at $x=0$.
$$ u(0) = \delta_s = \frac{F_s}{k_s}=\frac{P_L+\rho gAL}{k_s} $$In general we may have a density, area, or body force that vary with position so the force at the spring would not be known until the system is solved, and we can only prescribe a force‐displacement relationship at $x=0$.
$$F_s = k_su(0)$$We can obtain this result more formally using the general form of Robin Boundary Condition,
$$ \left(E(0)n(0)\frac{du}{dx}\Bigg|_{x=0}-\overline{t}\right) + \frac{k_s}{A(0)}\left(u(0) - \overline{u}\right)=0, \quad \text{at } x=0 $$The reference traction and reference displacement at $x=0$ are both zero, and the outward normal to the boundary is in the $-x$ direction. Substituting these values gives
$$ \left(EA(-1)\frac{du}{dx}\right)_{x=0}+k_s\left(u(x)\right)_{x=0} $$