I know this is an old thread, but it seems like the solution has not been found.
Today I encountered the same error while working with huge point cloud / mesh data sets. Actually this bug is hardly noticeable, because it’s rather random. I found out that rounding the coordinates of points helps and makes the results reliable.
I noticed that if points’ coordinates look like
new Point3d(1.237928032875061, 0.762071967124939, 1.237928032875061)
then the Mesh.IsPointInside() will give unpredictable results (I tested on simple mesh boxes). Rounding to numbers with six or less fractional digits like
new Point3d(1.237928, 0.762071, 1.237928)
makes the results reliable. Try this on your test definitions. It works for me and I think it’s fine since the guaranteed double precision is 6 fractional digits.
Btw writing own mesh inclusion (i.e. shooting semi-infinite rays from point and checking if the number of intersections is even or not) works fine, but is slower because Mesh.IsPointInside() uses fast unmanaged C++ libraries.
I guess the bug might have its roots in converting numerical types between managed and unmanaged code, but McNeel developers could probably elaborate a bit more about that.
I attach a test definition with internalized / hard-coded values to recreate the bug.
20190601_MeshPointInclusionBug.gh (14.5 KB)