ME751 Advanced Computational Multibody Dynamics November 23, 2016 - PowerPoint PPT Presentation

ME751 Advanced Computational Multibody Dynamics November 23, 2016 Antonio Recuero University of Wisconsin-Madison

Before we get started… Last time, we learned:  Plates and shells  Introduce Chrono’s bilinear ANCF shell element  Kinematics of ANCF shell element  Curvilinear initial configuration  Inertia forces  Internal forces  This lecture…  Homework  Convergence  Locking Issues  Assumed Natural Strain  Enhanced Assumed Strain  Implementation in Chrono  Examples  2

4. Shell strains: Orthotropic and curvilinear reference 3

4. Shell strains: Orthotropic and curvilinear reference 4

4. Shell strains: Orthotropic and curvilinear reference Current deformed reference 5

4. Shell strains: Orthotropic and curvilinear reference 6

4. Shell strains: Orthotropic and curvilinear reference 7

4. Shell strains: Orthotropic and curvilinear reference Final expression! These strains are related to constitutive equations 8

4. Shell strains: Orthotropic and curvilinear reference 9

1. Locking 10

1. Locking 11

1. Locking issues: Convergence 12

1.a Volumetric Locking 13

1.a Volumetric Locking A uniform strain hexahedron and quadrilateral with hourglass control: Flanagan and Belytschko, 1981 Rigid body Strain modes 14 Hourglass modes

1.b Shear Locking 15

2. ANCF bilinear shell element 16

2. ANCF bilinear shell. ANS 17

2. ANCF bilinear shell. ANS Substitution of these compatible strains in equations 18

2. ANCF bilinear shell. ANS 19

2. ANCF bilinear shell. EAS 20

2. ANCF bilinear shell. EAS 21

2. ANCF bilinear shell. EAS 22

2. ANCF bilinear shell. EAS Substitution of thickness strain 23

3. Implementation of bilinear ANCF shell element in Chrono void MyForce::Evaluate(ChMatrixNM<double, 54, 1>& result, // Transformation : Orthogonal transformation (A and J) const double x, const double y, const double z) { ChVector<double> G1xG2; // Cross product of first and // Element shape function Evaluation of 3D shape second column of ChMatrixNM<double, 1, 8> N; functions double G1dotG1; // Dot product of first column of m_element ‐ >ShapeFunctions(N, x, y, z); position vector gradient // Determinant of position vector gradient matrix: Initial G1xG2.x = Nx_d0[0][1] * Ny_d0[0][2] ‐ Nx_d0[0][2] * configuration Ny_d0[0][1]; Shape function ChMatrixNM<double, 1, 8> Nx; G1xG2.y = Nx_d0[0][2] * Ny_d0[0][0] ‐ Nx_d0[0][0] * derivatives ChMatrixNM<double, 1, 8> Ny; Ny_d0[0][2]; ChMatrixNM<double, 1, 8> Nz; G1xG2.z = Nx_d0[0][0] * Ny_d0[0][1] ‐ Nx_d0[0][1] * Vectors of initial ChMatrixNM<double, 1, 3> Nx_d0; Ny_d0[0][0]; deformation gradient ChMatrixNM<double, 1, 3> Ny_d0; G1dotG1 = Nx_d0[0][0] * Nx_d0[0][0] + Nx_d0[0][1] * ChMatrixNM<double, 1, 3> Nz_d0; Nx_d0[0][1] + Nx_d0[0][2] * Nx_d0[0][2]; double detJ0 = m_element ‐ >Calc_detJ0(x, y, z, Nx, Ny, Nz, Nx_d0, Ny_d0, Nz_d0); New interpolation for ANS // Tangent Frame Covariant vector base, initial config. ChVector<double> A1; // ANS shape function ChVector<double> A2; ChMatrixNM<double, 1, 4> S_ANS; // Shape function vector ChVector<double> A3; Normalization for Assumed Natural Strain A1.x = Nx_d0[0][0]; of such a frame ChMatrixNM<double, 6, 5> M; // Shape function vector A1.y = Nx_d0[0][1]; for Enhanced Assumed Strain A1.z = Nx_d0[0][2]; m_element ‐ >ShapeFunctionANSbilinearShell(S_ANS, x, y); A1 = A1 / sqrt(G1dotG1); m_element ‐ >Basis_M(M, x, y, z); A3 = G1xG2.GetNormalized(); Interpolation for 5 internal A2.Cross(A3, A1); parameters 24

3. Internal forces j01(0) = j0(0, 0); j02(0) = j0(1, 0); j03(0) = j0(2, 0); . // Direction for orthotropic material . double theta = m_element ‐ >GetLayer(m_kl).Get_theta(); // . Fiber angle j02(2) = j0(1, 2); ChVector<double> AA1; j03(2) = j0(2, 2); ChVector<double> AA2; ChVector<double> AA3; // Coefficients of contravariant transformation AA1 = A1 * cos(theta) + A2 * sin(theta); beta(0, 0) = Vdot(AA1, j01); AA2 = ‐ A1 * sin(theta) + A2 * cos(theta); beta(1, 0) = Vdot(AA2, j01); Calculate beta coefficients Transform orthogonal AA3 = A3; beta(2, 0) = Vdot(AA3, j01); system with direction of . fiber (orthotropic material) /// Beta . ChMatrixNM<double, 3, 3> j0; . Constant transformation to ChVector<double> j01; beta(7, 0) = Vdot(AA2, j03); account for initial ChVector<double> j02; beta(8, 0) = Vdot(AA3, j03); configuration ChVector<double> j03; ChMatrixNM<double, 9, 1> beta; // Transformation matrix, function of fiber angle // Calculates inverse of rd0 (j0) (position vector const ChMatrixNM<double, 6, 6>& T0 = m_element ‐ gradient: Initial Configuration) >GetLayer(m_kl).Get_T0(); j0(0, 0) = Ny_d0[0][1] * Nz_d0[0][2] ‐ Nz_d0[0][1] * // Determinant of the initial position vector gradient at Ny_d0[0][2]; the element center j0(0, 1) = Ny_d0[0][2] * Nz_d0[0][0] ‐ Ny_d0[0][0] * double detJ0C = m_element ‐ >GetLayer(m_kl).Get_detJ0C(); Nz_d0[0][2]; j0(0, 2) = Ny_d0[0][0] * Nz_d0[0][1] ‐ Nz_d0[0][0] * // Enhanced Assumed Strain Ny_d0[0][1]; ChMatrixNM<double, 6, 5> G = T0 * M * (detJ0C / detJ0); j0(1, 0) = Nz_d0[0][1] * Nx_d0[0][2] ‐ Nx_d0[0][1] * ChMatrixNM<double, 6, 1> strain_EAS = G * (*m_alpha_eas); . . ChMatrixNM<double, 8, 1> ddNx; EAS calculation . ChMatrixNM<double, 8, 1> ddNy; ChMatrixNM<double, 8, 1> ddNz; j0.MatrDivScale(detJ0); ddNx.MatrMultiplyT(m_element ‐ >m_ddT, Nx); ddNy.MatrMultiplyT(m_element ‐ >m_ddT, Ny); 25 ddNz.MatrMultiplyT(m_element ‐ >m_ddT, Nz);

3. Internal forces ChMatrixNM<double, 8, 1> ddNx; ChMatrixNM<double, 8, 1> ddNy; ChMatrixNM<double, 8, 1> ddNz; ddNx.MatrMultiplyT(m_element ‐ >m_ddT, Nx); ddNy.MatrMultiplyT(m_element ‐ >m_ddT, Ny); ddNz.MatrMultiplyT(m_element ‐ >m_ddT, Nz); ChMatrixNM<double, 8, 1> d0d0Nx; ChMatrixNM<double, 8, 1> d0d0Ny; ChMatrixNM<double, 8, 1> d0d0Nz; d0d0Nx.MatrMultiplyT(m_element ‐ >m_d0d0T, Nx); Derivatives w.r.t. coordinates d0d0Ny.MatrMultiplyT(m_element ‐ >m_d0d0T, Ny); d0d0Nz.MatrMultiplyT(m_element ‐ >m_d0d0T, Nz); 26

3. Internal forces Calculation of covariant compatible strains // Strain component ChMatrixNM<double, 6, 1> strain_til; strain_til(0, 0) = 0.5 * ((Nx * ddNx)(0, 0) ‐ (Nx * d0d0Nx)(0, 0)); strain_til(1, 0) = 0.5 * ((Ny * ddNy)(0, 0) ‐ (Ny * d0d0Ny)(0, 0)); strain_til(2, 0) = (Nx * ddNy)(0, 0) ‐ (Nx * d0d0Ny)(0, 0); strain_til(3, 0) = N(0, 0) * m_element ‐ >m_strainANS(0, 0) + N(0, 2) * m_element ‐ >m_strainANS(1, 0) + N(0, 4) * m_element ‐ >m_strainANS(2, 0) + N(0, 6) * m_element ‐ >m_strainANS(3, 0); strain_til(4, 0) = S_ANS(0, 2) * m_element ‐ >m_strainANS(6, 0) + S_ANS(0, 3) * m_element ‐ >m_strainANS(7, 0); strain_til(5, 0) = S_ANS(0, 0) * m_element ‐ >m_strainANS(4, 0) + S_ANS(0, 1) * m_element ‐ >m_strainANS(5, 0); ANS: Transverse shear and thickness strain // For orthotropic material ChMatrixNM<double, 6, 1> strain; strain(0, 0) = strain_til(0, 0) * beta(0) * beta(0) + strain_til(1, 0) * beta(3) * beta(3) + strain_til(2, 0) * beta(0) * beta(3) + strain_til(3, 0) * beta(6) * beta(6) + strain_til(4, 0) * beta(0) * beta(6) + strain_til(5, 0) * beta(3) * beta(6); . . . strain(5, 0) = strain_til(0, 0) * 2.0 * beta(1) * beta(2) + strain_til(1, 0) * 2.0 * beta(4) * beta(5) + strain_til(2, 0) * (beta(2) * beta(4) + beta(1) * beta(5)) + strain_til(3, 0) * 2.0 * beta(7) * beta(8) + Beta matrix/vector multiplication strain_til(4, 0) * (beta(2) * beta(7) + beta(1) * beta(8)) + strain_til(5, 0) * (beta(5) * beta(7) + beta(4) * beta(8)); 27