"anisotropic_points_to_implicit3" problem

[Immocat]

Hi,
It's me (#109) again 😸 . I find a problem with anisotropic point to implicit3.

What problem:

Here is a screenshot:

It seems not correct for a dambreak at frame 0.

How:

I just run a sph dambreak demo similar to sph_sim example 3, and here is the .pos file at frame 0 that I got.
frame_000000.pos

And then, I tried to reconstruct the surface using the anisotropic kernel.

./particles2obj -i frame_000000.pos -o frame_000000.obj -r 300,200,150 -g 0.01,0.01,0.01 -k 0.036

(kernal radius:0.036, resolution is 300 200 150, and the gridspacing is 0.01, 0.01, 0.01. sph example 3's domain is 3,2,1.5, so I think it is reasonable)

My goal is to get a surface as tight as possible to render. My sph particle radius is 0.02(and the kernel radius is 1.8 * radius) as default, so I chose 0.036 as reconstruction kernel radius. I'm not sure if it is appropriate. Finally, I get the .obj as the screenshot shows.

My opinion:

I looked into the source code and the original paper

I find that the code is a little bit different from what the paper describes. The main difference is in equation 15. Instead of using the magic number k_s = 1400, the code uses std::pow(relV, 1.0 / 3.0) at line 152 in "anisotropic_points_to_implicit3.cpp"

I do not understand why the code uses "pow 1/3", but I tried to use the magic number 1400 in the same way as the paper, and the problem is gone! Here is the screenshot using k_s = 1400.

(Sorry for the light spot in the middle. It is the UI of Blender)

I've not tested the 2D version of the anisotropic kernel, maybe it suffers from the same problem.

That's all I found... I hope I describe the problem clearly. I will maintain the picture and the .pos file till this issue is closed.

[doyubkim]

Thanks for finding this, @Immocat! Let me take a look at the problem. I’m not sure if the magic number is the real root cause, but still very interesting find.

[dokipen3d]

Hi all.

I've been implementing this paper at work using vdb pointDataGrids into a volume grid. I have been a bit confused by that magic number too. I've also encountered banding at the voxel splat stage. I think I'm close to figuring out what the issue is and it's related to the magic scaling.

So here's the problem (I think). I have a particle grid (which stores the particles in 'index' space ie the integer coordinates of the grid which is the way vdb accesses it's coordinates). I calculate the svd/matrix using those coordinates. The svd is fine as I can load the particles back into houdini and instance spheres onto them using the matrix and they are all oriented and squashed correctly.
Then I go to the voxel splat stage which is where I've been having the banding issues. The thing causing the issue is that I am writing to a voxel grid of arbitrary resolution. I do the neighbour search and splat the particles using a kernel (i've tried many) and mult the r vector by the matrix. (one tricky thing to keep in my mind is that I am doing the voxel splat in the particles index space to prevent me having to do a transform into a different 'space').

The thing I think is the problem is that I am using the matrix which is of a certain 'scale' but the voxel splat has a different 'support' radius. So the matrix scaling and the kernel radius are of different scales. So if the anisotropic matrix has eigenvalues of something like (2.0, 0.5, 2.0), and the voxel splat is only 1.0, then the r vectors will be stretched out of the kernel support and messes up the splatted value. I'm going to try and scale down the matrix by a factor that is the difference between the two kernels in the particle and voxel stages.

Also want to try a metaball approach with sdf blending still using the matrices.

I'll follow this thread with interest!

[doyubkim]

Thank you for chiming in, @dokipen3d! Yes, I also experienced similar issues with the kernel support not being fully captured by the grid points. Probably this is related to that issue based on the first screenshot @Immocat attached (SDF being not smooth which indicates kernel cut-off).

[Immocat]

Thanks for all your replies.
Although using the ks = 1400 solve the "hole" problem, I found a new problem:

I'm not sure if it is the banding problem @dokipen3d mentioned. I used the same frame particles to generate these two meshes, while the left one's grid spacing is 0.01 and the right is 0.02.

So for my current project, I choose to use zhu_bridson's method 😸. I personally prefer zhu_bridson's result for its sharp and thin surfaces, except it is still suffered from bumpy problem to some extent and the shortcoming described in their paper.

Comparing to zhu_bridson's method, anisotropic method is much slower because of the svd. So I also recommend @doyubkim to adopt the Implicit-shifted Symmetric QR Singular Value Decomposition of 3 × 3(2x2) Matrices.

Not only this method is widely used and critical in Disney's Material point method simulation(maybe in your future book?:smile_cat:), but also it is also fast and not that general as the code from "numerical recipe", which is currently used in the code base. QRSVD's result has some convention(mentioned in part 1 in its paper) which perfectly matches the requirement in the anisotropic paper(2x2,3x3 matrices, result is rotation matrix det = 1, sorted eigen values).

I will try to implement material point method and adopt the QRSVD as my future project and do some experiment with this anisotropic kernel but not in the near future.

[doyubkim]

That is great suggestion, @Immocat. Yes, I was also not very proud of the performance of the anisotropic implementation. The SDV is not very efficient, but also using the uniform grid also make things slow. Let me open a separate ticket to cover that issue.

[dokipen3d]

I've been using eigen's svd for the moment. At one point I had a slow access pattern for the vdb neighbour search and the svd wasn't an issue. Now I've got a faster search, the svd could probably now be looked at. There is a guy here at work who implemented the fast 3x3 svd so I'm going to look at that as well.

The banding does look familiar. It was definely tied to the voxel grid res. It varied with the angle too.

I wish I could share screenshots but here at Weta, we can't really get stuff out without getting it all approved.

[doyubkim]

Just pushed a branch with a potential fix. @Immocat, if you can checkout that branch and see if it fixes your problem, that would be great. It seems to work with most of my test cases including your attached file.

The actual bug was caused by wrong clamping of the singular values after SVD step.

Not related to the bug, but if I make additional comment on my magic parameter, I'm actually implementing the algorithm slightly different from the original paper in order to remove the original magic parameter, 1400. The scaling factor I'm using (pow(relV, 1/3) part) is trying to normalize the kernel volume so that either perfect spherical kernel or elongated kernel have same volume.

Also, about the aliasing, I'm using fixed 2x radius than the kernel size which can be improved I think.

[dokipen3d]

Why would absmax fix the holes? Can svd give negative scales?

[dokipen3d]

Yeah if we can fix this aliasing, I will be very happy!

[doyubkim]

@dokipen3d great question. I think my PR and comment were misleading. Apparently, the version of implementation I'm using for SVD does generate negative output, which is totally fine unless I carry the sign all the way to the G matrix. So the fix is cancelling out the sign part and use the abs scale to construct G. But as I type this comment here, I think I can make my code and comment more clear.