esa / polyhedral-gravity-model

Implementation of a polyhedral gravity model in C++17 with a Python Binding
https://esa.github.io/polyhedral-gravity-model/
GNU General Public License v3.0
21 stars 9 forks source link

Sign of the tensor #29

Closed JoachimSchwabe closed 7 months ago

JoachimSchwabe commented 8 months ago

This is also related to issue #28 I did a small test by comparing the tensor quantities with the numerical differences of the "gravity" components (derivative of the potential) for delta = 1 m, as follows ...

d2V/dx2(x,y,z) vs. [ dV/dx(x+1,y,z) - dV/dx(x,y,z) ] d2V/dy2(x,y,z) vs. [ dV/dy(x,y+1,z) - dV/dy(x,y,z) ] d2V/dz2(x,y,z) vs. [ dV/dz(x,y,z+1) - dV/dz(x,y,z) ]

... and found out that the signs of the tensor do not match the delta(Vx/Vy/Vz).

On the other hand, the same also holds for the comparison with the first-order derivatives and the potential ...

dV/dz(x,y,z) vs. [ V(x,y,z+1) - V(x,y,z) ]

... the sign of e.g. Vz does not match the potential difference for z+1.

So it appears that the Vx/Vy/Vz signs are indeed flipped compared to V and the tensor. But again, as mentioned in #28, the potential of a polyhedral body should be positive for a point outside the mesh hull.

Two ideas for reasons:

A) Does the definition of the x,y,z axes matter? I used UTM (a right-handed coordinate system) both for test body and computation points with x ... Easting, y ... Northing and z ... Up? But actually I think that should not make a difference...

B) Maybe the mesh of my test body is "flipped" regarding inside/outside. But then still, the first-order derivatives do not match the potential and the tensor.

UPDATE 2024-03-11 (copied from #28):

I checked once more with the example cube from the description. My computation points were located at (0, 0, 3) and (0, 0, 3.001) above the test cube, and the signs were as follows: V : positive (decreasing upwards) -> Vz > 0 acc[2] : positive acceleration (decreasing upwards) --> negative of Vz ! tensor[2] : positive --> correct sign of Vzz = d/dz(Vz), but flipped w.r.t. to output of acc[2]

From this I conclude that the mesh of my own test body is indeed "inverted", so that all quantities are additionally flipped.

schuhmaj commented 7 months ago

Many Thanks for reporting this issue. Please find my answer in #28 regarding the inverted acceleration sign.

Regarding the "inverted mesh" causing "double inversion" Have you already checked that the plane unit normals point outwards of the polyhedron? In case the plane unit normals point inwards, the signs are inverted.

You can do this via polyhedral_gravity.utility.check_mesh(vertices, faces)

JoachimSchwabe commented 7 months ago

Also many thanks for the hint with check_mesh!

Indeed, the check was "false" for my test body. Maybe because it is hollow on the bottom side (see attached screenshot). Dome

On the other hand, the check for the test cube was also "false". I realized that I had one flipped sign in the very last node that caused it.

That being said, I have to revoke my initial posting. With the correct test cube nodes, I now get the following signs, which are in fact according to the geodetic convention:

potential: positive (outside test cube, decreasing facing away from the cube) acceleration: positive (for x, y and z > 0 and outside test cube) => acc_z = -Tz negative etc. , magnitude decreasing outwards tensor: positive, i.e., Tzz = d/dz(Tz) = -d/dz(acc_z) => positive since magnitude of negative Tz becomes smaller

If e.g. x < 0 and outside the cube, the output acc_x then becomes negative, since Tx = d/dx(pot) > 0, etc.

Sorry for the confusion!