1 // Copyright (c) 2017-2026, 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 #include <ceed.h>
14 #include <petscsys.h>
15
16 // -----------------------------------------------------------------------------
17 // Boundary Functions
18 // -----------------------------------------------------------------------------
19 // Note: If additional boundary conditions are added, an update is needed in elasticity.h for the boundaryOptions variable.
20
21 // BCMMS - boundary function
22 // Values on all points of the mesh is set based on given solution below for u[0], u[1], u[2]
BCMMS(PetscInt dim,PetscReal load_increment,const PetscReal coords[],PetscInt num_comp_u,PetscScalar * u,void * ctx)23 PetscErrorCode BCMMS(PetscInt dim, PetscReal load_increment, const PetscReal coords[], PetscInt num_comp_u, PetscScalar *u, void *ctx) {
24 PetscScalar x = coords[0];
25 PetscScalar y = coords[1];
26 PetscScalar z = coords[2];
27
28 PetscFunctionBeginUser;
29
30 u[0] = exp(2 * x) * sin(3 * y) * cos(4 * z) / 1e8 * load_increment;
31 u[1] = exp(3 * y) * sin(4 * z) * cos(2 * x) / 1e8 * load_increment;
32 u[2] = exp(4 * z) * sin(2 * x) * cos(3 * y) / 1e8 * load_increment;
33
34 PetscFunctionReturn(PETSC_SUCCESS);
35 };
36
37 #ifndef M_PI
38 #define M_PI 3.14159265358979323846
39 #endif
40
41 // BCClamp - fix boundary values with affine transformation at fraction of load increment
BCClamp(PetscInt dim,PetscReal load_increment,const PetscReal coords[],PetscInt num_comp_u,PetscScalar * u,void * ctx)42 PetscErrorCode BCClamp(PetscInt dim, PetscReal load_increment, const PetscReal coords[], PetscInt num_comp_u, PetscScalar *u, void *ctx) {
43 PetscScalar x = coords[0];
44 PetscScalar y = coords[1];
45 PetscScalar z = coords[2];
46 PetscScalar(*clampMax) = (PetscScalar(*))ctx;
47
48 PetscFunctionBeginUser;
49 PetscScalar
50 // Translation vector
51 lx = clampMax[0] * load_increment,
52 ly = clampMax[1] * load_increment, lz = clampMax[2] * load_increment,
53 // Normalized rotation axis
54 kx = clampMax[3], ky = clampMax[4], kz = clampMax[5],
55 // Rotation polynomial
56 c_0 = clampMax[6] * M_PI, c_1 = clampMax[7] * M_PI, cx = kx * x + ky * y + kz * z,
57 // Rotation magnitude
58 theta = (c_0 + c_1 * cx) * load_increment;
59 PetscScalar c = cos(theta), s = sin(theta);
60
61 u[0] = lx + s * (-kz * y + ky * z) + (1 - c) * (-(ky * ky + kz * kz) * x + kx * ky * y + kx * kz * z);
62 u[1] = ly + s * (kz * x + -kx * z) + (1 - c) * (kx * ky * x + -(kx * kx + kz * kz) * y + ky * kz * z);
63 u[2] = lz + s * (-ky * x + kx * y) + (1 - c) * (kx * kz * x + ky * kz * y + -(kx * kx + ky * ky) * z);
64 PetscFunctionReturn(PETSC_SUCCESS);
65 };
66