This is an example of cs_user_extra_operations function which computes a local Nusselt number.
Local variables
Faces can be selected either using a selector (as is done here) or directly using boundary zone element ids if an appropriate zone has been defined.
FILE *file = nullptr;
const cs_real_t prandtl = 0.71, height = 1., qwall = 1.;
("normal[0, -1, 0, 0.1] and box[-1000, -1000, -1000, 1000, 0.01, 1000]",
&nlelt, lstelt);
Reconstruct value at selected boundary faces
Allocate a local array for the selected boundary faces.
General case (for non-orthogonal meshes)
Case of orthogonal meshes
Compute boundary value without reconstruction
for (
cs_lnum_t ielt = 0; ielt < nlelt; ielt++) {
treloc[ielt] = coefap[face_id] + coefbp[face_id] * f->
val[c_id];
}
- Note
- Here, we assign the non-reconstructed value.
Open file to print values and broadcast values to all parallel ranks. Values are ordered by the xabs
array values provided to the cs_parall_allgather_ordered_r, so this function is needed even when not running in parallel.
file = fopen("Nusselt.dat", "w");
}
cs_real_t *xabs =
nullptr, *xabsg =
nullptr, *xnusselt =
nullptr;
for (
cs_lnum_t ilelt = 0; ilelt < nlelt; ilelt ++) {
xabs[ilelt] =
cdgfbo[face_id][0];
}
Compute the bulk temperature, finalize the Nusselt number, print it and free memory not already freed before:
for (
cs_gnum_t ilelt = 0; ilelt < neltg; ilelt ++) {
int npoint = 200;
int irang1 = -1;
int irangv;
for (int ii = 0; ii < npoint; ii++) {
0.};
cell_cen,
xyz,
&c_id, &irangv);
if (c_id != c_id_prev || irangv != irang1) {
c_id_prev = c_id;
irang1 = irangv;
xtbulk += xtb;
xubulk += xub;
}
}
}
xtbulk = xbulk[0] / xbulk[1];
xubulk = xbulk[1];
xnusselt[ilelt] = qwall * 2. * height /
lambda / (tfac - xtbulk);
}
fprintf(file,
"%17.9e %17.9e\n",
xabsg[ii]*10.,
xnusselt[ii]/(0.023*pow(30000., 0.8)*pow(0.71, 0.4)));
if (file != nullptr)
fclose(file);