Program Listing for File change_neutral.hpp

↰ Return to documentation for file (fwdpp/sugar/change_neutral.hpp)

#ifndef FWDPP_SUGAR_CHANGE_NEUTRAL_HPP
#define FWDPP_SUGAR_CHANGE_NEUTRAL_HPP

#include <algorithm>
#include <exception>
#include <fwdpp/debug.hpp>

namespace fwdpp
{
    namespace sugar
    {
        template <typename MutationContainerType, typename mut_key_cont_t>
        void
        change_neutral_details(const MutationContainerType &mutations, const double pos,
                               const std::size_t mindex, mut_key_cont_t &a,
                               mut_key_cont_t &b)
        {
            /*
              Ask if haploid_genome has this mutation.  We do a linear
              search here in case of finite-site schemes and we
              also do not expect this function to be called a whole lot.
            */
            auto i = std::find(std::begin(a), std::end(a), mindex);
            // If it exists and is mindex, erase it
            if (i == std::end(a))
                return;

            a.erase(i);

            // Add mutation key into b
            b.insert(std::upper_bound(
                         std::begin(b), std::end(b), pos,
                         [&mutations](const double &p,
                                      const typename mut_key_cont_t::value_type &vt)

                         { return p < mutations[vt].pos; }),
                     mindex);
        }
    } // namespace sugar

    template <typename poptype>
    void
    change_neutral(poptype &p, const std::size_t mindex)
    {
        if (mindex >= p.mutations.size())
            throw std::out_of_range("mindex >= p.mutations.size()");

        bool is_neutral = p.mutations[mindex].neutral;

        // Change neutral flag
        p.mutations[mindex].neutral = !p.mutations[mindex].neutral;
        const auto pos = p.mutations[mindex].pos;

        for (auto &g : p.haploid_genomes)
            {
                if (g.n)
                    {
                        if (is_neutral)
                            {
                                sugar::change_neutral_details(p.mutations, pos, mindex,
                                                              g.mutations, g.smutations);
                            }
                        else
                            {
                                sugar::change_neutral_details(p.mutations, pos, mindex,
                                                              g.smutations, g.mutations);
                            }
                        debug::haploid_genome_data_valid(g, p.mutations, p.mcounts);
                    }
            }
    }
} // namespace fwdpp

#endif