1 // Copyright (c) 2017-2022, Lawrence Livermore National Security, LLC and other CEED contributors. 2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details. 3 // 4 // SPDX-License-Identifier: BSD-2-Clause 5 // 6 // This file is part of CEED: http://github.com/ceed 7 8 /// @file 9 /// Boundary condition functions for solid mechanics example using PETSc 10 11 #include "../include/boundary.h" 12 13 // ----------------------------------------------------------------------------- 14 // Boundary Functions 15 // ----------------------------------------------------------------------------- 16 // Note: If additional boundary conditions are added, an update is needed in 17 // elasticity.h for the boundaryOptions variable. 18 19 // BCMMS - boundary function 20 // Values on all points of the mesh is set based on given solution below 21 // for u[0], u[1], u[2] 22 PetscErrorCode BCMMS(PetscInt dim, PetscReal load_increment, 23 const PetscReal coords[], PetscInt num_comp_u, 24 PetscScalar *u, void *ctx) { 25 PetscScalar x = coords[0]; 26 PetscScalar y = coords[1]; 27 PetscScalar z = coords[2]; 28 29 PetscFunctionBeginUser; 30 31 u[0] = exp(2*x)*sin(3*y)*cos(4*z) / 1e8 * load_increment; 32 u[1] = exp(3*y)*sin(4*z)*cos(2*x) / 1e8 * load_increment; 33 u[2] = exp(4*z)*sin(2*x)*cos(3*y) / 1e8 * load_increment; 34 35 PetscFunctionReturn(0); 36 }; 37 38 #ifndef M_PI 39 # define M_PI 3.14159265358979323846 40 #endif 41 42 // BCClamp - fix boundary values with affine transformation at fraction of load 43 // increment 44 PetscErrorCode BCClamp(PetscInt dim, PetscReal load_increment, 45 const PetscReal coords[], PetscInt num_comp_u, 46 PetscScalar *u, void *ctx) { 47 PetscScalar x = coords[0]; 48 PetscScalar y = coords[1]; 49 PetscScalar z = coords[2]; 50 PetscScalar (*clampMax) = (PetscScalar(*))ctx; 51 52 PetscFunctionBeginUser; 53 PetscScalar 54 // Translation vector 55 lx = clampMax[0]*load_increment, 56 ly = clampMax[1]*load_increment, 57 lz = clampMax[2]*load_increment, 58 // Normalized rotation axis 59 kx = clampMax[3], 60 ky = clampMax[4], 61 kz = clampMax[5], 62 // Rotation polynomial 63 c_0 = clampMax[6] * M_PI, 64 c_1 = clampMax[7] * M_PI, 65 cx = kx * x + ky * y + kz * z, 66 // Rotation magnitude 67 theta = (c_0 + c_1 * cx) * load_increment; 68 PetscScalar c = cos(theta), s = sin(theta); 69 70 u[0] = lx + s*(-kz*y + ky*z) + (1-c)*(-(ky*ky+kz*kz)*x + kx*ky*y + kx*kz*z); 71 u[1] = ly + s*(kz*x + -kx*z) + (1-c)*(kx*ky*x + -(kx*kx+kz*kz)*y + ky*kz*z); 72 u[2] = lz + s*(-ky*x + kx*y) + (1-c)*(kx*kz*x + ky*kz*y + -(kx*kx+ky*ky)*z); 73 PetscFunctionReturn(0); 74 }; 75