Scaling factors

cpp_test/TestBPDecoder.cpp

@@ -631,6 +631,82 @@ TEST(BpDecoder, DefaultScheduleAndConstantWhenRandomSerialScheduleFalse) {
     ASSERT_EQ(first_schedule, second_schedule) << "Schedule changed despite random_serial_schedule being false.";
 }
 
+TEST(BpDecoder, DynamicScalingFactorSetup) {
+    int n = 5;
+    auto pcm = ldpc::gf2codes::rep_code<ldpc::bp::BpEntry>(n);
+    int maximum_iterations = 10;
+    auto channel_probabilities = vector<double>(pcm.n, 0.1);
+
+    // Initialize decoder with dynamic scaling factor damping
+    double min_sum_scaling_factor = 0.5;
+    double dynamic_scaling_factor_damping = 0.1;
+    auto decoder = ldpc::bp::BpDecoder(pcm, channel_probabilities, maximum_iterations, ldpc::bp::MINIMUM_SUM,
+                                       ldpc::bp::PARALLEL, min_sum_scaling_factor, 1, ldpc::bp::NULL_INT_VECTOR, 0, false, ldpc::bp::AUTO, dynamic_scaling_factor_damping);
+
+    // Verify that the scaling factors are set up correctly
+    ASSERT_EQ(decoder.ms_scaling_factor_vector.size(), maximum_iterations);
+    for (int i = 0; i < maximum_iterations; i++) {
+        double expected_factor = 1.0 - (1.0 - min_sum_scaling_factor) * std::pow(2.0, -1 * i * dynamic_scaling_factor_damping);
+        ASSERT_NEAR(decoder.ms_scaling_factor_vector[i], expected_factor, 1e-6);
+    }
+}
+
+TEST(BpDecoder, StaticScalingFactorSetup) {
+    int n = 5;
+    auto pcm = ldpc::gf2codes::rep_code<ldpc::bp::BpEntry>(n);
+    int maximum_iterations = 10;
+    auto channel_probabilities = vector<double>(pcm.n, 0.1);
+
+    // Initialize decoder with static scaling factor
+    double min_sum_scaling_factor = 0.5;
+    auto decoder = ldpc::bp::BpDecoder(pcm, channel_probabilities, maximum_iterations, ldpc::bp::MINIMUM_SUM,
+                                       ldpc::bp::PARALLEL, min_sum_scaling_factor);
+
+    // Verify that the scaling factors are set up correctly
+    ASSERT_EQ(decoder.ms_scaling_factor_vector.size(), maximum_iterations);
+    for (int i = 0; i < maximum_iterations; i++) {
+        ASSERT_NEAR(decoder.ms_scaling_factor_vector[i], min_sum_scaling_factor, 1e-6);
+    }
+}
+
+
+TEST(BpDecoder, MsConvergeValueSetup) {
+    int n = 5;
+    auto pcm = ldpc::gf2codes::rep_code<ldpc::bp::BpEntry>(n);
+    int maximum_iterations = 10;
+    auto channel_probabilities = vector<double>(pcm.n, 0.1);
+
+    // Initialize decoder with ms_converge value
+    double ms_converge = 0.01;
+    auto decoder = ldpc::bp::BpDecoder(pcm, channel_probabilities, maximum_iterations, ldpc::bp::MINIMUM_SUM,
+                                       ldpc::bp::PARALLEL, 0.5, 1, ldpc::bp::NULL_INT_VECTOR, 0, false, ldpc::bp::AUTO, 0.0, ms_converge);
+
+    // Verify that the ms_converge value is set correctly
+    ASSERT_EQ(decoder.ms_converge_value, 0.01);
+}
+
+TEST(BpDecoder, SetUpMsScalingFactorsTest) {
+    int n = 5;
+    auto pcm = ldpc::gf2codes::rep_code<ldpc::bp::BpEntry>(n);
+    int maximum_iterations = 10;
+    auto channel_probabilities = vector<double>(pcm.n, 0.1);
+
+    // Initialize decoder with ms_scaling_factor=1 and ms_converge_value=2.0
+    double ms_scaling_factor = 1.0;
+    double ms_converge_value = 2.0;
+    auto decoder = ldpc::bp::BpDecoder(pcm, channel_probabilities, maximum_iterations, ldpc::bp::MINIMUM_SUM,
+                                       ldpc::bp::PARALLEL, ms_scaling_factor, 1, ldpc::bp::NULL_INT_VECTOR, 0, false,
+                                       ldpc::bp::AUTO, 0.1, ms_converge_value);
+
+    // Verify that the scaling factors are set up correctly
+    ASSERT_EQ(decoder.ms_scaling_factor_vector.size(), maximum_iterations);
+    for (int i = 0; i < maximum_iterations; i++) {
+        double expected_factor = ms_converge_value - (ms_converge_value - ms_scaling_factor) * std::pow(2.0, -1 * i * 0.1);
+        ASSERT_NEAR(decoder.ms_scaling_factor_vector[i], expected_factor, 1e-6);
+    }
+}
+
+
 int main(int argc, char **argv) {
     ::testing::InitGoogleTest(&argc, argv);
     return RUN_ALL_TESTS();

pyproject.toml

@@ -25,7 +25,7 @@ dependencies = [
     "sinter>=1.12.0",
     "pymatching"
 ]
-version = "2.3.8"
+version = "2.3.9"
 
 [project.urls]
 Documentation = "https://software.roffe.eu/ldpc"

python_test/test_bp_dynamic_scaling_factor.py

@@ -0,0 +1,93 @@
+import pytest
+import numpy as np
+from ldpc.bp_decoder._bp_decoder import BpDecoder
+from ldpc import BpOsdDecoder
+from ldpc.bplsd_decoder import BpLsdDecoder
+from ldpc.belief_find_decoder import BeliefFindDecoder
+from ldpc.bp_flip import BpFlipDecoder
+
+def test_dynamic_scaling_factor_damping_initialization():
+    pcm = np.array([[1, 1, 0], [0, 1, 1]], dtype=np.uint8)
+    channel_probs = [0.1, 0.2, 0.3]
+    max_iter = 10
+    damping_factor = 0.1
+
+    decoder = BpDecoder(pcm, channel_probs=channel_probs, bp_method="ms", max_iter=max_iter, dynamic_scaling_factor_damping=damping_factor)
+
+    assert decoder.dynamic_scaling_factor_damping == damping_factor, "Dynamic scaling factor damping not set correctly."
+
+def test_dynamic_scaling_factor_damping_effect():
+    pcm = np.array([[1, 1, 0], [0, 1, 1]], dtype=np.uint8)
+    channel_probs = [0.1, 0.2, 0.3]
+    max_iter = 10
+    damping_factor = 0.1
+
+    decoder = BpDecoder(pcm, channel_probs=channel_probs, bp_method="ms", max_iter=max_iter, dynamic_scaling_factor_damping=damping_factor)
+
+    # Verify that the scaling factors are computed correctly
+    expected_factors = [
+        1.0 - (1.0 - decoder.ms_scaling_factor) * (2.0 ** (-1 * i * damping_factor))
+        for i in range(max_iter)
+    ]
+    for i, factor in enumerate(expected_factors):
+        assert pytest.approx(decoder.ms_scaling_factor_vector[i], rel=1e-6) == factor, f"Scaling factor mismatch at iteration {i}."
+
+def test_dynamic_scaling_factor_with_initial_and_converge_parameters():
+    pcm = np.array([[1, 1, 0], [0, 1, 1]], dtype=np.uint8)
+    channel_probs = [0.1, 0.2, 0.3]
+
+    max_iter = 10
+    damping_factor = 0.1 # Initial damping factor for dynamic scaling factor growth
+    ms_scaling_factor = 0.5  # Initial scaling factor for testing
+    ms_converge_value = 2.0  # Convergence value for dynamic scaling factor
+
+    decoder = BpDecoder(pcm, channel_probs=channel_probs, max_iter=max_iter,  bp_method="ms", ms_scaling_factor=ms_scaling_factor, dynamic_scaling_factor_damping=damping_factor, ms_converge_value=ms_converge_value)
+
+    print("ms_scaling_factor_start:", decoder.ms_scaling_factor)
+    print("damping_factor:", damping_factor)
+    print("ms_converge_value:", decoder.ms_converge_value)
+    print("Initial scaling factors:", decoder.ms_scaling_factor_vector)
+
+
+    # Verify that the scaling factors are recomputed correctly
+    expected_factors = [
+        ms_converge_value - (ms_converge_value - ms_scaling_factor) * (2.0 ** (-1 * i * damping_factor))
+        for i in range(max_iter)
+    ]
+    for i, factor in enumerate(expected_factors):
+        assert pytest.approx(decoder.ms_scaling_factor_vector[i], rel=1e-6) == factor, f"Scaling factor mismatch at iteration {i} after update."
+
+def test_dynamic_scaling_factor_damping_bplsd():
+    pcm = np.array([[1, 1, 0], [0, 1, 1]], dtype=np.uint8)
+    damping_factor = 0.2
+
+    decoder = BpLsdDecoder(pcm, error_rate = 0.1, bp_method="ms", dynamic_scaling_factor_damping=damping_factor)
+    assert decoder.dynamic_scaling_factor_damping == damping_factor, "Damping factor not set correctly."
+
+    updated_damping_factor = 0.5
+    decoder.dynamic_scaling_factor_damping = updated_damping_factor
+    assert decoder.dynamic_scaling_factor_damping == updated_damping_factor, "Damping factor update failed."
+
+
+def test_dynamic_scaling_factor_damping_bposd():
+    pcm = np.array([[1, 1, 0], [0, 1, 1]], dtype=np.uint8)
+    damping_factor = 0.3
+
+    decoder = BpOsdDecoder(pcm, error_rate=0.1, bp_method="ms", dynamic_scaling_factor_damping=damping_factor)
+    assert decoder.dynamic_scaling_factor_damping == damping_factor, "Damping factor not set correctly."
+
+    updated_damping_factor = 0.6
+    decoder.dynamic_scaling_factor_damping = updated_damping_factor
+    assert decoder.dynamic_scaling_factor_damping == updated_damping_factor, "Damping factor update failed."
+
+def test_dynamic_scaling_factor_damping_belief_find():
+    pcm = np.array([[1, 1, 0], [0, 1, 1]], dtype=np.uint8)
+    damping_factor = 0.4
+
+    decoder = BeliefFindDecoder(pcm, error_rate=0.1, bp_method="ms", dynamic_scaling_factor_damping=damping_factor)
+    assert decoder.dynamic_scaling_factor_damping == damping_factor, "Damping factor not set correctly."
+
+    updated_damping_factor = 0.7
+    decoder.dynamic_scaling_factor_damping = updated_damping_factor
+    assert decoder.dynamic_scaling_factor_damping == updated_damping_factor, "Damping factor update failed."
+

python_test/test_bplsd.py

@@ -151,6 +151,7 @@ def test_rep_code_ms():
     H = rep_code(3)
 
     lsd = BpLsdDecoder(H, error_rate=0.1, bp_method="min_sum", ms_scaling_factor=1.0)
+
     assert lsd is not None
     assert lsd.bp_method == "minimum_sum"
     assert lsd.schedule == "parallel"
@@ -190,3 +191,7 @@ def test_stats_reset():
     assert len(stats["bit_llrs"]) == 0
     assert len(stats["individual_cluster_stats"]) == 0
     assert len(stats["global_timestep_bit_history"]) == 0
+
+
+if __name__ == "__main__":
+    test_rep_code_ms()

python_test/test_qcodes.py

@@ -239,7 +239,8 @@ def test_400_16_6_hgp():
         error_rate=error_rate,
         max_iter=max_iter,
         bp_method="ms",
-        ms_scaling_factor=0,
+        ms_scaling_factor=0.5,
+        dynamic_scaling_factor_damping=1.0,
         schedule="parallel",
         bits_per_step=1,
         lsd_order=0,
@@ -260,7 +261,8 @@ def test_400_16_6_hgp():
         error_rate=error_rate,
         max_iter=max_iter,
         bp_method="ms",
-        ms_scaling_factor=0,
+        ms_scaling_factor=0.55,
+        dynamic_scaling_factor_damping=0.1,
         schedule="serial",
         bits_per_step=1,
         lsd_order=0,

src_cpp/bp.hpp

@@ -11,6 +11,7 @@
 #include <chrono>
 #include <stdexcept> // required for std::runtime_error
 #include <set>
+#include <iostream>
 
 #include "math.h"
 #include "sparse_matrix_base.hpp"
@@ -60,9 +61,11 @@ namespace ldpc {
             BpSchedule schedule;
             BpInputType bp_input_type;
             double ms_scaling_factor;
+            std::vector<double> ms_scaling_factor_vector;
+            double dynamic_scaling_factor_damping;
+            double ms_converge_value;
             std::vector<uint8_t> decoding;
             std::vector<uint8_t> candidate_syndrome;
-
             std::vector<double> log_prob_ratios;
             std::vector<double> initial_log_prob_ratios;
             std::vector<double> soft_syndrome;
@@ -85,10 +88,14 @@ namespace ldpc {
                     const std::vector<int> &serial_schedule = NULL_INT_VECTOR,
                     int random_schedule_seed = 0,
                     bool random_serial_schedule = false,
-                    BpInputType bp_input_type = AUTO) :
+                    BpInputType bp_input_type = AUTO,
+                    double dynamic_scaling_factor_damping = -1.0,
+                    double ms_converge_value = 1.0) : 
                     pcm(parity_check_matrix), channel_probabilities(std::move(channel_probabilities)),
                     check_count(pcm.m), bit_count(pcm.n), maximum_iterations(maximum_iterations), bp_method(bp_method),
                     schedule(schedule), ms_scaling_factor(min_sum_scaling_factor),
+                    dynamic_scaling_factor_damping(dynamic_scaling_factor_damping),
+                    ms_converge_value(ms_converge_value),
                     iterations(0) //the parity check matrix is passed in by reference
             {
 
@@ -126,13 +133,40 @@ namespace ldpc {
                     }
                 }
 
+
+                this->set_up_ms_scaling_factors();
+
                 //Initialise OMP thread pool
                 // this->omp_thread_count = omp_threads;
                 // this->set_omp_thread_count(this->omp_thread_count);
             }
 
             ~BpDecoder() = default;
 
+
+            void set_up_ms_scaling_factors(){
+
+                if(this->bp_method == MINIMUM_SUM){
+
+                    if(this->dynamic_scaling_factor_damping <= 0) {
+                        this->ms_scaling_factor_vector.resize(this->maximum_iterations);
+                        for (int i = 0; i < this->maximum_iterations; i++) {
+                            this->ms_scaling_factor_vector[i] = this->ms_scaling_factor;
+                        }
+                    } else {
+                        this->ms_scaling_factor_vector.resize(this->maximum_iterations);
+                        for (int i = 0; i < this->maximum_iterations; i++) {
+                            this->ms_scaling_factor_vector[i] = this->ms_converge_value - (this->ms_converge_value - this->ms_scaling_factor) * std::pow(2.0, -1*i*this->dynamic_scaling_factor_damping);
+                        }
+                    } 
+
+                } else {
+                    this->ms_scaling_factor_vector.clear();
+                }
+
+            }
+
+
             void set_omp_thread_count(int count) {
                 this->omp_thread_count = count;
                 // omp_set_num_threads(this->omp_thread_count);
@@ -191,6 +225,7 @@ namespace ldpc {
 
             std::vector<uint8_t> &bp_decode_parallel(std::vector<uint8_t> &syndrome) {
 
+
                 this->converge = 0;
 
                 this->initialise_log_domain_bp();
@@ -219,13 +254,9 @@ namespace ldpc {
                         }
                     } else if (this->bp_method == MINIMUM_SUM) {
 
-                        double alpha;
-                        if(this->ms_scaling_factor == 0.0) {
-                            alpha = 1.0 - std::pow(2.0, -1.0*it);
-                        }
-                        else {
-                            alpha = this->ms_scaling_factor;
-                        }
+
+                        double alpha = this->ms_scaling_factor_vector[it - 1];
+
 
                         //check to bit updates
                         for (int i = 0; i < check_count; i++) {
@@ -456,14 +487,6 @@ namespace ldpc {
 
                 for (int it = 1; it <= maximum_iterations; it++) {
 
-                    double alpha;
-                    if(this->ms_scaling_factor == 0.0) {
-                        alpha = 1.0 - std::pow(2.0, -1.0*it);
-                    }
-                    else {
-                        alpha = this->ms_scaling_factor;
-                    }
-
                     if (this->random_serial_schedule) {
                         this->rng_list_shuffle.shuffle(this->serial_schedule_order);
                     } else if (this->schedule == BpSchedule::SERIAL_RELATIVE) {
@@ -500,6 +523,7 @@ namespace ldpc {
                                 this->log_prob_ratios[bit_index] += e.check_to_bit_msg;
                             }
                         } else if (this->bp_method == 1) {
+                            double alpha = this->ms_scaling_factor_vector[it - 1];
                             for (auto &e: pcm.iterate_column(bit_index)) {
                                 check_index = e.row_index;
                                 int sgn = syndrome[check_index];

src_python/ldpc/belief_find_decoder/__init__.pyi

@@ -2,7 +2,6 @@ import numpy as np
 from scipy.sparse import spmatrix
 
 class BeliefFindDecoder(BpDecoderBase):
-
     """
     A class representing a decoder that combines Belief Propagation (BP) with the Union Find Decoder (UFD) algorithm.
 
@@ -41,18 +40,17 @@ class BeliefFindDecoder(BpDecoderBase):
         The inversion method can be applied to any parity check matrix.
     bits_per_step : int, optional
         Specifies the number of bits added to the cluster in each step of the UFD algorithm. If no value is provided, this is set the block length of the code.
-
-    Notes
-    -----
-    The `BeliefFindDecoder` class leverages soft information outputted by the BP decoder to guide the cluster growth
-    in the UFD algorithm. The number of bits added to the cluster in each step is controlled by the `bits_per_step` parameter.
-    The `uf_method` parameter activates a more general version of the UFD algorithm suitable for LDPC codes when set to True.
+    dynamic_scaling_factor_damping : Optional[float], optional
+        The damping factor for dynamic scaling in the minimum sum method, by default -1.0.
+    ms_converge_value (Optional[float]):
+        Convergence value for the minimum-sum method.
     """
 
     def __cinit__(self, pcm: Union[np.ndarray, scipy.sparse.spmatrix], error_rate: Optional[float] = None,
                  error_channel: Optional[List[float]] = None, max_iter: Optional[int] = 0, bp_method: Optional[str] = 'minimum_sum',
                  ms_scaling_factor: Optional[float] = 1.0, schedule: Optional[str] = 'parallel', omp_thread_count: Optional[int] = 1,
-                 random_schedule_seed: Optional[int] = 0, serial_schedule_order: Optional[List[int]] = None, uf_method: str = "peeling", bits_per_step:int = 0, input_vector_type: str = "syndrome"): ...
+                 random_schedule_seed: Optional[int] = 0, serial_schedule_order: Optional[List[int]] = None, uf_method: str = "peeling",
+                 bits_per_step: int = 0, input_vector_type: str = "syndrome", dynamic_scaling_factor_damping: Optional[float] = -1.0, ms_converge_value: Optional[float] = 1.0, **kwargs): ...
 
     def __del__(self): ...