1 // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 4 #include <ceed/types.h> 5 #include "utils.h" 6 7 /** 8 @brief Calculate CFL for 3D spatial dimensions 9 10 Defined as 11 @f[ 12 \mathrm{CFL} = \Delta t \sqrt{\overline{g}_{jk} u_j u_k} 13 @f] 14 15 where \f$ \Delta t \f$ is the `timestep`, \f$ \overline{g}_{jk} \f$ is the element metric tensor `gij`, and \f$ u_j \f$ is the `velocity`. 16 For best results, `gij` should be scaled to account for the reference element domain (often [-1,1] for Gaussian quadrature). 17 18 @param[in] velocity Velocity 19 @param[in] timestep Timestep 20 @param[in] gij Element metric tensor (should be scaled) 21 @return The CFL value for the element 22 **/ 23 CEED_QFUNCTION_HELPER CeedScalar CalculateCFL_3D(const CeedScalar velocity[3], CeedScalar timestep, const CeedScalar gij[3][3]) { 24 CeedScalar gij_uj[3] = {0.}; 25 26 MatVec3(gij, velocity, CEED_NOTRANSPOSE, gij_uj); 27 return sqrt(Dot3(velocity, gij_uj)) * timestep; 28 } 29 30 /** 31 @brief Calculate CFL for 2D spatial dimensions 32 33 Defined as 34 @f[ 35 \mathrm{CFL} = \Delta t \sqrt{\overline{g}_{jk} u_j u_k} 36 @f] 37 38 where \f$ \Delta t \f$ is the `timestep`, \f$ \overline{g}_{jk} \f$ is the element metric tensor `gij`, and \f$ u_j \f$ is the `velocity`. 39 For best results, `gij` should be scaled to account for the reference element domain (often [-1,1] for Gaussian quadrature). 40 41 @param[in] velocity Advection velocity 42 @param[in] timestep Timestep 43 @param[in] gij Element metric tensor (should be scaled) 44 @return The CFL value for the element 45 **/ 46 CEED_QFUNCTION_HELPER CeedScalar CalculateCFL_2D(const CeedScalar velocity[2], CeedScalar timestep, const CeedScalar gij[2][2]) { 47 CeedScalar gij_uj[2] = {0.}; 48 49 MatVec2(gij, velocity, CEED_NOTRANSPOSE, gij_uj); 50 return sqrt(Dot2(velocity, gij_uj)) * timestep; 51 } 52 53 /** 54 @brief Calculate Peclet number for 3D spatial dimensions 55 56 Defined as 57 @f[ 58 \mathrm{Pe} = \frac{\sqrt{\overline{g}_{jk} u_j u_k}}{\kappa} 59 @f] 60 61 where \f$ \kappa \f$ is the `diffusion_coeff`, \f$ \overline{g}_{jk} \f$ is the element metric tensor `gij`, and \f$ u_j \f$ is the `velocity`. 62 For best results, `gij` should be scaled to account for the reference element domain (often [-1,1] for Gaussian quadrature). 63 64 @param[in] velocity Advection velocity 65 @param[in] diffusion_coeff Diffusion coefficient 66 @param[in] gij Element metric tensor (should be scaled) 67 @return The Peclet value for the element 68 **/ 69 CEED_QFUNCTION_HELPER CeedScalar CalculatePe_3D(const CeedScalar velocity[3], CeedScalar diffusion_coeff, const CeedScalar gij[3][3]) { 70 CeedScalar gij_uj[3] = {0.}, gij_inv_mat[3][3] = {{0.}}; 71 72 MatInv3(gij, gij_inv_mat, NULL); 73 MatVec3(gij_inv_mat, velocity, CEED_NOTRANSPOSE, gij_uj); 74 return sqrt(Dot3(velocity, gij_uj)) / diffusion_coeff; 75 } 76 77 /** 78 @brief Calculate Peclet number for 2D spatial dimensions 79 80 Defined as 81 @f[ 82 \mathrm{Pe} = \frac{\sqrt{\overline{g}_{jk} u_j u_k}}{\kappa} 83 @f] 84 85 where \f$ \kappa \f$ is the `diffusion_coeff`, \f$ \overline{g}_{jk} \f$ is the element metric tensor `gij`, and \f$ u_j \f$ is the `velocity`. 86 For best results, `gij` should be scaled to account for the reference element domain (often [-1,1] for Gaussian quadrature). 87 88 @param[in] velocity Advection velocity 89 @param[in] diffusion_coeff Diffusion coefficient 90 @param[in] gij Element metric tensor (should be scaled) 91 @return The Peclet value for the element 92 **/ 93 CEED_QFUNCTION_HELPER CeedScalar CalculatePe_2D(const CeedScalar velocity[2], CeedScalar diffusion_coeff, const CeedScalar gij[2][2]) { 94 CeedScalar gij_uj[2] = {0.}, gij_inv_mat[2][2] = {{0.}}; 95 96 MatInv2(gij, gij_inv_mat, NULL); 97 MatVec2(gij_inv_mat, velocity, CEED_NOTRANSPOSE, gij_uj); 98 return sqrt(Dot2(velocity, gij_uj)) / diffusion_coeff; 99 } 100