Nonlinear finite elements/Homework 6/Solutions/Problem 1/Part 9

Problem 1: Part 9

Construct a laminar coordinate system at the blue point. (Use equations 9.3.16 and 9.3.17 from the book chapter). Assume that the blue point is at the center of element 5.

The $x$ and $y$ coordinates of the master and slave nodes are

${\begin{bmatrix}x_{1}\\y_{1}\\x_{2}\\y_{2}\end{bmatrix}}={\begin{bmatrix}r\cos \theta _{1}\\r\sin \theta _{1}\\r\cos \theta _{2}\\r\sin \theta _{2}\end{bmatrix}}={\begin{bmatrix}0.55000\\0.95263\\0.95263\\0.55000\end{bmatrix}}$

${\begin{bmatrix}x_{1-}\\y_{1-}\\x_{2-}\\y_{2-}\end{bmatrix}}={\begin{bmatrix}r_{1}\cos \theta _{1}\\r_{1}\sin \theta _{1}\\r_{1}\cos \theta _{2}\\r_{1}\sin \theta _{2}\end{bmatrix}}={\begin{bmatrix}0.50000\\0.86603\\0.86603\\0.50000\end{bmatrix}}$

${\begin{bmatrix}x_{1+}\\y_{1+}\\x_{2+}\\y_{2+}\end{bmatrix}}={\begin{bmatrix}r_{2}\cos \theta _{1}\\r_{2}\sin \theta _{1}\\r_{2}\cos \theta _{2}\\r_{2}\sin \theta _{2}\end{bmatrix}}={\begin{bmatrix}0.60000\\1.0392\\1.0392\\0.60000\end{bmatrix}}$

Since there are two master nodes in an element, the parent element is a four-noded quad with shape functions

${\begin{aligned}N_{1-}(\xi ,\eta )&={\cfrac {1}{4}}(1-\xi )(1-\eta )&N_{2-}(\xi ,\eta )&={\cfrac {1}{4}}(1+\xi )(1-\eta )\\N_{1+}(\xi ,\eta )&={\cfrac {1}{4}}(1-\xi )(1+\eta )&N_{2+}(\xi ,\eta )&={\cfrac {1}{4}}(1+\xi )(1+\eta )~.\end{aligned}}$

The isoparametric map is

${\begin{aligned}x(\xi ,\eta )&=x_{1-}~N_{1^{-}}(\xi ,\eta )+x_{2-}~N_{2^{-}}(\xi ,\eta )+x_{1+}~N_{1^{+}}(\xi ,\eta )+x_{2+}~N_{2^{+}}(\xi ,\eta )\\y(\xi ,\eta )&=y_{1-}~N_{1^{-}}(\xi ,\eta )+y_{2-}~N_{2^{-}}(\xi ,\eta )+y_{1+}~N_{1^{+}}(\xi ,\eta )+y_{2+}~N_{2^{+}}(\xi ,\eta )~.\end{aligned}}$

Plugging in the numbers, we get

${\begin{aligned}x&=0.12500(1-\xi )(1-\eta )+0.21651(1+\xi )(1-\eta )\\&+0.15000(1-\xi )(1+\eta )+0.25981(1+\xi )(1+\eta )\\y&=0.21651(1-\xi )(1-\eta )+0.12500(1+\xi )(1-\eta )\\&+0.25981(1-\xi )(1+\eta )+0.15000(1+\xi )(1+\eta )~.\end{aligned}}$

Therefore, the derivatives with respect to $\xi$ are

${\begin{aligned}x_{,\xi }={\frac {\partial x}{\partial \xi }}&=0.20131+0.018301\eta \\y_{,\xi }={\frac {\partial y}{\partial \xi }}&=-0.20131-0.018301\eta ~.\end{aligned}}$

Since the blue point is at the center of the element, the values of $\xi$ and $\eta$ at that point are zero. Therefore, the values of the derivatives at that point are

${\begin{aligned}x_{,\xi }={\frac {\partial x}{\partial \xi }}&=0.20131\\y_{,\xi }={\frac {\partial y}{\partial \xi }}&=-0.20131~.\end{aligned}}$

Therefore,

$\mathbf {x} _{,\xi }=0.20131\mathbf {e} _{x}-0.20131\mathbf {e} _{y},~~\lVert \mathbf {x} _{,\xi }\rVert _{}={\sqrt {0.20131^{2}+0.20131^{2}}}=0.28470~.$

The local laminar basis vector ${\widehat {\mathbf {e} }}_{x}$ is given by

${\widehat {\mathbf {e} }}_{x}={\cfrac {\mathbf {x} _{,\xi }}{\lVert \mathbf {x} _{,\xi }\rVert _{}}}=0.70711\mathbf {e} _{x}-0.70711\mathbf {e} _{y}~.$

The laminar basis vector ${\widehat {\mathbf {e} }}_{y}$ is given by

${\widehat {\mathbf {e} }}_{y}=\mathbf {e} _{z}\times {\widehat {\mathbf {e} _{x}}}=0.70711\mathbf {e} _{x}+0.70711\mathbf {e} _{y}~.$

A plot of the vectors is shown in Figure 13.

Figure 13. Laminar base vectors.

The Maple code for this calculation is shown below.

> #
> # Shape functions
> #
> N1m := 1/4*(1-xi)*(1-eta):
> N2m := 1/4*(1+xi)*(1-eta):
> N1p := 1/4*(1-xi)*(1+eta):
> N2p := 1/4*(1+xi)*(1+eta):
> N := linalg[matrix](1,4,[N1m,N2m,N1p,N2p]);
> #
> # Compute local laminar basis vectors at the blue point
> #
> # Isoparametric map
> #
> x := x1m*N1m + x2m*N2m + x1p*N1p + x2p*N2p;
> y := y1m*N1m + y2m*N2m + y1p*N1p + y2p*N2p;
> #
> # Derivative of Isoparametric map
> #
> dx_dxi := diff(x, xi);
> dy_dxi := diff(y, xi);
> #
> # Jacobian evluated at blue point
> #
> dx_dxi_cen := subs(eta=0, dx_dxi);
> dy_dxi_cen := subs(eta=0, dy_dxi);
> #
> # Local laminar basis vectors at the blue point
> #
> norm_dx_dxi_cen := simplify(sqrt(dx_dxi_cen^2 + dy_dxi_cen^2));
> ehat_x := vector([simplify(dx_dxi_cen/norm_dx_dxi_cen),
> simplify(dy_dxi_cen/norm_dx_dxi_cen), 0]);
> ehat_z := vector([0, 0, 1]);
> ehat_y := crossprod(ehat_z, ehat_x);

Nonlinear finite elements/Homework 6/Solutions/Problem 1