protogalaxy_prover_internal.hpp

Go to the documentation of this file.

// === AUDIT STATUS ===
// internal:    { status: not started, auditors: [], date: YYYY-MM-DD }
// external_1:  { status: not started, auditors: [], date: YYYY-MM-DD }
// external_2:  { status: not started, auditors: [], date: YYYY-MM-DD }
// =====================
 
#pragma once
#include "barretenberg/common/bb_bench.hpp"
#include "barretenberg/common/container.hpp"
#include "barretenberg/common/thread.hpp"
#include "barretenberg/honk/execution_trace/execution_trace_usage_tracker.hpp"
#include "barretenberg/protogalaxy/constants.hpp"
#include "barretenberg/protogalaxy/prover_verifier_shared.hpp"
#include "barretenberg/relations/relation_parameters.hpp"
#include "barretenberg/relations/relation_types.hpp"
#include "barretenberg/relations/utils.hpp"
#include "barretenberg/ultra_honk/oink_prover.hpp"
 
namespace bb {
 
template <class ProverInstance> class ProtogalaxyProverInternal {
  public:
    using Flavor = typename ProverInstance::Flavor;
    using FF = typename Flavor::FF;
    using RelationUtils = bb::RelationUtils<Flavor>;
    using ProverPolynomials = typename Flavor::ProverPolynomials;
    using Relations = typename Flavor::Relations;
    using AllValues = typename Flavor::AllValues;
    using SubrelationSeparators = typename Flavor::SubrelationSeparators;
    using ProverInstances = std::array<std::shared_ptr<ProverInstance>, NUM_INSTANCES>;
 
    static constexpr size_t EXTENDED_LENGTH = computed_extended_length<Flavor>();
    static constexpr size_t BATCHED_EXTENDED_LENGTH = computed_batched_extended_length<Flavor>();
    static constexpr size_t NUM_SUBRELATIONS = Flavor::NUM_SUBRELATIONS;
 
    using UnivariateRelationParametersNoOptimisticSkipping = bb::RelationParameters<Univariate<FF, EXTENDED_LENGTH>>;
    using UnivariateRelationParameters =
        bb::RelationParameters<Univariate<FF, EXTENDED_LENGTH, 0, /*skip_count=*/SKIP_COUNT>>;
    using UnivariateSubrelationSeparators = std::array<Univariate<FF, BATCHED_EXTENDED_LENGTH>, NUM_SUBRELATIONS - 1>;
 
    // Univariates that interpolate polynomial evaluations at a given vertex across two instances
    using ExtendedUnivariate = Univariate<FF, EXTENDED_LENGTH>;
    // Combiner univariate
    using ExtendedUnivariateWithRandomization = Univariate<FF, BATCHED_EXTENDED_LENGTH>;
 
    using ShortUnivariates = typename Flavor::template ProverUnivariates<NUM_INSTANCES>;
 
    using ExtendedUnivariates =
        typename Flavor::template ProverUnivariatesWithOptimisticSkipping<ExtendedUnivariate::LENGTH,
                                                                          /*skip_count=*/SKIP_COUNT>;
 
    using ExtendedUnivariatesType =
        std::conditional_t<Flavor::USE_SHORT_MONOMIALS, ShortUnivariates, ExtendedUnivariates>;
 
    using TupleOfTuplesOfUnivariates = typename Flavor::template ProtogalaxyTupleOfTuplesOfUnivariates<NUM_INSTANCES>;
    using TupleOfTuplesOfUnivariatesNoOptimisticSkipping =
        typename Flavor::template ProtogalaxyTupleOfTuplesOfUnivariatesNoOptimisticSkipping<NUM_INSTANCES>;
 
    using RelationEvaluations = decltype(create_tuple_of_arrays_of_values<typename Flavor::Relations>());
 
    ExecutionTraceUsageTracker trace_usage_tracker;
 
    ProtogalaxyProverInternal(ExecutionTraceUsageTracker trace_usage_tracker = ExecutionTraceUsageTracker{})
        : trace_usage_tracker(std::move(trace_usage_tracker))
    {}
 
    template <size_t skip_count = 0> static auto row_to_univariates(const ProverInstances& instances, size_t row_idx)
    {
        using ContainerType =
            std::conditional_t<Flavor::USE_SHORT_MONOMIALS,
                               typename Flavor::template ProverUnivariates<2>,
                               std::array<Univariate<FF, ExtendedUnivariate::LENGTH, 0, skip_count>, NUM_INSTANCES>>;
        // As a practical measure, get the first prover instance's view to deduce the array type
        std::array<decltype(instances[0]->polynomials.get_all()), NUM_INSTANCES> views;
        views[0] = instances[0]->polynomials.get_all();
        views[1] = instances[1]->polynomials.get_all();
 
        ContainerType results;
        // Set the size corresponding to the number of rows in the execution trace
        // Iterate over the prover polynomials' views corresponding to each prover instance
        for (size_t inst_idx = 0; auto& get_all : views) {
            // Iterate over all columns in the trace execution of an prover instance and extract their value at row_idx.
            if constexpr (Flavor::USE_SHORT_MONOMIALS) {
                // In this case, the elements of the ContainerType are AllEntities, so we need to get the underlying
                // polynomials via get_all()
                for (auto [result, poly_ptr] : zip_view(results.get_all(), get_all)) {
                    result.evaluations[inst_idx] = poly_ptr[row_idx];
                }
            } else {
                for (auto [result, poly_ptr] : zip_view(results, get_all)) {
                    result.evaluations[inst_idx] = poly_ptr[row_idx];
                }
            }
            inst_idx++;
        }
        return results;
    }
 
    inline static FF process_subrelation_evaluations(const RelationEvaluations& evals,
                                                     const SubrelationSeparators& challenges,
                                                     FF& linearly_dependent_contribution)
    {
        // Initialize result with the contribution from the first subrelation
        FF linearly_independent_contribution = std::get<0>(evals)[0];
        size_t idx = 0;
 
        auto scale_by_challenge_and_accumulate =
            [&]<size_t relation_idx, size_t subrelation_idx, typename Element>(Element& element) {
                if constexpr (!(relation_idx == 0 && subrelation_idx == 0)) {
                    using Relation = typename std::tuple_element_t<relation_idx, Relations>;
                    // Accumulate scaled subrelation contribution
                    const Element contribution = element * challenges[idx++];
                    if constexpr (subrelation_is_linearly_independent<Relation, subrelation_idx>()) {
                        linearly_independent_contribution += contribution;
                    } else {
                        linearly_dependent_contribution += contribution;
                    }
                }
            };
        RelationUtils::apply_to_tuple_of_arrays_elements(scale_by_challenge_and_accumulate, evals);
        return linearly_independent_contribution;
    }
 
    Polynomial<FF> compute_row_evaluations(const ProverPolynomials& polynomials,
                                           const SubrelationSeparators& alphas,
                                           const RelationParameters<FF>& relation_parameters)
 
    {
 
        BB_BENCH_NAME("ProtogalaxyProver_::compute_row_evaluations");
 
        const size_t polynomial_size = polynomials.get_polynomial_size();
        Polynomial<FF> aggregated_relation_evaluations(polynomial_size);
 
        // Determine the number of threads over which to distribute the work
        const size_t num_threads = compute_num_threads(polynomial_size);
 
        std::vector<FF> linearly_dependent_contribution_accumulators(num_threads);
 
        // Distribute the execution trace rows across threads so that each handles an equal number of active rows
        trace_usage_tracker.construct_thread_ranges(num_threads, polynomial_size, /*use_prev_accumulator=*/true);
 
        parallel_for(num_threads, [&](size_t thread_idx) {
            for (const ExecutionTraceUsageTracker::Range& range : trace_usage_tracker.thread_ranges[thread_idx]) {
                for (size_t idx = range.first; idx < range.second; idx++) {
                    const AllValues row = polynomials.get_row(idx);
                    // Evaluate all subrelations on given row. Separator is 1 since we are not summing across rows here.
                    const RelationEvaluations evals =
                        RelationUtils::accumulate_relation_evaluations(row, relation_parameters, FF(1));
 
                    // Sum against challenges alpha
                    aggregated_relation_evaluations.at(idx) = process_subrelation_evaluations(
                        evals, alphas, linearly_dependent_contribution_accumulators[thread_idx]);
                }
            }
        });
 
        aggregated_relation_evaluations.at(0) += sum(linearly_dependent_contribution_accumulators);
 
        return aggregated_relation_evaluations;
    }
    static std::vector<FF> construct_coefficients_tree(std::span<const FF> betas,
                                                       std::span<const FF> deltas,
                                                       std::vector<FF>& current_level_buffer,
                                                       size_t current_width,
                                                       size_t current_degree)
    {
        if (betas.size() == 1) {
            return std::vector<FF>(current_level_buffer.begin(),
                                   current_level_buffer.begin() + static_cast<std::ptrdiff_t>(current_degree));
        }
 
        size_t num_levels = betas.size();
 
        // Pre-allocate buffers with maximum size needed
        // At each level i, we have width/(2^i) nodes, each with degree i+1
        // Calculate total size needed for the next level buffer
        size_t max_next_buffer_size = 0;
        size_t temp_width = current_width;
        for (size_t level = 1; level < num_levels; ++level) {
            temp_width /= 2;
            size_t degree = level + 1;
            max_next_buffer_size = std::max(max_next_buffer_size, temp_width * (degree + 1));
        }
 
        std::vector<FF> next_level_buffer(max_next_buffer_size);
 
        // Iterate through levels 1 to num_levels-1
        for (size_t level = 1; level < num_levels; ++level) {
            size_t next_degree = level + 1;
            size_t next_width = current_width / 2;
            size_t next_stride = next_degree + 1;
            size_t current_stride = current_degree;
 
            // Parallel initialization and computation
            parallel_for_heuristic(
                next_width,
                [&](size_t parent) {
                    size_t node = parent * 2;
                    size_t parent_offset = parent * next_stride;
                    size_t left_child_offset = node * current_stride;
                    size_t right_child_offset = (node + 1) * current_stride;
 
                    // Copy left child coefficients and zero-initialize remaining
                    for (size_t d = 0; d < current_stride; d++) {
                        next_level_buffer[parent_offset + d] = current_level_buffer[left_child_offset + d];
                    }
                    for (size_t d = current_stride; d < next_stride; d++) {
                        next_level_buffer[parent_offset + d] = FF(0);
                    }
 
                    // Add right child contribution: n_r * (β_i + δ_i X)
                    for (size_t d = 0; d < current_stride; d++) {
                        next_level_buffer[parent_offset + d] +=
                            current_level_buffer[right_child_offset + d] * betas[level];
                        next_level_buffer[parent_offset + d + 1] +=
                            current_level_buffer[right_child_offset + d] * deltas[level];
                    }
                },
                /* overestimate */ thread_heuristics::FF_MULTIPLICATION_COST * next_degree * 3);
 
            // Swap buffers: next becomes current
            std::swap(current_level_buffer, next_level_buffer);
            current_width = next_width;
            current_degree = next_stride;
        }
 
        // Return the final result (the root node)
        return std::vector<FF>(current_level_buffer.begin(),
                               current_level_buffer.begin() + static_cast<std::ptrdiff_t>(current_degree));
    }
 
    static std::vector<FF> construct_perturbator_coefficients(std::span<const FF> betas,
                                                              std::span<const FF> deltas,
                                                              const Polynomial<FF>& full_honk_evaluations)
    {
        BB_BENCH();
 
        auto width = full_honk_evaluations.size();
 
        // Construct first level coefficients (degree 1 polynomials from leaf pairs)
        // Use flat buffer: each node has 2 coefficients (degree 1)
        size_t first_level_width = width / 2;
        std::vector<FF> first_level_buffer(first_level_width * 2);
 
        parallel_for_heuristic(
            first_level_width,
            [&](size_t parent) {
                size_t node = parent * 2;
                size_t offset = parent * 2;
                first_level_buffer[offset] = full_honk_evaluations[node] + full_honk_evaluations[node + 1] * betas[0];
                first_level_buffer[offset + 1] = full_honk_evaluations[node + 1] * deltas[0];
            },
            /* overestimate */ thread_heuristics::FF_MULTIPLICATION_COST * 3);
 
        // Build the tree iteratively in-place
        return construct_coefficients_tree(betas, deltas, first_level_buffer, first_level_width, 2);
    }
 
    Polynomial<FF> compute_perturbator(const std::shared_ptr<const ProverInstance>& accumulator,
                                       const std::vector<FF>& deltas)
    {
        BB_BENCH();
        auto full_honk_evaluations =
            compute_row_evaluations(accumulator->polynomials, accumulator->alphas, accumulator->relation_parameters);
        const auto betas = accumulator->gate_challenges;
        BB_ASSERT_EQ(betas.size(), deltas.size());
        const size_t log_circuit_size = accumulator->log_dyadic_size();
 
        // Compute the perturbator using only the first log_circuit_size-many betas/deltas
        std::vector<FF> perturbator = construct_perturbator_coefficients(std::span{ betas.data(), log_circuit_size },
                                                                         std::span{ deltas.data(), log_circuit_size },
                                                                         full_honk_evaluations);
 
        // Populate the remaining coefficients with zeros to reach the required constant size
        for (size_t idx = log_circuit_size; idx < CONST_PG_LOG_N; ++idx) {
            perturbator.emplace_back(FF(0));
        }
 
        // Check that the perturbator zeroth coefficient is equal to the target sum stored in the accumulator
        BB_ASSERT_EQ(perturbator[0],
                     accumulator->target_sum,
                     "ProtogalaxyProver: the zeroth coefficient of the perturbator is different from the target sum "
                     "stored in the accumulator.");
 
        return Polynomial<FF>{ perturbator };
    }
 
    template <size_t skip_count = 0>
    BB_INLINE static void extend_univariates(ExtendedUnivariatesType& extended_univariates,
                                             const ProverInstances& instances,
                                             const size_t row_idx)
    {
        if constexpr (Flavor::USE_SHORT_MONOMIALS) {
            extended_univariates = std::move(row_to_univariates(instances, row_idx));
        } else {
            auto incoming_univariates = row_to_univariates<skip_count>(instances, row_idx);
            for (auto [extended_univariate, incoming_univariate] :
                 zip_view(extended_univariates.get_all(), incoming_univariates)) {
                incoming_univariate.template self_extend_from<NUM_INSTANCES>();
                extended_univariate = std::move(incoming_univariate);
            }
        }
    }
 
    template <size_t relation_idx = 0>
    BB_INLINE static void accumulate_relation_univariates(TupleOfTuplesOfUnivariates& univariate_accumulators,
                                                          const ExtendedUnivariatesType& extended_univariates,
                                                          const UnivariateRelationParameters& relation_parameters,
                                                          const FF& scaling_factor)
    {
        using Relation = std::tuple_element_t<relation_idx, Relations>;
 
        //  Check if the relation is skippable to speed up accumulation
        if constexpr (!isSkippable<Relation, decltype(extended_univariates)>) {
            // If not, accumulate normally
            Relation::accumulate(std::get<relation_idx>(univariate_accumulators),
                                 extended_univariates,
                                 relation_parameters,
                                 scaling_factor);
        } else {
            // If so, only compute the contribution if the relation is active
            if (!Relation::skip(extended_univariates)) {
                Relation::accumulate(std::get<relation_idx>(univariate_accumulators),
                                     extended_univariates,
                                     relation_parameters,
                                     scaling_factor);
            }
        }
 
        // Repeat for the next relation.
        if constexpr (relation_idx + 1 < Flavor::NUM_RELATIONS) {
            accumulate_relation_univariates<relation_idx + 1>(
                univariate_accumulators, extended_univariates, relation_parameters, scaling_factor);
        }
    }
 
    ExtendedUnivariateWithRandomization compute_combiner(const ProverInstances& instances,
                                                         const GateSeparatorPolynomial<FF>& gate_separators,
                                                         const UnivariateRelationParameters& relation_parameters,
                                                         const UnivariateSubrelationSeparators& alphas,
                                                         TupleOfTuplesOfUnivariates& univariate_accumulators)
    {
        BB_BENCH();
 
        // Determine the number of threads over which to distribute the work
        // The polynomial size is given by the virtual size since the computation includes
        // the incoming key which could have nontrivial values on the larger domain in case of overflow.
        const size_t common_polynomial_size = instances[0]->polynomials.w_l.virtual_size();
        const size_t num_threads = compute_num_threads(common_polynomial_size);
 
        using ThreadAccumulators = TupleOfTuplesOfUnivariates;
 
        // Construct univariate accumulator containers; one per thread
        // Note: std::vector will trigger {}-initialization of the contents. Therefore no need to zero the univariates.
        std::vector<ThreadAccumulators> thread_univariate_accumulators(num_threads);
 
        // Distribute the execution trace rows across threads so that each handles an equal number of active rows
        trace_usage_tracker.construct_thread_ranges(num_threads, common_polynomial_size);
 
        // Accumulate the contribution from each sub-relation
        parallel_for(num_threads, [&](size_t thread_idx) {
            // Construct extended univariates containers; one per thread
            ExtendedUnivariatesType extended_univariates;
 
            for (const ExecutionTraceUsageTracker::Range& range : trace_usage_tracker.thread_ranges[thread_idx]) {
                for (size_t idx = range.first; idx < range.second; idx++) {
                    // Instantiate extended univariates: they contain the evaluations of the prover polynomials from the
                    // instances being folded and are extended to length EXTENDED_LENGTH, which is the maximum length of
                    // a folded subrelation polynomial. We set skip_count to NUM_INSTANCES - 1 = 1 because the
                    // relations evaluate to zero on valid keys, and the evaluations of the combiner on the point 1,
                    // is exactly the relation contributions coming from the key being folded (note that the Univariate
                    // class skips computing evaluations between 1 and skip_count).
                    //
                    // The values of the evaluations are still available for skipping relations, but all
                    // derived univariates (coming from addition, multiplication, etc) will ignore those evaluations.
                    //
                    // No need to initialize extended_univariates to 0, as it's assigned to a value in the following
                    // function.
                    extend_univariates<SKIP_COUNT>(extended_univariates, instances, idx);
 
                    const FF pow_challenge = gate_separators[idx];
 
                    // Accumulate the i-th row's univariate contribution. Note that the relation parameters passed
                    // to this function have already been folded. Moreover, linear-dependent relations that act over
                    // the entire execution trace rather than on rows will not be multiplied by the pow challenge.
                    accumulate_relation_univariates(thread_univariate_accumulators[thread_idx],
                                                    extended_univariates,
                                                    relation_parameters, // these parameters have already been folded
                                                    pow_challenge);
                }
            }
        });
 
        RelationUtils::zero_univariates(univariate_accumulators);
        // Accumulate the per-thread univariate accumulators into a single set of accumulators
        for (auto& accumulators : thread_univariate_accumulators) {
            RelationUtils::add_nested_tuples(univariate_accumulators, accumulators);
        }
        // This does nothing if TupleOfTuples is TupleOfTuplesOfUnivariates
        TupleOfTuplesOfUnivariatesNoOptimisticSkipping deoptimized_univariates =
            deoptimize_univariates(univariate_accumulators);
        // Batch the univariate contributions from each folded sub-relation using the folded batching challenges to
        // obtain the combiner
        return batch_over_relations(deoptimized_univariates, alphas);
    }
 
    template <typename TupleOfTuplesOfUnivariatePossiblyOptimistic>
    static TupleOfTuplesOfUnivariatesNoOptimisticSkipping deoptimize_univariates(
        const TupleOfTuplesOfUnivariatePossiblyOptimistic& tup)
    {
        // If input does not have optimized operators, return the input
        if constexpr (std::same_as<TupleOfTuplesOfUnivariatePossiblyOptimistic,
                                   TupleOfTuplesOfUnivariatesNoOptimisticSkipping>) {
            return tup;
        }
 
        const auto deoptimize = [&]<size_t outer_idx, size_t inner_idx>(auto& element) {
            auto& element_with_skipping = std::get<inner_idx>(std::get<outer_idx>(tup));
            element = element_with_skipping.convert();
        };
 
        // Note: {} is required to initialize the tuple contents. Otherwise the univariates contain garbage.
        TupleOfTuplesOfUnivariatesNoOptimisticSkipping result{};
        RelationUtils::apply_to_tuple_of_tuples(result, deoptimize);
        return result;
    }
 
    static ExtendedUnivariateWithRandomization batch_over_relations(
        TupleOfTuplesOfUnivariatesNoOptimisticSkipping& univariate_accumulators,
        const UnivariateSubrelationSeparators& alphas)
    {
        auto result = std::get<0>(std::get<0>(univariate_accumulators)).template extend_to<BATCHED_EXTENDED_LENGTH>();
 
        size_t idx = 0;
        const auto scale_and_sum = [&]<size_t outer_idx, size_t inner_idx>(auto& element) {
            if constexpr (outer_idx == 0 && inner_idx == 0) {
                return;
            }
 
            // Extend folded subrelation polynomials to the length of the combiner
            auto extended = element.template extend_to<BATCHED_EXTENDED_LENGTH>();
            extended *= alphas[idx];
            result += extended;
            idx++;
        };
 
        RelationUtils::apply_to_tuple_of_tuples(univariate_accumulators, scale_and_sum);
        RelationUtils::zero_univariates(univariate_accumulators);
 
        return result;
    }
 
    static std::pair<typename ProverInstance::FF, std::array<typename ProverInstance::FF, NUM_INSTANCES>>
    compute_vanishing_polynomial_and_lagranges(const FF& challenge)
    {
        FF vanishing_polynomial_at_challenge;
        std::array<FF, NUM_INSTANCES> lagranges;
        vanishing_polynomial_at_challenge = challenge * (challenge - FF(1));
        lagranges = { FF(1) - challenge, challenge };
 
        return { vanishing_polynomial_at_challenge, lagranges };
    }
 
    static Univariate<FF, BATCHED_EXTENDED_LENGTH, NUM_INSTANCES> compute_combiner_quotient(
        FF perturbator_evaluation, ExtendedUnivariateWithRandomization combiner)
    {
        std::array<FF, BATCHED_EXTENDED_LENGTH - NUM_INSTANCES> combiner_quotient_evals = {};
 
        for (size_t point = NUM_INSTANCES; point < combiner.size(); point++) {
            auto idx = point - NUM_INSTANCES;
            FF lagrange_0 = FF(1) - FF(point);
            FF vanishing_polynomial = FF(point) * (FF(point) - 1);
            combiner_quotient_evals[idx] =
                (combiner.value_at(point) - perturbator_evaluation * lagrange_0) * vanishing_polynomial.invert();
        }
 
        return Univariate<FF, BATCHED_EXTENDED_LENGTH, NUM_INSTANCES>(combiner_quotient_evals);
    }
 
    template <typename ExtendedRelationParameters>
    static ExtendedRelationParameters compute_extended_relation_parameters(const ProverInstances& instances)
    {
        using UnivariateParameter = typename ExtendedRelationParameters::DataType;
        ExtendedRelationParameters result;
        size_t param_idx = 0;
        for (auto& param : result.get_to_fold()) {
            Univariate<FF, NUM_INSTANCES> tmp(0);
            tmp.value_at(0) = instances[0]->relation_parameters.get_to_fold()[param_idx];
            tmp.value_at(1) = instances[1]->relation_parameters.get_to_fold()[param_idx];
            param = tmp.template extend_to<UnivariateParameter::LENGTH, UnivariateParameter::SKIP_COUNT>();
            param_idx++;
        }
        return result;
    }
 
    static UnivariateSubrelationSeparators compute_and_extend_alphas(const ProverInstances& instances)
    {
        UnivariateSubrelationSeparators result;
        size_t alpha_idx = 0;
        for (auto& alpha : result) {
            Univariate<FF, NUM_INSTANCES> tmp;
            tmp.value_at(0) = instances[0]->alphas[alpha_idx];
            tmp.value_at(1) = instances[1]->alphas[alpha_idx];
            alpha = tmp.template extend_to<BATCHED_EXTENDED_LENGTH>();
            alpha_idx++;
        }
        return result;
    }
 
    static size_t compute_num_threads(const size_t domain_size)
    {
        const size_t max_num_threads = get_num_cpus_pow2(); // number of available threads (power of 2)
        const size_t min_iterations_per_thread =
            1 << 6; // min number of iterations for which we'll spin up a unique thread
        const size_t desired_num_threads = domain_size / min_iterations_per_thread;
        size_t num_threads = std::min(desired_num_threads, max_num_threads); // fewer than max if justified
        num_threads = num_threads > 0 ? num_threads : 1;                     // ensure num threads is >= 1
 
        return num_threads;
    }
};
} // namespace bb