diff --git a/doc/python_api_reference_vDev.md b/doc/python_api_reference_vDev.md
index ee5b43b49..9905cd76c 100644
--- a/doc/python_api_reference_vDev.md
+++ b/doc/python_api_reference_vDev.md
@@ -173,6 +173,7 @@ API references for stable versions are kept on the [stim github wiki](https://gi
 - [`stim.FlipSimulator`](#stim.FlipSimulator)
     - [`stim.FlipSimulator.__init__`](#stim.FlipSimulator.__init__)
     - [`stim.FlipSimulator.batch_size`](#stim.FlipSimulator.batch_size)
+    - [`stim.FlipSimulator.broadcast_pauli_errors`](#stim.FlipSimulator.broadcast_pauli_errors)
     - [`stim.FlipSimulator.do`](#stim.FlipSimulator.do)
     - [`stim.FlipSimulator.get_detector_flips`](#stim.FlipSimulator.get_detector_flips)
     - [`stim.FlipSimulator.get_measurement_flips`](#stim.FlipSimulator.get_measurement_flips)
@@ -6153,6 +6154,61 @@ def batch_size(
     """
 ```
 
+<a name="stim.FlipSimulator.broadcast_pauli_errors"></a>
+```python
+# stim.FlipSimulator.broadcast_pauli_errors
+
+# (in class stim.FlipSimulator)
+def broadcast_pauli_errors(
+    self,
+    *,
+    pauli: Union[str, int],
+    mask: np.ndarray,
+) -> None:
+    """Applies a pauli error to all qubits in all instances, filtered by a mask.
+
+    Args:
+        pauli: The pauli, specified as an integer or string.
+            Uses the convention 0=I, 1=X, 2=Y, 3=Z.
+            Any value from [0, 1, 2, 3, 'X', 'Y', 'Z', 'I', '_'] is allowed.
+        mask: A 2d numpy array specifying where to apply errors. The first axis
+            is qubits, the second axis is simulation instances. The first axis
+            can have a length less than the current number of qubits (or more,
+            which adds qubits to the simulation). The length of the second axis
+            must match the simulator's `batch_size`. The array must satisfy
+
+                mask.dtype == np.bool_
+                len(mask.shape) == 2
+                mask.shape[1] == flip_sim.batch_size
+
+            The error is only applied to qubit q in instance k when
+
+                mask[q, k] == True.
+
+    Examples:
+        >>> import stim
+        >>> import numpy as np
+        >>> sim = stim.FlipSimulator(
+        ...     batch_size=2,
+        ...     num_qubits=3,
+        ...     disable_stabilizer_randomization=True,
+        ... )
+        >>> sim.broadcast_pauli_errors(
+        ...     pauli='X',
+        ...     mask=np.asarray([[True, False],[False, False],[True, True]]),
+        ... )
+        >>> sim.peek_pauli_flips()
+        [stim.PauliString("+X_X"), stim.PauliString("+__X")]
+
+        >>> sim.broadcast_pauli_errors(
+        ...     pauli='Z',
+        ...     mask=np.asarray([[False, True],[False, False],[True, True]]),
+        ... )
+        >>> sim.peek_pauli_flips()
+        [stim.PauliString("+X_Y"), stim.PauliString("+Z_Y")]
+    """
+```
+
 <a name="stim.FlipSimulator.do"></a>
 ```python
 # stim.FlipSimulator.do
diff --git a/doc/stim.pyi b/doc/stim.pyi
index 697608bdc..6a7ff42f8 100644
--- a/doc/stim.pyi
+++ b/doc/stim.pyi
@@ -4670,6 +4670,54 @@ class FlipSimulator:
             >>> sim.batch_size
             42
         """
+    def broadcast_pauli_errors(
+        self,
+        *,
+        pauli: Union[str, int],
+        mask: np.ndarray,
+    ) -> None:
+        """Applies a pauli error to all qubits in all instances, filtered by a mask.
+
+        Args:
+            pauli: The pauli, specified as an integer or string.
+                Uses the convention 0=I, 1=X, 2=Y, 3=Z.
+                Any value from [0, 1, 2, 3, 'X', 'Y', 'Z', 'I', '_'] is allowed.
+            mask: A 2d numpy array specifying where to apply errors. The first axis
+                is qubits, the second axis is simulation instances. The first axis
+                can have a length less than the current number of qubits (or more,
+                which adds qubits to the simulation). The length of the second axis
+                must match the simulator's `batch_size`. The array must satisfy
+
+                    mask.dtype == np.bool_
+                    len(mask.shape) == 2
+                    mask.shape[1] == flip_sim.batch_size
+
+                The error is only applied to qubit q in instance k when
+
+                    mask[q, k] == True.
+
+        Examples:
+            >>> import stim
+            >>> import numpy as np
+            >>> sim = stim.FlipSimulator(
+            ...     batch_size=2,
+            ...     num_qubits=3,
+            ...     disable_stabilizer_randomization=True,
+            ... )
+            >>> sim.broadcast_pauli_errors(
+            ...     pauli='X',
+            ...     mask=np.asarray([[True, False],[False, False],[True, True]]),
+            ... )
+            >>> sim.peek_pauli_flips()
+            [stim.PauliString("+X_X"), stim.PauliString("+__X")]
+
+            >>> sim.broadcast_pauli_errors(
+            ...     pauli='Z',
+            ...     mask=np.asarray([[False, True],[False, False],[True, True]]),
+            ... )
+            >>> sim.peek_pauli_flips()
+            [stim.PauliString("+X_Y"), stim.PauliString("+Z_Y")]
+        """
     def do(
         self,
         obj: Union[stim.Circuit, stim.CircuitInstruction, stim.CircuitRepeatBlock],
diff --git a/glue/python/src/stim/__init__.pyi b/glue/python/src/stim/__init__.pyi
index 697608bdc..6a7ff42f8 100644
--- a/glue/python/src/stim/__init__.pyi
+++ b/glue/python/src/stim/__init__.pyi
@@ -4670,6 +4670,54 @@ class FlipSimulator:
             >>> sim.batch_size
             42
         """
+    def broadcast_pauli_errors(
+        self,
+        *,
+        pauli: Union[str, int],
+        mask: np.ndarray,
+    ) -> None:
+        """Applies a pauli error to all qubits in all instances, filtered by a mask.
+
+        Args:
+            pauli: The pauli, specified as an integer or string.
+                Uses the convention 0=I, 1=X, 2=Y, 3=Z.
+                Any value from [0, 1, 2, 3, 'X', 'Y', 'Z', 'I', '_'] is allowed.
+            mask: A 2d numpy array specifying where to apply errors. The first axis
+                is qubits, the second axis is simulation instances. The first axis
+                can have a length less than the current number of qubits (or more,
+                which adds qubits to the simulation). The length of the second axis
+                must match the simulator's `batch_size`. The array must satisfy
+
+                    mask.dtype == np.bool_
+                    len(mask.shape) == 2
+                    mask.shape[1] == flip_sim.batch_size
+
+                The error is only applied to qubit q in instance k when
+
+                    mask[q, k] == True.
+
+        Examples:
+            >>> import stim
+            >>> import numpy as np
+            >>> sim = stim.FlipSimulator(
+            ...     batch_size=2,
+            ...     num_qubits=3,
+            ...     disable_stabilizer_randomization=True,
+            ... )
+            >>> sim.broadcast_pauli_errors(
+            ...     pauli='X',
+            ...     mask=np.asarray([[True, False],[False, False],[True, True]]),
+            ... )
+            >>> sim.peek_pauli_flips()
+            [stim.PauliString("+X_X"), stim.PauliString("+__X")]
+
+            >>> sim.broadcast_pauli_errors(
+            ...     pauli='Z',
+            ...     mask=np.asarray([[False, True],[False, False],[True, True]]),
+            ... )
+            >>> sim.peek_pauli_flips()
+            [stim.PauliString("+X_Y"), stim.PauliString("+Z_Y")]
+        """
     def do(
         self,
         obj: Union[stim.Circuit, stim.CircuitInstruction, stim.CircuitRepeatBlock],
diff --git a/src/stim/simulators/frame_simulator.pybind.cc b/src/stim/simulators/frame_simulator.pybind.cc
index 56e6795a8..6d503756f 100644
--- a/src/stim/simulators/frame_simulator.pybind.cc
+++ b/src/stim/simulators/frame_simulator.pybind.cc
@@ -33,6 +33,32 @@ std::optional<size_t> py_index_to_optional_size_t(
     return (size_t)i;
 }
 
+uint8_t pybind11_object_to_pauli_ixyz(const pybind11::object &obj) {
+    if (pybind11::isinstance<pybind11::str>(obj)) {
+        std::string s = pybind11::cast<std::string>(obj);
+        if (s == "X") {
+            return 1;
+        } else if (s == "Y") {
+            return 2;
+        } else if (s == "Z") {
+            return 3;
+        } else if (s == "I" || s == "_") {
+            return 0;
+        }
+    } else if (pybind11::isinstance<pybind11::int_>(obj)) {
+        uint8_t v = 255;
+        try {
+            v = pybind11::cast<uint8_t>(obj);
+        } catch (const pybind11::cast_error &) {
+        }
+        if (v < 4) {
+            return (uint8_t)v;
+        }
+    }
+
+    throw std::invalid_argument("Need pauli in ['I', 'X', 'Y', 'Z', 0, 1, 2, 3, '_'].");
+}
+
 pybind11::class_<FrameSimulator<MAX_BITWORD_WIDTH>> stim_pybind::pybind_frame_simulator(pybind11::module &m) {
     return pybind11::class_<FrameSimulator<MAX_BITWORD_WIDTH>>(
         m,
@@ -374,27 +400,7 @@ void stim_pybind::pybind_frame_simulator_methods(
            const pybind11::object &pauli,
            int64_t qubit_index,
            int64_t instance_index) {
-            uint8_t p = 255;
-            try {
-                p = pybind11::cast<uint8_t>(pauli);
-            } catch (const pybind11::cast_error &) {
-                try {
-                    std::string s = pybind11::cast<std::string>(pauli);
-                    if (s == "X") {
-                        p = 1;
-                    } else if (s == "Y") {
-                        p = 2;
-                    } else if (s == "Z") {
-                        p = 3;
-                    } else if (s == "I" || s == "_") {
-                        p = 0;
-                    }
-                } catch (const pybind11::cast_error &) {
-                }
-            }
-            if (p > 3) {
-                throw std::invalid_argument("Expected pauli in [0, 1, 2, 3, '_', 'I', 'X', 'Y', 'Z']");
-            }
+            uint8_t p = pybind11_object_to_pauli_ixyz(pauli);
             if (instance_index < 0) {
                 instance_index += self.batch_size;
             }
@@ -409,6 +415,7 @@ void stim_pybind::pybind_frame_simulator_methods(
                 stats.num_qubits = qubit_index + 1;
                 self.ensure_safe_to_do_circuit_with_stats(stats);
             }
+
             p ^= p >> 1;
             self.x_table[qubit_index][instance_index] = (p & 1) != 0;
             self.z_table[qubit_index][instance_index] = (p & 2) != 0;
@@ -770,4 +777,93 @@ void stim_pybind::pybind_frame_simulator_methods(
                 [stim.PauliString("+YX__")]
         )DOC")
             .data());
+
+    c.def(
+        "broadcast_pauli_errors",
+        [](FrameSimulator<MAX_BITWORD_WIDTH> &self, const pybind11::object &pauli, const pybind11::object &mask) {
+            uint8_t p = pybind11_object_to_pauli_ixyz(pauli);
+
+            if (!pybind11::isinstance<pybind11::array_t<bool>>(mask)) {
+                throw std::invalid_argument("Need isinstance(mask, np.ndarray) and mask.dtype == np.bool_");
+            }
+            const pybind11::array_t<bool> &arr = pybind11::cast<pybind11::array_t<bool>>(mask);
+
+            if (arr.ndim() != 2) {
+                throw std::invalid_argument(
+                    "Need a 2d mask (first axis is qubits, second axis is simulation instances). Need len(mask.shape) "
+                    "== 2.");
+            }
+
+            pybind11::ssize_t s_mask_num_qubits = arr.shape(0);
+            pybind11::ssize_t s_mask_batch_size = arr.shape(1);
+            if ((uint64_t)s_mask_batch_size != self.batch_size) {
+                throw std::invalid_argument("Need mask.shape[1] == flip_sim.batch_size");
+            }
+            if (s_mask_num_qubits > UINT32_MAX) {
+                throw std::invalid_argument("Mask exceeds maximum number of simulated qubits.");
+            }
+            uint32_t mask_num_qubits = (uint32_t)s_mask_num_qubits;
+            uint32_t mask_batch_size = (uint32_t)s_mask_batch_size;
+
+            self.ensure_safe_to_do_circuit_with_stats(CircuitStats{.num_qubits = mask_num_qubits});
+            auto u = arr.unchecked<2>();
+            bool p_x = (0b0110 >> p) & 1;  // parity of 2 bit number
+            bool p_z = p & 2;
+            for (size_t i = 0; i < mask_num_qubits; i++) {
+                for (size_t j = 0; j < mask_batch_size; j++) {
+                    bool b = *u.data(i, j);
+                    self.x_table[i][j] ^= b & p_x;
+                    self.z_table[i][j] ^= b & p_z;
+                }
+            }
+        },
+        pybind11::kw_only(),
+        pybind11::arg("pauli"),
+        pybind11::arg("mask"),
+        clean_doc_string(R"DOC(
+            @signature def broadcast_pauli_errors(self, *, pauli: Union[str, int], mask: np.ndarray) -> None:
+            Applies a pauli error to all qubits in all instances, filtered by a mask.
+
+            Args:
+                pauli: The pauli, specified as an integer or string.
+                    Uses the convention 0=I, 1=X, 2=Y, 3=Z.
+                    Any value from [0, 1, 2, 3, 'X', 'Y', 'Z', 'I', '_'] is allowed.
+                mask: A 2d numpy array specifying where to apply errors. The first axis
+                    is qubits, the second axis is simulation instances. The first axis
+                    can have a length less than the current number of qubits (or more,
+                    which adds qubits to the simulation). The length of the second axis
+                    must match the simulator's `batch_size`. The array must satisfy
+
+                        mask.dtype == np.bool_
+                        len(mask.shape) == 2
+                        mask.shape[1] == flip_sim.batch_size
+
+                    The error is only applied to qubit q in instance k when
+
+                        mask[q, k] == True.
+
+            Examples:
+                >>> import stim
+                >>> import numpy as np
+                >>> sim = stim.FlipSimulator(
+                ...     batch_size=2,
+                ...     num_qubits=3,
+                ...     disable_stabilizer_randomization=True,
+                ... )
+                >>> sim.broadcast_pauli_errors(
+                ...     pauli='X',
+                ...     mask=np.asarray([[True, False],[False, False],[True, True]]),
+                ... )
+                >>> sim.peek_pauli_flips()
+                [stim.PauliString("+X_X"), stim.PauliString("+__X")]
+
+                >>> sim.broadcast_pauli_errors(
+                ...     pauli='Z',
+                ...     mask=np.asarray([[False, True],[False, False],[True, True]]),
+                ... )
+                >>> sim.peek_pauli_flips()
+                [stim.PauliString("+X_Y"), stim.PauliString("+Z_Y")]
+
+        )DOC")
+            .data());
 }
diff --git a/src/stim/simulators/frame_simulator_pybind_test.py b/src/stim/simulators/frame_simulator_pybind_test.py
index 945d7acfc..0078a8d8e 100644
--- a/src/stim/simulators/frame_simulator_pybind_test.py
+++ b/src/stim/simulators/frame_simulator_pybind_test.py
@@ -191,15 +191,15 @@ def test_set_pauli_flip():
         stim.PauliString('XZ_'),
     ]
 
-    with pytest.raises(ValueError, match='Expected pauli'):
+    with pytest.raises(ValueError, match='pauli'):
         sim.set_pauli_flip(-1, qubit_index=0, instance_index=0)
-    with pytest.raises(ValueError, match='Expected pauli'):
+    with pytest.raises(ValueError, match='pauli'):
         sim.set_pauli_flip(4, qubit_index=0, instance_index=0)
-    with pytest.raises(ValueError, match='Expected pauli'):
+    with pytest.raises(ValueError, match='pauli'):
         sim.set_pauli_flip('R', qubit_index=0, instance_index=0)
-    with pytest.raises(ValueError, match='Expected pauli'):
+    with pytest.raises(ValueError, match='pauli'):
         sim.set_pauli_flip('XY', qubit_index=0, instance_index=0)
-    with pytest.raises(ValueError, match='Expected pauli'):
+    with pytest.raises(ValueError, match='pauli'):
         sim.set_pauli_flip(object(), qubit_index=0, instance_index=0)
 
     with pytest.raises(IndexError, match='instance_index'):
@@ -216,6 +216,193 @@ def test_set_pauli_flip():
         stim.PauliString('XZ___'),
     ]
 
+def test_broadcast_pauli_errors():
+    sim = stim.FlipSimulator(
+        batch_size=2,
+        num_qubits=3,
+        disable_stabilizer_randomization=True,
+    )
+    sim.broadcast_pauli_errors(
+        pauli='X',
+        mask=np.asarray([
+            [True, False],
+            [False, False],
+            [True, True]]
+        ),
+    )
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+X_X"),
+        stim.PauliString("+__X")
+    ]
+    sim.broadcast_pauli_errors(
+        pauli='Z',
+        mask=np.asarray([
+            [True, True],
+            [True, False],
+            [False, False]]
+        ),
+    )
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+YZX"),
+        stim.PauliString("+Z_X")
+    ]
+    sim.broadcast_pauli_errors(
+        pauli='Y',
+        mask=np.asarray([
+            [True, False],
+            [False, True],
+            [False, True]]
+        ),
+    )
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+_ZX"),
+        stim.PauliString("+ZYZ")
+    ]
+    sim.broadcast_pauli_errors(
+        pauli='I',
+        mask=np.asarray([
+            [True, True],
+            [False, True],
+            [True, True]]
+        ),
+    )
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+_ZX"),
+        stim.PauliString("+ZYZ")
+    ]
+
+    # do it again with ints
+    sim = stim.FlipSimulator(
+        batch_size=2,
+        num_qubits=3,
+        disable_stabilizer_randomization=True,
+    )
+    sim.broadcast_pauli_errors(
+        pauli=1,
+        mask=np.asarray([
+            [True, False],
+            [False, False],
+            [True, True]]
+        ),
+    )
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+X_X"),
+        stim.PauliString("+__X")
+    ]
+    sim.broadcast_pauli_errors(
+        pauli=3,
+        mask=np.asarray([
+            [True, True],
+            [True, False],
+            [False, False]]
+        ),
+    )
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+YZX"),
+        stim.PauliString("+Z_X")
+    ]
+    sim.broadcast_pauli_errors(
+        pauli=2,
+        mask=np.asarray([
+            [True, False],
+            [False, True],
+            [False, True]]
+        ),
+    )
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+_ZX"),
+        stim.PauliString("+ZYZ")
+    ]
+    sim.broadcast_pauli_errors(
+        pauli=0,
+        mask=np.asarray([
+            [True, True],
+            [False, True],
+            [True, True]]
+        ),
+    )
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+_ZX"),
+        stim.PauliString("+ZYZ")
+    ]
+
+    with pytest.raises(ValueError, match='pauli'):
+        sim.broadcast_pauli_errors(
+            pauli='whoops',
+            mask=np.asarray([
+                [True, True],
+                [False, True],
+                [True, True]]
+            ),
+        )
+    with pytest.raises(ValueError, match='pauli'):
+        sim.broadcast_pauli_errors(
+            pauli=4,
+            mask=np.asarray([
+                [True, True],
+                [False, True],
+                [True, True]]
+            ),
+        )
+    with pytest.raises(ValueError, match='batch_size'):
+        sim.broadcast_pauli_errors(
+            pauli='X',
+            mask=np.asarray([
+                [True, True,True],
+                [False, True, True],
+                [True, True, True]]
+            ),
+        )
+    with pytest.raises(ValueError, match='batch_size'):
+        sim.broadcast_pauli_errors(
+            pauli='X',
+            mask=np.asarray([
+                [True],
+                [False],
+                [True]]
+            ),
+        )
+    sim = stim.FlipSimulator(
+        batch_size=2,
+        num_qubits=3,
+        disable_stabilizer_randomization=True,
+    )
+    sim.broadcast_pauli_errors(
+        pauli='X',
+        mask=np.asarray([
+            [True, False],
+            [False, False],
+            [True, True],
+            [True, True]]
+        ),
+    ) # automatically expands the qubit basis
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+X_XX"),
+        stim.PauliString("+__XX")
+    ]
+    sim.broadcast_pauli_errors(
+        pauli='X',
+        mask=np.asarray([
+            [True, False],
+            [False, False],
+            ]
+        ),
+    )  # tolerates fewer qubits in mask than in simulator
+    peek = sim.peek_pauli_flips()
+    assert peek == [
+        stim.PauliString("+__XX"),
+        stim.PauliString("+__XX")
+    ]
+
 
 def test_repro_heralded_pauli_channel_1_bug():
     circuit = stim.Circuit("""