C++ Subset类代码示例

void debugDisplay(const Subset& SS)
{
	std::cout << " " << SS.count() << "\t" << SS.size() << "\t";
	for(unsigned int i = 0; i < SS.count(); ++i)
		std::cout << SS.pop() << " ";
	std::cout << std::endl;
}

bool sameSubset(const Subset &a, const Subset &b)
{
  if (a.size() != b.size()) return false;
  for (size_t i = 0; i < a.size(); i++)
    if (a[i] != b[i]) return false;
  return true;
}

/*
	Rule II: If B contains more elements than C then S(B) > S(C).
*/
bool validate2(const Subset &x, const Subset &y) {
	if(x.size() > y.size())
		return (S(x) > S(y));
	else if(x.size() < y.size())
		return (S(x) < S(y));
	else
		return false;
}

void load_particle_states(const Subset &s,
                          const Assignment &ss,
                          const ParticleStatesTable *pst) {
  load_particle_states(s.begin(), s.end(), ss, pst);
  if (s.size()!=0) {
    s[0]->get_model()->update();
  }
}

/**
 * Compute the Kramer-Mesner matrix.
 *
 * Precondition: t < k
 */
KramerMesnerMatrix computeKMMatrix(const Group& G, unsigned int t, unsigned int k) {
	std::vector<KMBuilderOutput> builderOutputs;			// stores the relevant data for k = 2 onwards
	Matrix A;
	
	for (int i = 2; i <= k; i++) {
		boost::scoped_ptr<KMBuilder> builder;
		
		// Get orbit representatives of (i - 1)-subsets
		std::vector<Subset> orbitReps;
		if (i == 2) {
			// As it is the first time, we have to compute the orbit representatives of singleton subsets...
			Subset pointsRemaining = generateX(G.getNumPoints());
			while (!pointsRemaining.empty()) {
				typedef OrbitSet<Permutation, unsigned long>::const_iterator OrbitSetIterator;
				
				unsigned long point = *(pointsRemaining.begin());	// A new representative...
				Subset singletonSet = makeSingletonSet(point);		// ...as a set
				orbitReps.push_back(singletonSet);
				
				// Compute the orbit of the point, and remove it from the remaining points
				OrbitSet<Permutation, unsigned long> orbit = G.orbit(point, Transversal::TrivialAction());
				for (OrbitSetIterator it = orbit.begin(); it != orbit.end(); it++) {
					pointsRemaining.erase(pointsRemaining.find(*it));
				}
			}
			
			builder.reset(new KMBuilder(G, i, orbitReps));
		} else {
			// We get them from what we computed earlier
			KMBuilderOutput& input = builderOutputs[i - 3];
			orbitReps = input.getNewReps();
			
			builder.reset(new KMBuilder(G, i, orbitReps, input.getPrunerData()));
		}
		
		KMBuilderOutput builderOutput = builder->build();
		builderOutputs.push_back(builderOutput);
		
		// From the identity that A[t][k] = A[t][s] * A[s][k] / combinat(k - t, k - s) for any s
		// between t and k, we get the following:
		if (i == t + 1) {
			assignMatrix(A, builderOutput.getNewMatrix());
		} else if (i > t + 1) {
			assignMatrix(A, matrixMultiply(A, builderOutput.getNewMatrix()));
			scalarDivide(i - t, A);
		}
	}
	
	return KramerMesnerMatrix(G, builderOutputs[t - 2].getNewReps(), builderOutputs[k - 2].getNewReps(), A);
}

float TSPProblem::solveTSP()
{
    allocateTSPCache_();

    //base cases
    {
        for (unsigned int i=0; i<rows_; i++) {
            cache_[i][0] = std::numeric_limits<float>::max();
        }
        Subset s; s.add(0);
        cache_[s.repr()][0] = 0; // from node to itself using exactly itself
    }

    //meat of the algorithm
    for (unsigned int m=2; m <= cities_.size(); m++) {      //subset size (including 0)
        vector<Subset> subsets;
        getSubsets(subsets, m, cities_.size());

        for (const Subset& s : subsets) {

            for (unsigned int j=1; j<cities_.size(); j++) { //omit start node 0
                if (!(s.contains(j))) {                     // omit stuff not in subset
                    continue;
                }

                Subset sprime(s);
                sprime.remove(j);

                float minval = std::numeric_limits<float>::max();
                for (unsigned int k=0; k<cities_.size(); k++) {
                    if (!(s.contains(k)) || k == j) {
                        continue;
                    }
                    float localCost = cache_[sprime.repr()][k] + distance(cities_[k], cities_[j]);
                    if (localCost < minval) {
                        minval = localCost;
                    }
                }
                cache_[s.repr()][j] = minval;
            }
        }
    }

    //finally find minimum of all paths computed so far
    float minval = std::numeric_limits<float>::max();
    {
        Subset s;
        for (unsigned int i=0; i<cities_.size() ; i++) {
            s.add(i);
        }
        for (unsigned int j=1; j < cities_.size(); j++) {
            float localCost = cache_[s.repr()][j] + distance(cities_[j], cities_[0]);
            if (localCost < minval) {
                minval = localCost;
            }
        }
    }

    return minval;
}

void draw(const Subset &s) {
	int i;
	
	cout << "{";
	
	for(i = 0; i < s.size(); i++) {
		cout << s[i];
		
		if(i + 1 != s.size())
			cout << ", ";
	}
	
	cout << "}";
}

Assignment get_merged_assignment(const Subset &s, const Assignment &ss0,
                                 const Ints &i0, const Assignment &ss1,
                                 const Ints &i1) {
  Ints ret(s.size(), -1);
  IMP_USAGE_CHECK(ss0.size() == i0.size(),
                  "The size of the subset and "
                      << "the index don't match: " << ss0.size() << " vs "
                      << i0.size());
  IMP_USAGE_CHECK(ss1.size() == i1.size(),
                  "The size of the subset and "
                      << "the index don't match: " << ss1.size() << " vs "
                      << i1.size());
  for (unsigned int i = 0; i < i0.size(); ++i) {
    ret[i0[i]] = ss0[i];
  }
  for (unsigned int i = 0; i < i1.size(); ++i) {
    ret[i1[i]] = ss1[i];
  }
  IMP_IF_CHECK(USAGE) {
    for (unsigned int i = 0; i < ret.size(); ++i) {
      IMP_USAGE_CHECK(ret[i] >= 0, "Not all set");
    }
  }
  return Assignment(ret);
}

Subset MinRepPruner::minRep(const Subset& candidate) {
	typedef boost::dynamic_bitset<unsigned long> dset;
	
	dset candidateBitset(G->getNumPoints());
	for (Subset::const_iterator it = candidate.begin(); it != candidate.end(); ++it) {
		candidateBitset[*it] = true;
	}
	
	dset resultBitset = search.lexMin(candidateBitset);
	
	Subset result;
	for (dset::size_type i = resultBitset.find_first(); i != dset::npos; i = resultBitset.find_next(i)) {
		result.insert(i);
	}
	return result;
}

Assignment
get_nearest_assignment(const Subset &s,
                       const algebra::VectorKD &embedding,
                       ParticleStatesTable *pst) {

  Ints ret(s.size());
  unsigned int cur=0;
  // kind of a hack to get size
  for (unsigned int i=0; i< s.size(); ++i) {
    unsigned int sz
        = pst->get_particle_states(s[i])->get_embedding(0).get_dimension();
    algebra::VectorKD cpt(embedding.coordinates_begin()+cur,
                          embedding.coordinates_begin()+cur+sz);
    cur+=sz;
    ret[i]= pst->get_particle_states(s[i])->get_nearest_state(cpt);
  }
  return Assignment(ret);
}

/*
	Example:
		f({1, 2, 3}, {4}) = true because 4 is not in the other list
		f({1, 2, 3}, {2}) = false because 2 is in the other list
*/
bool uniek(const Subset &a, const Subset &b) {
	int x, y;
	
	//Init
	x = 0;
	y = 0;
	
	while(x < a.size() && y < b.size()) {
		if(a[x] == b[y])
			return false;
		
		if(a[x] < b[y])
			x++;
		else
			y++;
	}
	
	return true;
}

display::Geometries get_subset_graph_geometry(const SubsetGraph &ig) {
  display::Geometries ret;
  SubsetGraphConstVertexName vm = boost::get(boost::vertex_name, ig);
  for (std::pair<InteractionGraphTraits::vertex_iterator,
                 InteractionGraphTraits::vertex_iterator> be =
           boost::vertices(ig);
       be.first != be.second; ++be.first) {
    Subset s = vm[*be.first];
    display::Color c = display::get_display_color(*be.first);
    for (unsigned int i = 0; i < s.size(); ++i) {
      core::XYZ pi(s[i]);
      IMP_NEW(display::SphereGeometry, cg,
              (algebra::Sphere3D(pi.get_coordinates(), 1)));
      cg->set_color(c);
      cg->set_name(s.get_name());
      ret.push_back(cg.get());
    }
  }
  return ret;
}

 algebra::VectorKD get_embedding(const Subset &s,
                                 const Assignment &a,
                                 ParticleStatesTable *pst){
   Floats embed;
  for (unsigned int i=0; i< s.size(); ++i) {
    algebra::VectorKD cur
        = pst->get_particle_states(s[i])->get_embedding(a[i]);
    embed.insert(embed.end(), cur.coordinates_begin(),
                 cur.coordinates_end());
  }
  return embed;
 }

RestraintsTemp get_restraints(const Subset &s,
                              const ParticleStatesTable *pst,
                              const DependencyGraph &dg,
                              RestraintSet *rs) {
  RestraintsTemp rw= get_restraints(RestraintsTemp(1, rs));
  Subset other=pst->get_subset();
  ParticlesTemp oms;
  std::set_difference(other.begin(), other.end(),
                      s.begin(), s.end(),
                      std::back_inserter(oms));
  IMP::compatibility::map<Restraint*, int> index
    = IMP::base::internal::get_graph_index<Restraint>(dg);
  Ints to_remove;
  for (unsigned int i=0; i< rw.size(); ++i) {
    if (IMP::internal::get_has_ancestor(dg, index[rw[i]], oms)) {
      to_remove.push_back(i);
    }
  }
  for (int i=to_remove.size()-1; i >=0; --i) {
    rw.erase(rw.begin()+to_remove[i]);
  }
  return rw;
}

Assignments get_state_clusters(const Subset &subset, const Assignments &states,
                               ParticleStatesTable *pst, double resolution) {
  Vector<Ints> rotated(subset.size(), Ints(states.size(), -1));
  for (unsigned int i = 0; i < states.size(); ++i) {
    for (unsigned int j = 0; j < states[i].size(); ++j) {
      rotated[j][i] = states[i][j];
    }
  }
  for (unsigned int i = 0; i < rotated.size(); ++i) {
    std::sort(rotated[i].begin(), rotated[i].end());
    rotated[i].erase(std::unique(rotated[i].begin(), rotated[i].end()),
                     rotated[i].end());
  }
  Vector<Ints> clustering(states.size());
  for (unsigned int i = 0; i < subset.size(); ++i) {
    IMP::PointerMember<ParticleStates> ps =
        pst->get_particle_states(subset[i]);
    Ints c = get_state_clusters(subset[i], ps, rotated[i], resolution);
    clustering[i] = c;
  }

  return filter_states(subset, states, clustering);
}