LinearRegression with zero sample_weights is not the same as excluding those rows

### Describe the bug

Excluding rows having `sample_weight == 0` in `LinearRegression` does not give the same results.

### Steps/Code to Reproduce

```python
from collections import Counter
import numpy as np
from sklearn.linear_model import LinearRegression

results = []
for i in range(100):
    rng = np.random.RandomState(i)
    n_samples, n_features = 10, 5
    X = rng.rand(n_samples, n_features)
    y = rng.rand(n_samples)
    reg = LinearRegression()
    sample_weight = rng.uniform(low=0.01, high=2, size=X.shape[0])
    sample_weight_0 = sample_weight.copy()
    sample_weight_0[-5:] = 0
    y[-5:] *= 1000  # to make excluding those samples important
    reg.fit(X, y, sample_weight=sample_weight_0)
    coef_0, intercept_0 = reg.coef_.copy(), reg.intercept_

    reg.fit(X[:-5], y[:-5], sample_weight=sample_weight[:-5])
    print(f"{coef_0=}")
    print(f"{reg.coef_=}")
    results.append(np.allclose(reg.coef_, coef_0, rtol=1e-6))
print(Counter(results))
```

### Expected Results

Always `True`.

### Actual Results

```python
Counter({True: 79, False: 21})  # it fails 20% of the time
```

The print statement gives:
```
coef_0 =    array([ 1.43516166, -1.78826443,  0.15365526,  1.82233166, -1.6       ])
reg.coef_ = array([ 2.24022351, -1.04917851,  0.45341088,  0.80315086, -0.25726798])
```


### Versions

```shell
System:
    python: 3.9.15 (main, Nov 15 2022, 05:24:15)  [Clang 14.0.0 (clang-1400.0.29.202)]
   machine: macOS-13.3.1-x86_64-i386-64bit

Python dependencies:
      sklearn: 1.2.2
        numpy: 1.23.5
        scipy: 1.10.1
       Cython: 0.29.33

       user_api: blas
   internal_api: openblas
         prefix: libopenblas
```
```

COMMENT:
My take at this issue: It's the minimum norm issue with `LinearRegression`, see #22947 and #22855. The effective dimensions/shape of `X` are `(5, 6)` including the intercept, so an underdetermined system.
```python
import numpy as np
from scipy.linalg import lstsq
import scipy.sparse as sparse
from sklearn.linear_model import LinearRegression


rng = np.random.RandomState(2)  # Seed dependent, np.allclose should be False.
n_samples, n_features = 10, 5
X = rng.rand(n_samples, n_features)
y = rng.rand(n_samples)
reg = LinearRegression()
sample_weight = rng.uniform(low=0.01, high=2, size=X.shape[0])
sample_weight_0 = sample_weight.copy()
sample_weight_0[-5:] = 0
y[-5:] *= 1000  # to make excluding those samples important

reg.fit(X, y, sample_weight=sample_weight_0)
coef_0 = np.r_[reg.coef_.copy(), reg.intercept_]

reg.fit(X[:-5], y[:-5], sample_weight=sample_weight[:-5])
coef_1 = np.r_[reg.coef_.copy(), reg.intercept_]

# For fun, we also use the lsqr solve by converting so sparse X
reg.fit(sparse.csr_array(X), y, sample_weight=sample_weight_0)
coef_2 = np.r_[reg.coef_.copy(), reg.intercept_] # same as coef_1

np.allclose(coef_0, coef_1, rtol=1e-6), np.allclose(coef_1, coef_2, rtol=1e-6)
```
Should be `False, True`.
```python
# We solve with plain lstsq and intercept term.
X_with_intercept = np.c_[X, np.ones(10)]
coef_lstsq = lstsq(np.sqrt(sample_weight_0)[:, None] * X_with_intercept, np.sqrt(sample_weight_0) * y)[0]

# Norm of residues
[
    np.linalg.norm(sample_weight_0 * (X_with_intercept @ coef_0 - y)),
    np.linalg.norm(sample_weight_0 * (X_with_intercept @ coef_1 - y)),
    np.linalg.norm(sample_weight_0 * (X_with_intercept @ coef_lstsq - y)),
]
```
All around 1e-16. So all are valid solutions.
```python
# Norm of coefficients
[
    np.linalg.norm(coef_0),
    np.linalg.norm(coef_1),
    np.linalg.norm(coef_lstsq),
]
```
`[3.3496024918115093, 2.7998306762269527, 2.7862672453159094]`
So the manual `lstsq` has the smallest L2-norm.

Labels: Bug, module:linear_model