Update perturbation.py for NumPy 2.0 compatibility (#656)

### Description
This changes `perturbation.py` to work with NumPy 2.0.

### Related issues
Fixes #653.

### Checklist
_Before this pull request can be reviewed, all of these tasks should be
completed. Denote completed tasks with an `x` inside the square brackets
`[ ]` in the Markdown source below:_
- [x] I have added a description (see above).
- [x] I have added a link to any related issues see (see above).
- [x] I have read the [Contributing
Guide](https://github.com/quokka-astro/quokka/blob/development/CONTRIBUTING.md).
- [ ] I have added tests for any new physics that this PR adds to the
code.
- [ ] I have tested this PR on my local computer and all tests pass.
- [x] I have manually triggered the GPU tests with the magic comment
`/azp run`.
- [x] I have requested a reviewer for this PR.

src/perturbation.py

@@ -69,7 +69,7 @@ def init_perturbations(n, kmin, kmax, dtype):
 
 
 def normalize(fx, fy, fz):
-    norm = np.sqrt(np.sum(fx**2 + fy**2 + fz**2)/np.product(n))
+    norm = np.sqrt(np.sum(fx**2 + fy**2 + fz**2)/np.prod(n))
     fx = fx/norm
     fy = fy/norm
     fz = fz/norm
@@ -162,7 +162,7 @@ def get_erot_ke_ratio(pertx, perty, pertz, rad=-1.0):
         r2[idx0] = 0.0
     erot_ke_ratio = (np.sum(y*pertz-z*perty)**2 +
                      np.sum(z*pertx-x*pertz)**2 +
-                     np.sum(x*perty-y*pertx)**2)/(np.sum(r2)*np.product(n)) 
+                     np.sum(x*perty-y*pertx)**2)/(np.sum(r2)*np.prod(n)) 
     return erot_ke_ratio