`lax.rem` has wrong gradient for `y` when `x<0`
`lax.rem` has wrong gradient when `x<0`. For example,

```py
import jax

def fn(y):
    x = jax.numpy.array(-19., jax.numpy.float64)
    return jax.lax.rem(x, y)

y = jax.numpy.array(3., jax.numpy.float64)
res = fn(y)
print(res) # -1.0

a = jax.jacfwd(fn)(y)
print(a) # 7.0
```

Actually, the gradient should be 6 instead of 7. Because `-19 = (-6)*3 + (-1)`. When the `y` is `3.1`, the results will be `-0.4`. Thus, the gradient should be 6.

By the way, another two similar APIs `numpy.remainder` and `numpy.mod` also have this issue.

```py
def fn(y):
    x = jax.numpy.array(-19., jax.numpy.float64)
    return jax.numpy.remainder(x, y)

y = jax.numpy.array(3., jax.numpy.float64)
res = fn(y)
print(res) # 2.0

a = jax.jacfwd(fn)(y)
print(a) # 8.0
```

In these two APIs, `-19 = (-7)*3 + 2`, so their gradient should be 7 instead of 8.
diff --git a/jax/_src/lax/lax.py b/jax/_src/lax/lax.py
index aeb80847c..3980d79bc 100644
--- a/jax/_src/lax/lax.py
+++ b/jax/_src/lax/lax.py
@@ -2161,7 +2161,7 @@ rem_p = standard_naryop([_int | _float, _int | _float], 'rem')
 ad.defjvp(
     rem_p,
     lambda g, x, y: _maybe_broadcast(broadcast_shapes(np.shape(x), np.shape(y)), g),
-    lambda g, x, y: mul(neg(g), floor(div(x, y))))
+    lambda g, x, y: mul(neg(g), mul(sign(div(x, y)), floor(abs(div(x, y))))))
 mlir.register_lowering(rem_p, partial(_nary_lower_mhlo, mhlo.RemOp))
 
 def _minmax_complex_lowering(x, y, *, lax_cmp_pick_x):
diff --git a/tests/lax_autodiff_test.py b/tests/lax_autodiff_test.py
index 155fab36e..41346ae45 100644
--- a/tests/lax_autodiff_test.py
+++ b/tests/lax_autodiff_test.py
@@ -1060,19 +1060,24 @@ class LaxAutodiffTest(jtu.JaxTestCase):
 
   # TODO(mattjj): make this a more systematic test
   def testRemainder(self):
+    def gen_x(rng, size):
+      return rng.uniform(-9, 9, size=size)
+
+    def gen_y(rng, size):
+      # avoid values near zero because gradients diverge
+      return rng.uniform(0.1, 5, size=size) * rng.choice([-1, 1], size=size)
+
     rng = self.rng()
-    x = rng.uniform(-0.9, 9, size=(3, 4))
-    y = rng.uniform(0.7, 1.9, size=(3, 1))
+    x = gen_x(rng, (5, 8))
+    y = gen_y(rng, (1, 8))
     assert not set(np.unique(x)) & set(np.unique(y))
-    # TODO(jakevdp) try to make these tolerances tighter.
-    tol = 1e-1
-    check_grads(lax.rem, (x, y), 2, ["fwd", "rev"], tol, tol)
+    check_grads(lax.rem, (x, y), 2, ["fwd", "rev"])
 
     rng = self.rng()
-    x = rng.uniform(-0.9, 9, size=(1, 4))
-    y = rng.uniform(0.7, 1.9, size=(3, 4))
+    x = gen_x(rng, (1, 8))
+    y = gen_y(rng, (5, 8))
     assert not set(np.unique(x)) & set(np.unique(y))
-    check_grads(lax.rem, (x, y), 2, ["fwd", "rev"], tol, tol)
+    check_grads(lax.rem, (x, y), 2, ["fwd", "rev"])
 
   def testHigherOrderGradientOfReciprocal(self):
     # Regression test for https://github.com/google/jax/issues/3136
['jax/_src/lax/lax.py', 'tests/lax_autodiff_test.py']