# SPH to Eulerian grid question

1. Jun 9, 2014

### TimeFall

Hello! I'm currently working on an algorithm that converts SPH (smoothed particle hydrodynamics) data to an Eulerian grid. Basically, the problem is this: each sph particle has an associated smoothing length, which is effectively the radius of a sphere centered at the sph particle's position. The properties of the particle as "smeared" out by a so-called smoothing kernel W(r, h), where r is the distance from the particle and h is the smoothing length. W is a piecewise polynomial. The idea is to get the points on the cell where the sphere overlaps it, and then integrate the smoothing kernel over this overlap volume. My problem is that I do not know of an efficient, effective way of getting these intersection points so that I can do the kernel integration, and was wondering if anyone out there did?

In a nutshell, I have a sphere that overlaps with a cube (not necessarily fully), and I would like to know the points of intersection in the x, y, and z directions. Or, if anyone knows of a better way to approach this problem, I am open to suggestions. Any help would be greatly appreciated! Thanks!

2. Jun 10, 2014

### .Scott

The intersection problem isn't that bad, but there are a lot of combinations to consider.
Basically, you start with the intersection of the sphere with each of the 6 planes - but I don't you need to do that.

You haven't said this explicitly, but I believe that your cube is already lined up with the x, y, and z axis. If not, rotate to a coordinate system where it is. Assuming this is required, it will involve changing the sph center coordinates.

I'm assuming that you only have one cell - I'm a little foggy on exactly what you're doing in that regard.

My first thought was that you wanted to scan through x,y,z coordinates that are shared by the cube and the sph. I'm now thinking that's not exactly what you want, but this might be useful anyway. When you respond, I can be more exact:
So start as if you were going to scan the entire cure - with three nested loops, one each for x, y, and z. Then limit the min x and max x value to those that are shared by the cube and the sphere. If the center of the sphere is Xs, Ys, Zs; the radius R; and the cube limits Xcmin, Xcmax, Ycmin, Ycmax, Zcmin, Zcmax, then you'll be limiting the range of "X" in the outer for loop to (Xcmin to Xcmax) and (Xs +/- R).
Then inside the X loop, compute the range of Y values - limited by (Ycmin to Ycmax) and (Ys +/- sqrt(R*R-(Xs-X)^2)). Then inside the Y loop, compute the range of Z values - limited by (Zcmin to Zcmax) and (Zs +/- sqrt(R*R-(Xs-X)^2-(Ys-Y)^2)).

Do the optimizations explicitly - as opposed to hoping the compiler will do them for you. For example, precompute R*R.

3. Jun 10, 2014

### .Scott

When you say, "the points on the cell where the sphere overlaps it", are these points along a grid in the interior of the cell and sphere, or the interiors without any grid, or points along the surface of the cell?

I don't know what the integration is that you're performing. The details may reveal significant possibilities for optimization.

4. Jun 10, 2014

### TimeFall

Hey .Scott, thanks for replying! There is more than one cell. I probably should have been more explicit in my description of the problem. What I'm doing is dividing up the simulation volume (which has periodic boundary conditions, so if one sph smoothing sphere goes off the edge of the simulation volume, it has to wraparound to the opposite side) into a regular grid of a user-defined number of cells in each dimension. Yes, the grid is lined up with the x,y, and z axis. The function I need to integrate is:
$$W(r, h) = \frac{8}{\pi h^{3}} \begin{cases} 1 - 6(\frac{r}{h})^{2} + 6(\frac{r}{h})^{3}, &0 \leq \frac{r}{h} \leq \frac{1}{2}\\ 2(1 - \frac{r}{h})^{3}, & \frac{1}{2} < \frac{r}{h} \leq 1 \\ 0, & \frac{r}{h} > 1\end{cases}$$
To answer your last question, the region of overlap is defined by points along the surface of the cell (basically, where the sphere "passes through" the surface of the cell).
Another thing that I tried is a sort of quasi Monte Carlo integration that found a normalization weight for each sph particle and then filled each cell with a random field. I then checked to see how many of the random field particles were inside the smoothing sphere and then applied the normalization constant. This worked decently well, but it takes, effectively, forever to run on the larger simulation output files. However, the reason I tried doing it this way was, as you said, there are a ton of different cases you have to worry about when trying to find the integration limits.