commit e9865d6e6227c828c750502e7473d4dae5759232
parent fccb72868c926dec8b6458e80c328fb8dad9c3c5
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Tue, 17 Mar 2026 09:17:44 +0100
Revert "Update how the upper mesh of a line is calculated"
This reverts commit fe3d4473bae17697db5c6d87aa84b21dc5c2d8ee.
In reality, this change no longer guaranteed that the generated mesh was
indeed an upper boundary of the line. It merely ensured that its
_vertices_ lay above the line, but not necessarily the polyline itself;
in other words, some of its edges could intersect the line.
This reversal therefore aims to restore the previous algorithm, which
specifically ensured that the polyline was an upper bound of the line by
assigning a vertex the value of the one preceding it. As a reminder,
this works because only the upper half of the line is actually meshed,
and this upper half is a strictly decreasing function.
Diffstat:
| M | src/sln_line.c | | | 107 | +++++++++++++++++++------------------------------------------------------------ |
1 file changed, 26 insertions(+), 81 deletions(-)
diff --git a/src/sln_line.c b/src/sln_line.c
@@ -246,7 +246,7 @@ error:
/* Calculate line values for a set of wave numbers */
static res_T
-eval_mesh_fit
+eval_mesh
(const struct sln_tree* tree,
const struct line* line,
const struct darray_double* wavenumbers,
@@ -254,8 +254,7 @@ eval_mesh_fit
{
const double* nu = NULL;
double* ha = NULL;
- size_t ivertex = 0;
- size_t nvertices = 0;
+ size_t ivertex, nvertices;
res_T res = RES_OK;
ASSERT(tree && line && wavenumbers && values);
@@ -278,88 +277,34 @@ error:
goto exit;
}
-/* Calculate an upper bound for the line values for a set of wave numbers */
-static res_T
-eval_mesh_upper
- (const struct sln_tree* tree,
- const struct line* line,
- const struct darray_double* wavenumbers,
- /* #vertices used to mesh the line on its first interval, from the center of
- * the line to the width of the line at mid-height. This parameter is used to
- * calculate an upper bound mesh for the line in accordance with its width.
- * This user-defined discretization parameter controls the desired mesh
- * approximation. Thus, the mesh will be upper bounding “but not too much”
- * with respect to this numerical parameter */
- const size_t nvertices_around_line_center,
+static void
+snap_mesh_to_upper_bound
+ (const struct darray_double* wavenumbers,
struct darray_double* values)
{
- const double* nu = NULL;
double* ha = NULL;
- double nu_offset = 0;
- size_t ivertex = 0;
- size_t nvertices = 0;
- res_T res = RES_OK;
- ASSERT(tree && line && wavenumbers && values && nvertices_around_line_center);
-
- /* The line is meshed on its upper half, which is a strictly decreasing
- * function. To ensure that the mesh is an upper bound of the line, it is
- * therefore sufficient to evaluate its value at a wave number slightly lower
- * than the current wave number.
- *
- * This refers to the behavior of the following nu_offset. It corresponds to
- * an offset to be applied to the current nu so as to evaluate the line
- * slightly before the actual nu, and thus have a line value greater than its
- * actual value. */
- nu_offset = line->gamma_l / (double)nvertices_around_line_center;
-
- nvertices = darray_double_size_get(wavenumbers);
- ASSERT(nvertices);
-
- res = darray_double_resize(values, nvertices);
- if(res != RES_OK) goto error;
+ size_t ivertex, nvertices;
+ ASSERT(wavenumbers && values);
+ ASSERT(darray_double_size_get(wavenumbers) == darray_double_size_get(values));
+ (void)wavenumbers;
- nu = darray_double_cdata_get(wavenumbers);
ha = darray_double_data_get(values);
+ nvertices = darray_double_size_get(wavenumbers);
- FOR_EACH(ivertex, 0, nvertices) {
- double nu_adjusted = nu[ivertex];
-
- /* The first node is actually the center of the line, so it is the absolute
- * upper limit of the line. No offset should therefore be applied to its wave
- * number to ensure that the resulting mesh is an upper bound of the line */
- if(ivertex == 0) {
- ASSERT(line->wavenumber == nu_adjusted);
- } else {
- nu_adjusted -= nu_offset;
-
- if(nu_adjusted == nu[ivertex]) {
- /* Offset was too small regarding the numeric precision */
- nu_adjusted = nextafter(-DBL_MAX, nu_adjusted);
- nu_adjusted = MMIN(nu_adjusted, line->wavenumber);
- }
- }
-
- ha[ivertex] = line_value(tree, line, nu_adjusted);
-
-#ifndef NDEBUG
- if(ivertex != 0) {
- double ha_ref = line_value(tree, line, nu[ivertex]);
- ASSERT(ha_ref < ha[ivertex]);
- }
-#endif
-
- /* Ensure that the value of the vertex, stored in single precision, is
- * indeed an exclusive upper bound of the original value. */
- if(ha[ivertex] != (float)ha[ivertex]) {
- ha[ivertex] = nextafterf((float)ha[ivertex], FLT_MAX);
- }
+ /* Ensure that the stored vertex value is an exclusive upper bound of the
+ * original value. We do this by storing a value in single precision that is
+ * strictly greater than its encoding in double precision */
+ if(ha[0] != (float)ha[0]) {
+ ha[0] = nextafterf((float)ha[0], FLT_MAX);
}
-exit:
- return res;
-error:
- ERROR(tree->sln, "Error evaluating the line mesh -- %s.\n", res_to_cstr(res));
- goto exit;
+ /* We have meshed the upper half of the line which is a strictly decreasing
+ * function. To ensure that the mesh is an upper limit of this function,
+ * simply align the value of each vertex with the value of the preceding
+ * vertex */
+ FOR_EACH_REVERSE(ivertex, nvertices-1, 0) {
+ ha[ivertex] = ha[ivertex-1];
+ }
}
static INLINE int
@@ -589,13 +534,13 @@ line_mesh
/* Evaluate the mesh vertices, i.e. define the line value for the list of
* wavenumbers */
+ eval_mesh(tree, &line, &wavenumbers, &values);
+
switch(tree->args.mesh_type) {
case SLN_MESH_UPPER:
- eval_mesh_upper(tree, &line, &wavenumbers, nvertices_adjusted, &values);
- break;
- case SLN_MESH_FIT:
- eval_mesh_fit(tree, &line, &wavenumbers, &values);
+ snap_mesh_to_upper_bound(&wavenumbers, &values);
break;
+ case SLN_MESH_FIT: /* Do nothing */ break;
default: FATAL("Unreachable code.\n"); break;
}
