sadist.h (5222B)
1 /* Copyright (C) 2024 |Méso|Star> (contact@meso-star.com) 2 * 3 * This program is free software: you can redistribute it and/or modify 4 * it under the terms of the GNU General Public License as published by 5 * the Free Software Foundation, either version 3 of the License, or 6 * (at your option) any later version. 7 * 8 * This program is distributed in the hope that it will be useful, 9 * but WITHOUT ANY WARRANTY; without even the implied warranty of 10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11 * GNU General Public License for more details. 12 * 13 * You should have received a copy of the GNU General Public License 14 * along with this program. If not, see <http://www.gnu.org/licenses/>. */ 15 16 #ifndef SADIST_H 17 #define SADIST_H 18 19 #include <rsys/rsys.h> 20 #include <math.h> 21 #include <errno.h> /* EINVAL */ 22 23 /******************************************************************************* 24 * Trilinear profile 25 ******************************************************************************/ 26 struct sadist_trilinear_profile { 27 /* Spatial range in which the trilinear profile is defined */ 28 double lower[3]; 29 double upper[3]; 30 31 double a[2]; /* Interpolated values along X */ 32 double b[2]; /* Interpolated values along Y */ 33 double c[2]; /* Interpolated values along Z */ 34 }; 35 #define SADIST_TRILINEAR_PROFILE_NULL__ {{0,0,0}, {0,0,0}, {0,0}, {0,0}, {0,0}} 36 static const struct sadist_trilinear_profile SADIST_TRILINEAR_PROFILE_NULL = 37 SADIST_TRILINEAR_PROFILE_NULL__; 38 39 static INLINE int 40 sadist_trilinear_profile_check(const struct sadist_trilinear_profile* profile) 41 { 42 ASSERT(profile); 43 if(profile->lower[0] >= profile->upper[0] 44 || profile->lower[1] >= profile->upper[1] 45 || profile->lower[2] >= profile->upper[2]) 46 return EINVAL; 47 48 return 0; 49 } 50 51 static INLINE double 52 sadist_trilinear_profile_temperature 53 (const struct sadist_trilinear_profile* profile, 54 const double p[3]) 55 { 56 double u, v, w; 57 ASSERT(profile && p); 58 u = (p[0] - profile->lower[0]) / (profile->upper[0] - profile->lower[0]); 59 v = (p[1] - profile->lower[1]) / (profile->upper[1] - profile->lower[1]); 60 w = (p[2] - profile->lower[2]) / (profile->upper[2] - profile->lower[2]); 61 return u * (profile->a[1] - profile->a[0]) + profile->a[0] 62 + v * (profile->b[1] - profile->b[0]) + profile->b[0] 63 + w * (profile->c[1] - profile->c[0]) + profile->c[0]; 64 } 65 66 /******************************************************************************* 67 * Unsteady profile 68 ******************************************************************************/ 69 struct sadist_unsteady_profile { 70 double A; /* Influence of the thermal source */ 71 double B1; /* Influence of space variation */ 72 double B2; /* Influence of spatio-temporal variations */ 73 74 double kx; /* Spatio-temporal periodicity in X */ 75 double ky; /* Spatio-temporal periodicity in Y */ 76 double kz; /* Spatio-temporal periodicity in Z */ 77 78 double lambda; /* Thermal conductivity [W/m/K] */ 79 double rho; /* Volumic mass [kg/m^3] */ 80 double cp; /* Calorific capacity [J/K/kg] */ 81 }; 82 #define SADIST_UNSTEADY_PROFILE_NULL__ {0,0,0,0,0,0,0,0,0} 83 static const struct sadist_unsteady_profile SADIST_UNSTEADY_PROFILE_NULL = 84 SADIST_UNSTEADY_PROFILE_NULL__; 85 86 static INLINE double 87 sadist_unsteady_profile_temperature 88 (const struct sadist_unsteady_profile* profile, 89 const double pos[3], 90 const double time) 91 { 92 double alpha; /* Diffusivity */ 93 double kx, ky, kz; 94 double x, y, z, t; 95 double a, b, c; 96 double temp; 97 98 ASSERT(profile && pos); 99 100 x = pos[0]; 101 y = pos[1]; 102 z = pos[2]; 103 t = time; 104 105 alpha = profile->lambda / (profile->rho*profile->cp); 106 107 kx = profile->kx; 108 ky = profile->ky; 109 kz = profile->kz; 110 111 a = (x*x*x*z-3*x*y*y*z); 112 b = sin(kx*x)*sin(ky*y)*sin(kz*z)*exp(-alpha*(kx*kx + ky*ky + kz*kz)*t); 113 c = x*x*x*x * y*y*y * z*z; 114 115 temp = (profile->B1*a + profile->B2*b - profile->A*c) / profile->lambda; 116 return temp; 117 } 118 119 static INLINE double 120 sadist_unsteady_profile_power 121 (const struct sadist_unsteady_profile* profile, 122 const double pos[3]) 123 { 124 double x, y, z; 125 double p; 126 ASSERT(profile && pos); 127 128 x = pos[0]; 129 y = pos[1]; 130 z = pos[2]; 131 132 p = 12*x*x*y*y*y*z*z + 6*x*x*x*x*y*z*z + 2*x*x*x*x*y*y*y; 133 p = profile->A * p; 134 return p; 135 } 136 137 /******************************************************************************* 138 * Miscellaneous 139 ******************************************************************************/ 140 static INLINE void 141 sadist_write_stl 142 (FILE* stream, 143 const double* positions, 144 const size_t nvertices, 145 const size_t* indices, 146 const size_t ntriangles) 147 { 148 size_t itri; 149 (void)nvertices; 150 151 fprintf(stream, "solid ascii\n"); 152 FOR_EACH(itri, 0, ntriangles) { 153 const double* v0 = positions + indices[itri*3+0]*3; 154 const double* v1 = positions + indices[itri*3+1]*3; 155 const double* v2 = positions + indices[itri*3+2]*3; 156 fprintf(stream, "facet normal 0 1 0\n"); 157 fprintf(stream, "\touter loop\n"); 158 fprintf(stream, "\t\tvertex %g %g %g\n", SPLIT3(v0)); 159 fprintf(stream, "\t\tvertex %g %g %g\n", SPLIT3(v1)); 160 fprintf(stream, "\t\tvertex %g %g %g\n", SPLIT3(v2)); 161 fprintf(stream, "\tendloop\n"); 162 fprintf(stream, "endfacet\n"); 163 } 164 fprintf(stream, "endsolid\n"); 165 } 166 167 #endif /* SADIST_H */