BMesh: simplify normal calculation, resolve partial update error

Simplify vertex normal calculation by moving the main normal
accumulation function to operate on vertices instead of faces.

Using faces had the down side that it needed to zero, accumulate and
normalize the vertex normals in 3 separate passes, accumulating also
needed a spin-lock for thread since the face would write it's normal
to all of it's vertices which could be shared with other faces.

Now a single loop over vertices is performed without locking.
This gives 5-6% speedup calculating all normals.

This also simplifies partial updates, fixing a problem where
all connected faces were being read from when calculating normals.
While this could have been resolved separately,
it's simpler to operate on vertices directly.
This commit is contained in:
Campbell Barton 2021-06-08 15:54:21 +10:00
parent f651cc6c4e
commit 496045fc30
Notes: blender-bot 2023-02-14 10:18:56 +01:00
Referenced by issue #88550, Mesh Optimization Project Progress
3 changed files with 167 additions and 211 deletions

View File

@ -35,8 +35,6 @@
#include "BKE_global.h"
#include "BKE_mesh.h"
#include "atomic_ops.h"
#include "intern/bmesh_private.h"
/* -------------------------------------------------------------------- */
@ -59,19 +57,28 @@ typedef struct BMEdgesCalcVectorsData {
float (*edgevec)[3];
} BMEdgesCalcVectorsData;
static void mesh_edges_calc_vectors_cb(void *userdata, MempoolIterData *mp_e)
static void bm_edge_calc_vectors_cb(void *userdata, MempoolIterData *mp_e)
{
BMEdgesCalcVectorsData *data = userdata;
BMEdge *e = (BMEdge *)mp_e;
if (e->l) {
const float *v1_co = data->vcos ? data->vcos[BM_elem_index_get(e->v1)] : e->v1->co;
const float *v2_co = data->vcos ? data->vcos[BM_elem_index_get(e->v2)] : e->v2->co;
sub_v3_v3v3(data->edgevec[BM_elem_index_get(e)], v2_co, v1_co);
normalize_v3(data->edgevec[BM_elem_index_get(e)]);
/* The edge vector will not be needed when the edge has no radial. */
if (e->l != NULL) {
float(*edgevec)[3] = userdata;
float *e_diff = edgevec[BM_elem_index_get(e)];
sub_v3_v3v3(e_diff, e->v2->co, e->v1->co);
normalize_v3(e_diff);
}
else {
/* the edge vector will not be needed when the edge has no radial */
}
static void bm_edge_calc_vectors_with_coords_cb(void *userdata, MempoolIterData *mp_e)
{
BMEdge *e = (BMEdge *)mp_e;
/* The edge vector will not be needed when the edge has no radial. */
if (e->l != NULL) {
BMEdgesCalcVectorsData *data = userdata;
float *e_diff = data->edgevec[BM_elem_index_get(e)];
sub_v3_v3v3(
e_diff, data->vcos[BM_elem_index_get(e->v2)], data->vcos[BM_elem_index_get(e->v1)]);
normalize_v3(e_diff);
}
}
@ -79,118 +86,123 @@ static void bm_mesh_edges_calc_vectors(BMesh *bm, float (*edgevec)[3], const flo
{
BM_mesh_elem_index_ensure(bm, BM_EDGE | (vcos ? BM_VERT : 0));
BMEdgesCalcVectorsData data = {
.vcos = vcos,
.edgevec = edgevec,
};
BM_iter_parallel(
bm, BM_EDGES_OF_MESH, mesh_edges_calc_vectors_cb, &data, bm->totedge >= BM_OMP_LIMIT);
if (vcos == NULL) {
BM_iter_parallel(
bm, BM_EDGES_OF_MESH, bm_edge_calc_vectors_cb, edgevec, bm->totedge >= BM_OMP_LIMIT);
}
else {
BMEdgesCalcVectorsData data = {
.edgevec = edgevec,
.vcos = vcos,
};
BM_iter_parallel(bm,
BM_EDGES_OF_MESH,
bm_edge_calc_vectors_with_coords_cb,
&data,
bm->totedge >= BM_OMP_LIMIT);
}
}
typedef struct BMVertsCalcNormalsData {
typedef struct BMVertsCalcNormalsWithCoordsData {
/* Read-only data. */
const float (*fnos)[3];
const float (*edgevec)[3];
const float (*vcos)[3];
/* Read-write data, protected by an atomic-based fake spin-lock like system. */
/* Write data. */
float (*vnos)[3];
} BMVertsCalcNormalsData;
} BMVertsCalcNormalsWithCoordsData;
static void mesh_verts_calc_normals_accum(
BMFace *f,
const float *f_no,
const float (*edgevec)[3],
/* Read-write data, protected by an atomic-based fake spin-lock like system. */
float (*vnos)[3])
BLI_INLINE void bm_vert_calc_normals_accum_loop(const BMLoop *l_iter,
const float (*edgevec)[3],
const float f_no[3],
float v_no[3])
{
#define FLT_EQ_NONAN(_fa, _fb) (*((const uint32_t *)&_fa) == *((const uint32_t *)&_fb))
BMLoop *l_first, *l_iter;
l_iter = l_first = BM_FACE_FIRST_LOOP(f);
do {
const float *e1diff, *e2diff;
float dotprod;
float fac;
/* calculate the dot product of the two edges that
* meet at the loop's vertex */
e1diff = edgevec[BM_elem_index_get(l_iter->prev->e)];
e2diff = edgevec[BM_elem_index_get(l_iter->e)];
dotprod = dot_v3v3(e1diff, e2diff);
/* edge vectors are calculated from e->v1 to e->v2, so
* adjust the dot product if one but not both loops
* actually runs from from e->v2 to e->v1 */
if ((l_iter->prev->e->v1 == l_iter->prev->v) ^ (l_iter->e->v1 == l_iter->v)) {
dotprod = -dotprod;
}
fac = saacos(-dotprod);
if (fac != fac) { /* NAN detection. */
/* Degenerated case, nothing to do here, just ignore that vertex. */
continue;
}
/* accumulate weighted face normal into the vertex's normal */
float *v_no = vnos ? vnos[BM_elem_index_get(l_iter->v)] : l_iter->v->no;
/* This block is a lockless threadsafe madd_v3_v3fl.
* It uses the first float of the vector as a sort of cheap spin-lock,
* assuming FLT_MAX is a safe 'illegal' value that cannot be set here otherwise.
* It also assumes that collisions between threads are highly unlikely,
* else performances would be quite bad here. */
float virtual_lock = v_no[0];
while (true) {
/* This loops until following conditions are met:
* - v_no[0] has same value as virtual_lock (i.e. it did not change since last try).
* - v_no[0] was not FLT_MAX, i.e. it was not locked by another thread.
*/
const float vl = atomic_cas_float(&v_no[0], virtual_lock, FLT_MAX);
if (FLT_EQ_NONAN(vl, virtual_lock) && vl != FLT_MAX) {
break;
}
virtual_lock = vl;
}
BLI_assert(v_no[0] == FLT_MAX);
/* Now we own that normal value, and can change it.
* But first scalar of the vector must not be changed yet, it's our lock! */
virtual_lock += f_no[0] * fac;
v_no[1] += f_no[1] * fac;
v_no[2] += f_no[2] * fac;
/* Second atomic operation to 'release'
* our lock on that vector and set its first scalar value. */
/* Note that we do not need to loop here, since we 'locked' v_no[0],
* nobody should have changed it in the mean time. */
virtual_lock = atomic_cas_float(&v_no[0], FLT_MAX, virtual_lock);
BLI_assert(virtual_lock == FLT_MAX);
} while ((l_iter = l_iter->next) != l_first);
#undef FLT_EQ_NONAN
}
static void mesh_verts_calc_normals_accum_cb(void *userdata, MempoolIterData *mp_f)
{
BMVertsCalcNormalsData *data = userdata;
BMFace *f = (BMFace *)mp_f;
const float *f_no = data->fnos ? data->fnos[BM_elem_index_get(f)] : f->no;
mesh_verts_calc_normals_accum(f, f_no, data->edgevec, data->vnos);
}
static void mesh_verts_calc_normals_normalize_cb(void *userdata, MempoolIterData *mp_v)
{
BMVertsCalcNormalsData *data = userdata;
BMVert *v = (BMVert *)mp_v;
float *v_no = data->vnos ? data->vnos[BM_elem_index_get(v)] : v->no;
if (UNLIKELY(normalize_v3(v_no) == 0.0f)) {
const float *v_co = data->vcos ? data->vcos[BM_elem_index_get(v)] : v->co;
normalize_v3_v3(v_no, v_co);
/* Calculate the dot product of the two edges that meet at the loop's vertex. */
const float *e1diff = edgevec[BM_elem_index_get(l_iter->prev->e)];
const float *e2diff = edgevec[BM_elem_index_get(l_iter->e)];
/* Edge vectors are calculated from e->v1 to e->v2, so adjust the dot product if one but not
* both loops actually runs from from e->v2 to e->v1. */
float dotprod = dot_v3v3(e1diff, e2diff);
if ((l_iter->prev->e->v1 == l_iter->prev->v) ^ (l_iter->e->v1 == l_iter->v)) {
dotprod = -dotprod;
}
const float fac = saacos(-dotprod);
/* NAN detection, otherwise this is a degenerated case, ignore that vertex in this case. */
if (fac == fac) { /* NAN detection. */
madd_v3_v3fl(v_no, f_no, fac);
}
}
static void bm_vert_calc_normals_impl(const float (*edgevec)[3], BMVert *v)
{
float *v_no = v->no;
zero_v3(v_no);
BMEdge *e_first = v->e;
if (e_first != NULL) {
BMEdge *e_iter = e_first;
do {
BMLoop *l_first = e_iter->l;
if (l_first != NULL) {
BMLoop *l_iter = l_first;
do {
if (l_iter->v == v) {
bm_vert_calc_normals_accum_loop(l_iter, edgevec, l_iter->f->no, v_no);
}
} while ((l_iter = l_iter->radial_next) != l_first);
}
} while ((e_iter = BM_DISK_EDGE_NEXT(e_iter, v)) != e_first);
if (LIKELY(normalize_v3(v_no) != 0.0f)) {
return;
}
}
/* Fallback normal. */
normalize_v3_v3(v_no, v->co);
}
static void bm_vert_calc_normals_cb(void *userdata, MempoolIterData *mp_v)
{
const float(*edgevec)[3] = userdata;
BMVert *v = (BMVert *)mp_v;
bm_vert_calc_normals_impl(edgevec, v);
}
static void bm_vert_calc_normals_with_coords(BMVert *v, BMVertsCalcNormalsWithCoordsData *data)
{
float *v_no = data->vnos[BM_elem_index_get(v)];
zero_v3(v_no);
/* Loop over edges. */
BMEdge *e_first = v->e;
if (e_first != NULL) {
BMEdge *e_iter = e_first;
do {
BMLoop *l_first = e_iter->l;
if (l_first != NULL) {
BMLoop *l_iter = l_first;
do {
if (l_iter->v == v) {
bm_vert_calc_normals_accum_loop(
l_iter, data->edgevec, data->fnos[BM_elem_index_get(l_iter->f)], v_no);
}
} while ((l_iter = l_iter->radial_next) != l_first);
}
} while ((e_iter = BM_DISK_EDGE_NEXT(e_iter, v)) != e_first);
if (LIKELY(normalize_v3(v_no) != 0.0f)) {
return;
}
}
/* Fallback normal. */
normalize_v3_v3(v_no, data->vcos[BM_elem_index_get(v)]);
}
static void bm_vert_calc_normals_with_coords_cb(void *userdata, MempoolIterData *mp_v)
{
BMVertsCalcNormalsWithCoordsData *data = userdata;
BMVert *v = (BMVert *)mp_v;
bm_vert_calc_normals_with_coords(v, data);
}
static void bm_mesh_verts_calc_normals(BMesh *bm,
@ -201,29 +213,34 @@ static void bm_mesh_verts_calc_normals(BMesh *bm,
{
BM_mesh_elem_index_ensure(bm, (BM_EDGE | BM_FACE) | ((vnos || vcos) ? BM_VERT : 0));
BMVertsCalcNormalsData data = {
.fnos = fnos,
.edgevec = edgevec,
.vcos = vcos,
.vnos = vnos,
};
BM_iter_parallel(
bm, BM_FACES_OF_MESH, mesh_verts_calc_normals_accum_cb, &data, bm->totface >= BM_OMP_LIMIT);
/* normalize the accumulated vertex normals */
BM_iter_parallel(bm,
BM_VERTS_OF_MESH,
mesh_verts_calc_normals_normalize_cb,
&data,
bm->totvert >= BM_OMP_LIMIT);
if (vcos == NULL) {
BM_iter_parallel(bm,
BM_VERTS_OF_MESH,
bm_vert_calc_normals_cb,
(void *)edgevec,
bm->totvert >= BM_OMP_LIMIT);
}
else {
BLI_assert(!ELEM(NULL, fnos, vnos));
BMVertsCalcNormalsWithCoordsData data = {
.edgevec = edgevec,
.fnos = fnos,
.vcos = vcos,
.vnos = vnos,
};
BM_iter_parallel(bm,
BM_VERTS_OF_MESH,
bm_vert_calc_normals_with_coords_cb,
&data,
bm->totvert >= BM_OMP_LIMIT);
}
}
static void mesh_faces_calc_normals_cb(void *UNUSED(userdata), MempoolIterData *mp_f)
static void bm_face_calc_normals_cb(void *UNUSED(userdata), MempoolIterData *mp_f)
{
BMFace *f = (BMFace *)mp_f;
BM_face_normal_update(f);
BM_face_calc_normal(f, f->no);
}
/**
@ -235,27 +252,13 @@ void BM_mesh_normals_update(BMesh *bm)
{
float(*edgevec)[3] = MEM_mallocN(sizeof(*edgevec) * bm->totedge, __func__);
/* Parallel mempool iteration does not allow generating indices inline anymore... */
/* Parallel mempool iteration does not allow generating indices inline anymore. */
BM_mesh_elem_index_ensure(bm, (BM_EDGE | BM_FACE));
/* calculate all face normals */
/* Calculate all face normals. */
BM_iter_parallel(
bm, BM_FACES_OF_MESH, mesh_faces_calc_normals_cb, NULL, bm->totface >= BM_OMP_LIMIT);
bm, BM_FACES_OF_MESH, bm_face_calc_normals_cb, NULL, bm->totface >= BM_OMP_LIMIT);
/* Zero out vertex normals */
BMIter viter;
BMVert *v;
int i;
BM_ITER_MESH_INDEX (v, &viter, bm, BM_VERTS_OF_MESH, i) {
BM_elem_index_set(v, i); /* set_inline */
zero_v3(v->no);
}
bm->elem_index_dirty &= ~BM_VERT;
/* Compute normalized direction vectors for each edge.
* Directions will be used for calculating the weights of the face normals on the vertex normals.
*/
bm_mesh_edges_calc_vectors(bm, edgevec, NULL);
/* Add weighted face normals to vertices, and normalize vert normals. */
@ -269,14 +272,14 @@ void BM_mesh_normals_update(BMesh *bm)
/** \name Update Vertex & Face Normals (Partial Updates)
* \{ */
static void mesh_faces_parallel_range_calc_normals_cb(
static void bm_partial_faces_parallel_range_calc_normals_cb(
void *userdata, const int iter, const TaskParallelTLS *__restrict UNUSED(tls))
{
BMFace *f = ((BMFace **)userdata)[iter];
BM_face_normal_update(f);
BM_face_calc_normal(f, f->no);
}
static void mesh_edges_parallel_range_calc_vectors_cb(
static void bm_partial_edges_parallel_range_calc_vectors_cb(
void *userdata, const int iter, const TaskParallelTLS *__restrict UNUSED(tls))
{
BMEdge *e = ((BMEdge **)((void **)userdata)[0])[iter];
@ -285,21 +288,12 @@ static void mesh_edges_parallel_range_calc_vectors_cb(
normalize_v3(r_edgevec);
}
static void mesh_verts_parallel_range_calc_normals_accum_cb(
static void bm_partial_verts_parallel_range_calc_normal_cb(
void *userdata, const int iter, const TaskParallelTLS *__restrict UNUSED(tls))
{
BMFace *f = ((BMFace **)((void **)userdata)[0])[iter];
const float(*edgevec)[3] = (float(*)[3])((void **)userdata)[1];
mesh_verts_calc_normals_accum(f, f->no, edgevec, NULL);
}
static void mesh_verts_parallel_range_calc_normals_normalize_cb(
void *userdata, const int iter, const TaskParallelTLS *__restrict UNUSED(tls))
{
BMVert *v = ((BMVert **)userdata)[iter];
if (UNLIKELY(normalize_v3(v->no) == 0.0f)) {
normalize_v3_v3(v->no, v->co);
}
BMVert *v = ((BMVert **)((void **)userdata)[0])[iter];
const float(*edgevec)[3] = (const float(*)[3])((void **)userdata)[1];
bm_vert_calc_normals_impl(edgevec, v);
}
/**
@ -316,21 +310,15 @@ void BM_mesh_normals_update_with_partial(BMesh *bm, const BMPartialUpdate *bmpin
const int verts_len = bmpinfo->verts_len;
const int edges_len = bmpinfo->edges_len;
const int faces_len = bmpinfo->faces_len;
const int faces_len_normal_calc_accumulate = bmpinfo->faces_len_normal_calc_accumulate;
float(*edgevec)[3] = MEM_mallocN(sizeof(*edgevec) * edges_len, __func__);
for (int i = 0; i < verts_len; i++) {
zero_v3(verts[i]->no);
}
TaskParallelSettings settings;
BLI_parallel_range_settings_defaults(&settings);
{
/* Faces. */
BLI_task_parallel_range(
0, faces_len, faces, mesh_faces_parallel_range_calc_normals_cb, &settings);
}
/* Faces. */
BLI_task_parallel_range(
0, faces_len, faces, bm_partial_faces_parallel_range_calc_normals_cb, &settings);
/* Temporarily override the edge indices,
* storing the correct indices in the case they're not dirty.
@ -366,19 +354,12 @@ void BM_mesh_normals_update_with_partial(BMesh *bm, const BMPartialUpdate *bmpin
* normals. */
void *data[2] = {edges, edgevec};
BLI_task_parallel_range(
0, edges_len, data, mesh_edges_parallel_range_calc_vectors_cb, &settings);
0, edges_len, data, bm_partial_edges_parallel_range_calc_vectors_cb, &settings);
/* Add weighted face normals to vertices. */
data[0] = faces;
BLI_task_parallel_range(0,
faces_len_normal_calc_accumulate,
data,
mesh_verts_parallel_range_calc_normals_accum_cb,
&settings);
/* Normalize the accumulated vertex normals. */
/* Calculate vertex normals. */
data[0] = verts;
BLI_task_parallel_range(
0, verts_len, verts, mesh_verts_parallel_range_calc_normals_normalize_cb, &settings);
0, verts_len, data, bm_partial_verts_parallel_range_calc_normal_cb, &settings);
}
if (edge_index_value != NULL) {

View File

@ -182,8 +182,6 @@ BMPartialUpdate *BM_mesh_partial_create_from_verts(BMesh *bm,
}
}
const int faces_len_direct = bmpinfo->faces_len;
if (params->do_normals) {
/* - Extend to all faces vertices:
* Any changes to the faces normal needs to update all surrounding vertices.
@ -219,31 +217,13 @@ BMPartialUpdate *BM_mesh_partial_create_from_verts(BMesh *bm,
BMEdge *e_iter = e_first;
do {
if (e_iter->l) {
if (!partial_elem_edge_ensure(bmpinfo, edges_tag, e_iter)) {
continue;
}
/* These faces need to be taken into account when weighting vertex normals
* but aren't needed for tessellation nor do their normals need to be recalculated.
* These faces end up between `faces_len` and `faces_len_normal_calc_accumulate`
* in the faces array. */
BMLoop *l_first_radial = e_iter->l;
BMLoop *l_iter_radial = l_first_radial;
/* Loop over radial loops. */
do {
if (l_iter_radial->v == v) {
partial_elem_face_ensure(bmpinfo, faces_tag, l_iter_radial->f);
}
} while ((l_iter_radial = l_iter_radial->radial_next) != l_first_radial);
partial_elem_edge_ensure(bmpinfo, edges_tag, e_iter);
}
} while ((e_iter = BM_DISK_EDGE_NEXT(e_iter, v)) != e_first);
} while ((l_iter = l_iter->next) != l_first);
}
}
bmpinfo->faces_len_normal_calc_accumulate = bmpinfo->faces_len;
bmpinfo->faces_len = faces_len_direct;
if (verts_tag) {
MEM_freeN(verts_tag);
}

View File

@ -49,11 +49,6 @@ typedef struct BMPartialUpdate {
int verts_len, verts_len_alloc;
int edges_len, edges_len_alloc;
int faces_len, faces_len_alloc;
/**
* Faces at the end of `faces` that don't need to have the normals recalculated
* but must be included when waiting the vertex normals.
*/
int faces_len_normal_calc_accumulate;
/** Store the parameters used in creation so invalid use can be asserted. */
BMPartialUpdate_Params params;