Skip to content

Commit

Permalink
review feedback
Browse files Browse the repository at this point in the history
  • Loading branch information
VeckoTheGecko committed Feb 3, 2025
1 parent 3d5513e commit ecf9d80
Show file tree
Hide file tree
Showing 3 changed files with 12 additions and 12 deletions.
10 changes: 5 additions & 5 deletions parcels/application_kernels/advection.py
Original file line number Diff line number Diff line change
Expand Up @@ -232,14 +232,14 @@ def AdvectionAnalytical(particle, fieldset, time): # pragma: no cover
else:
dz = 1.0

c1 = i_u.geodesic_distance(py[0], py[1], px[0], px[1], grid.mesh, np.dot(i_u.phi2D_lin(0.0, xsi), py))
c2 = i_u.geodesic_distance(py[1], py[2], px[1], px[2], grid.mesh, np.dot(i_u.phi2D_lin(eta, 1.0), py))
c3 = i_u.geodesic_distance(py[2], py[3], px[2], px[3], grid.mesh, np.dot(i_u.phi2D_lin(1.0, xsi), py))
c4 = i_u.geodesic_distance(py[3], py[0], px[3], px[0], grid.mesh, np.dot(i_u.phi2D_lin(eta, 0.0), py))
c1 = i_u._geodetic_distance(py[0], py[1], px[0], px[1], grid.mesh, np.dot(i_u.phi2D_lin(0.0, xsi), py))
c2 = i_u._geodetic_distance(py[1], py[2], px[1], px[2], grid.mesh, np.dot(i_u.phi2D_lin(eta, 1.0), py))
c3 = i_u._geodetic_distance(py[2], py[3], px[2], px[3], grid.mesh, np.dot(i_u.phi2D_lin(1.0, xsi), py))
c4 = i_u._geodetic_distance(py[3], py[0], px[3], px[0], grid.mesh, np.dot(i_u.phi2D_lin(eta, 0.0), py))
rad = np.pi / 180.0
deg2m = 1852 * 60.0
meshJac = (deg2m * deg2m * math.cos(rad * particle.lat)) if grid.mesh == "spherical" else 1
dxdy = i_u.compute_jacobian_determinant(py, px, eta, xsi) * meshJac
dxdy = i_u._compute_jacobian_determinant(py, px, eta, xsi) * meshJac

if withW:
U0 = direction * fieldset.U.data[ti, zi + 1, yi + 1, xi] * c4 * dz
Expand Down
10 changes: 5 additions & 5 deletions parcels/field.py
Original file line number Diff line number Diff line change
Expand Up @@ -2003,10 +2003,10 @@ def spatial_c_grid_interpolation2D(self, ti, z, y, x, time, particle=None, apply
px[1:] = np.where(-px[1:] + px[0] > 180, px[1:] + 360, px[1:])
xx = (1 - xsi) * (1 - eta) * px[0] + xsi * (1 - eta) * px[1] + xsi * eta * px[2] + (1 - xsi) * eta * px[3]
assert abs(xx - x) < 1e-4
c1 = i_u.geodesic_distance(py[0], py[1], px[0], px[1], grid.mesh, np.dot(i_u.phi2D_lin(0.0, xsi), py))
c2 = i_u.geodesic_distance(py[1], py[2], px[1], px[2], grid.mesh, np.dot(i_u.phi2D_lin(eta, 1.0), py))
c3 = i_u.geodesic_distance(py[2], py[3], px[2], px[3], grid.mesh, np.dot(i_u.phi2D_lin(1.0, xsi), py))
c4 = i_u.geodesic_distance(py[3], py[0], px[3], px[0], grid.mesh, np.dot(i_u.phi2D_lin(eta, 0.0), py))
c1 = i_u._geodetic_distance(py[0], py[1], px[0], px[1], grid.mesh, np.dot(i_u.phi2D_lin(0.0, xsi), py))
c2 = i_u._geodetic_distance(py[1], py[2], px[1], px[2], grid.mesh, np.dot(i_u.phi2D_lin(eta, 1.0), py))
c3 = i_u._geodetic_distance(py[2], py[3], px[2], px[3], grid.mesh, np.dot(i_u.phi2D_lin(1.0, xsi), py))
c4 = i_u._geodetic_distance(py[3], py[0], px[3], px[0], grid.mesh, np.dot(i_u.phi2D_lin(eta, 0.0), py))
if grid.zdim == 1:
if self.gridindexingtype == "nemo":
U0 = self.U.data[ti, yi + 1, xi] * c4
Expand Down Expand Up @@ -2038,7 +2038,7 @@ def spatial_c_grid_interpolation2D(self, ti, z, y, x, time, particle=None, apply
else:
meshJac = deg2m if grid.mesh == "spherical" else 1

jac = i_u.compute_jacobian_determinant(py, px, eta, xsi) * meshJac
jac = i_u._compute_jacobian_determinant(py, px, eta, xsi) * meshJac

u = (
(-(1 - eta) * U - (1 - xsi) * V) * px[0]
Expand Down
4 changes: 2 additions & 2 deletions parcels/tools/interpolation_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -175,7 +175,7 @@ def interpolate(phi: Callable[[float], list[float]], f: list[float], xsi: float)
return np.dot(phi(xsi), f)


def geodesic_distance(lat1: float, lat2: float, lon1: float, lon2: float, mesh: Mesh, lat: float) -> float:
def _geodetic_distance(lat1: float, lat2: float, lon1: float, lon2: float, mesh: Mesh, lat: float) -> float:
if mesh == "spherical":
rad = np.pi / 180.0
deg2m = 1852 * 60.0
Expand All @@ -184,7 +184,7 @@ def geodesic_distance(lat1: float, lat2: float, lon1: float, lon2: float, mesh:
return np.sqrt((lon2 - lon1) ** 2 + (lat2 - lat1) ** 2)


def compute_jacobian_determinant(py: np.ndarray, px: np.ndarray, eta: float, xsi: float) -> float:
def _compute_jacobian_determinant(py: np.ndarray, px: np.ndarray, eta: float, xsi: float) -> float:
dphidxsi = [eta - 1, 1 - eta, eta, -eta]
dphideta = [xsi - 1, -xsi, xsi, 1 - xsi]

Expand Down

0 comments on commit ecf9d80

Please sign in to comment.