29 std::vector<int> &sizes,
31 std::vector<KDT *> &correspondence_kdt_map,
32 std::vector<double> *min_diag_price,
33 std::vector<double> *min_price,
34 std::vector<std::vector<MatchingType>> *all_matchings,
36 bool actual_distance) {
37 Timer const time_matchings;
39 double local_cost = *total_cost;
40#ifdef TTK_ENABLE_OPENMP
41#pragma omp parallel for num_threads(threadNumber_) schedule(dynamic, 1) reduction(+:local_cost)
43 for(
int i = 0; i < numberOfInputs_; i++) {
44 double const delta_lim = delta_lim_;
46 current_bidder_diagrams_[i], barycenter_goods_[i], wasserstein_,
47 geometrical_factor_, lambda_, delta_lim, kdt, correspondence_kdt_map,
48 epsilon, min_diag_price->at(i), use_kdt, nonMatchingWeight_);
57 min_price->at(i) = getMinimalPrice(i);
58 std::vector<MatchingType> matchings;
60 all_matchings->at(i) = matchings;
62 local_cost += sqrt(cost);
69 precision_[i] = quotient < 1 ? 1. / sqrt(1 - quotient) - 1 : 10;
75 current_bidder_diagrams_[i].resize(sizes[i]);
77 *total_cost = local_cost;
82 std::vector<int> &sizes,
84 std::vector<KDT *> &correspondence_kdt_map,
85 std::vector<double> *min_diag_price,
86 std::vector<std::vector<MatchingType>> *all_matchings,
88 bool actual_distance) {
89 double local_cost = *total_cost;
90#ifdef TTK_ENABLE_OPENMP
91#pragma omp parallel for num_threads(threadNumber_) schedule(dynamic, 1) reduction(+:local_cost)
93 for(
int i = 0; i < numberOfInputs_; i++) {
95 current_bidder_diagrams_[i], barycenter_goods_[i], wasserstein_,
96 geometrical_factor_, lambda_, delta_lim_, kdt, correspondence_kdt_map, 0,
97 (*min_diag_price)[i], use_kdt, nonMatchingWeight_);
98 std::vector<MatchingType> matchings;
99 double const cost = auction.
run(matchings, i);
100 all_matchings->at(i) = matchings;
101 if(actual_distance) {
102 local_cost += sqrt(cost);
110 current_bidder_diagrams_[i].resize(sizes[i]);
112 *total_cost = local_cost;
137 std::vector<std::vector<MatchingType>> &previous_matchings) {
139 std::vector<std::vector<MatchingType>> corrected_matchings(numberOfInputs_);
140 for(
int i = 0; i < numberOfInputs_; i++) {
142 std::vector<int> new_to_old_id(current_bidder_diagrams_[i].size());
143 for(
size_t j = 0; j < current_bidder_ids_[i].size(); j++) {
144 int const new_id = current_bidder_ids_[i][j];
146 new_to_old_id[new_id] = j;
150 std::vector<MatchingType> matchings_diagram_i;
151 for(
size_t j = 0; j < previous_matchings[i].size(); j++) {
153 int const new_id = std::get<0>(m);
154 if(new_id >= 0 && std::get<1>(m) >= 0) {
155 std::get<0>(m) = new_to_old_id[new_id];
156 matchings_diagram_i.push_back(m);
159 corrected_matchings[i] = matchings_diagram_i;
161 return corrected_matchings;
165 std::vector<std::vector<MatchingType>> &matchings) {
167 Timer const t_update;
168 size_t const n_goods = barycenter_goods_[0].size();
170 size_t const n_diagrams = current_bidder_diagrams_.size();
173 double max_shift = 0;
175 std::vector<size_t> count_diag_matchings(
177 std::vector<double> x(n_goods);
178 std::vector<double> y(n_goods);
179 std::vector<double> crit_coords_x(n_goods);
180 std::vector<double> crit_coords_y(n_goods);
181 std::vector<double> crit_coords_z(n_goods);
182 for(
size_t i = 0; i < n_goods; i++) {
183 count_diag_matchings[i] = 0;
186 crit_coords_x[i] = 0;
187 crit_coords_y[i] = 0;
188 crit_coords_z[i] = 0;
190 std::vector<double> min_prices(n_diagrams);
191 for(
size_t j = 0; j < n_diagrams; j++) {
192 min_prices[j] = std::numeric_limits<double>::max();
195 std::vector<Bidder *>
198 for(
size_t j = 0; j < matchings.size(); j++) {
199 for(
size_t i = 0; i < matchings[j].size(); i++) {
200 int const bidder_id = std::get<0>(matchings[j][i]);
201 int const good_id = std::get<1>(matchings[j][i]);
202 if(good_id < 0 && bidder_id >= 0) {
204 points_to_append.push_back(¤t_bidder_diagrams_[j].at(bidder_id));
207 else if(good_id >= 0 && bidder_id >= 0) {
209 x[good_id] += current_bidder_diagrams_[j].at(bidder_id).x_;
210 y[good_id] += current_bidder_diagrams_[j].at(bidder_id).y_;
211 if(geometrical_factor_ < 1) {
212 const auto &critical_coordinates = current_bidder_diagrams_[j]
214 .GetCriticalCoordinates();
215 crit_coords_x[good_id] += critical_coordinates[0];
216 crit_coords_y[good_id] += critical_coordinates[1];
217 crit_coords_z[good_id] += critical_coordinates[2];
219 }
else if(good_id >= 0 && bidder_id < 0) {
222 count_diag_matchings[good_id] = count_diag_matchings[good_id] + 1;
228 for(
size_t i = 0; i < n_goods; i++) {
229 if(count_diag_matchings[i] < n_diagrams) {
233 = x[i] / (double)(n_diagrams - count_diag_matchings[i]);
235 = y[i] / (double)(n_diagrams - count_diag_matchings[i]);
239 = ((double)(n_diagrams - count_diag_matchings[i]) * x_bar
240 + (double)count_diag_matchings[i] * (x_bar + y_bar) / 2.)
241 / (
double)n_diagrams;
243 = ((double)(n_diagrams - count_diag_matchings[i]) * y_bar
244 + (double)count_diag_matchings[i] * (x_bar + y_bar) / 2.)
245 / (
double)n_diagrams;
247 double const new_crit_coord_x
248 = crit_coords_x[i] / (double)(n_diagrams - count_diag_matchings[i]);
249 double const new_crit_coord_y
250 = crit_coords_y[i] / (double)(n_diagrams - count_diag_matchings[i]);
251 double const new_crit_coord_z
252 = crit_coords_z[i] / (double)(n_diagrams - count_diag_matchings[i]);
256 double const dx = barycenter_goods_[0].at(i).x_ - new_x;
257 double const dy = barycenter_goods_[0].at(i).y_ - new_y;
258 double const shift =
Geometry::pow(std::abs(dx), wasserstein_)
260 if(shift > max_shift) {
264 for(
size_t j = 0; j < n_diagrams; j++) {
265 barycenter_goods_[j].at(i).SetCoordinates(new_x, new_y);
266 if(geometrical_factor_ < 1) {
267 barycenter_goods_[j].at(i).SetCriticalCoordinates(
268 new_crit_coord_x, new_crit_coord_y, new_crit_coord_z);
270 if(barycenter_goods_[j].at(i).getPrice() < min_prices[j]) {
271 min_prices[j] = barycenter_goods_[j].at(i).getPrice();
277 for(
size_t j = 0; j < n_diagrams; j++) {
278 if(min_prices[j] >= std::numeric_limits<double>::max() / 2.) {
285 for(
size_t i = 0; i < n_goods; i++) {
286 if(count_diag_matchings[i] == n_diagrams) {
287 points_deleted_ += 1;
291 barycenter_goods_[0].at(i).getPersistence() / 2., wasserstein_);
292 if(shift > max_shift) {
295 for(
size_t j = 0; j < n_diagrams; j++) {
296 barycenter_goods_[j].at(i).id_ = -1;
302 for(
size_t k = 0; k < points_to_append.size(); k++) {
304 Bidder *b = points_to_append[k];
306 = (b->
x_ + (n_diagrams - 1) * (b->
x_ + b->
y_) / 2.) / (n_diagrams);
308 = (b->
y_ + (n_diagrams - 1) * (b->
x_ + b->
y_) / 2.) / (n_diagrams);
310 for(
size_t j = 0; j < n_diagrams; j++) {
311 Good g =
Good(gx, gy,
false, barycenter_goods_[j].size());
313 if(geometrical_factor_ < 1) {
315 std::get<1>(critical_coordinates),
316 std::get<2>(critical_coordinates));
318 barycenter_goods_[j].emplace_back(g);
322 barycenter_goods_[j].at(g.
id_).getPersistence() / 2., wasserstein_);
323 if(shift > max_shift) {
330 for(
size_t j = 0; j < n_diagrams; j++) {
333 for(
size_t i = 0; i < barycenter_goods_[j].size(); i++) {
334 Good g = barycenter_goods_[j].at(i);
337 new_barycenter.emplace_back(g);
341 barycenter_goods_[j] = new_barycenter;
382 double previous_min_persistence,
383 double min_persistence,
384 std::vector<double> &initial_diagonal_prices,
385 std::vector<double> &initial_off_diagonal_prices,
386 int min_points_to_add,
387 bool add_points_to_barycenter) {
389 double new_min_persistence = min_persistence;
393 size_t max_diagram_size{};
394 for(
int i = 0; i < numberOfInputs_; i++) {
395 if(current_bidder_diagrams_[i].size() > max_diagram_size) {
396 max_diagram_size = current_bidder_diagrams_[i].size();
399 int const max_points_to_add = std::max(
400 min_points_to_add, min_points_to_add + (
int)(max_diagram_size / 10));
403 std::vector<std::vector<int>> candidates_to_be_added(numberOfInputs_);
404 std::vector<std::vector<int>> idx(numberOfInputs_);
405 for(
int i = 0; i < numberOfInputs_; i++) {
407 std::vector<double> persistences;
408 for(
size_t j = 0; j < bidder_diagrams_[i].size(); j++) {
409 Bidder const b = bidder_diagrams_[i].at(j);
411 if(persistence >= min_persistence
412 && persistence < previous_min_persistence) {
413 candidates_to_be_added[i].push_back(j);
414 idx[i].push_back(idx[i].size());
415 persistences.push_back(persistence);
418 sort(idx[i].
begin(), idx[i].
end(), [&persistences](
int &a,
int &b) {
419 return ((persistences[a] > persistences[b])
420 || ((persistences[a] == persistences[b]) && (a > b)));
422 int const size = candidates_to_be_added[i].size();
423 if(size >= max_points_to_add) {
424 double const last_persistence_added
425 = persistences[idx[i][max_points_to_add - 1]];
426 if(last_persistence_added > new_min_persistence) {
427 new_min_persistence = last_persistence_added;
435 int compteur_for_adding_points = 0;
437 for(
int i = 0; i < numberOfInputs_; i++) {
438 int const size = candidates_to_be_added[i].size();
439 for(
int j = 0; j < std::min(max_points_to_add, size); j++) {
440 Bidder b = bidder_diagrams_[i].at(candidates_to_be_added[i][idx[i][j]]);
442 b.
id_ = current_bidder_diagrams_[i].size();
445 current_bidder_diagrams_[i].emplace_back(b);
447 current_bidder_ids_[i][candidates_to_be_added[i][idx[i][j]]]
448 = current_bidder_diagrams_[i].size() - 1;
450 int const to_be_added_to_barycenter
451 = deterministic_ ? compteur_for_adding_points % numberOfInputs_
452 : rand() % numberOfInputs_;
454 if(to_be_added_to_barycenter == 0 && add_points_to_barycenter) {
455 for(
int k = 0; k < numberOfInputs_; k++) {
456 Good g =
Good(b.
x_, b.
y_,
false, barycenter_goods_[k].size());
457 g.
setPrice(initial_off_diagonal_prices[k]);
459 barycenter_goods_[k].emplace_back(g);
463 compteur_for_adding_points++;
466 return new_min_persistence;
557 auto kdt = std::make_unique<KDT>(
true, wasserstein_);
559 const int dimension = geometrical_factor_ >= 1 ? 2 : 5;
561 std::vector<double> coordinates;
562 std::vector<std::vector<double>> weights;
564 for(
size_t i = 0; i < barycenter_goods_[0].size(); i++) {
565 const Good &g = barycenter_goods_[0].at(i);
566 coordinates.push_back(geometrical_factor_ * g.
x_);
567 coordinates.push_back(geometrical_factor_ * g.
y_);
568 if(geometrical_factor_ < 1) {
569 coordinates.push_back((1 - geometrical_factor_) * g.
coords_[0]);
570 coordinates.push_back((1 - geometrical_factor_) * g.
coords_[1]);
571 coordinates.push_back((1 - geometrical_factor_) * g.
coords_[2]);
575 for(
size_t idx = 0; idx < barycenter_goods_.size(); idx++) {
576 std::vector<double>
const empty_weights;
577 weights.push_back(empty_weights);
578 for(
size_t i = 0; i < barycenter_goods_[idx].size(); i++) {
579 const Good &g = barycenter_goods_[idx].at(i);
580 weights[idx].push_back(g.
getPrice());
585 auto correspondence_kdt_map
586 = kdt->build(coordinates.data(), barycenter_goods_[0].size(), dimension,
587 weights, barycenter_goods_.size());
590 return std::make_pair(std::move(kdt), correspondence_kdt_map);
596 std::vector<std::vector<MatchingType>> previous_matchings;
597 double const min_persistence = 0;
598 double min_cost = std::numeric_limits<double>::max();
599 int last_min_cost_obtained = 0;
601 this->setBidderDiagrams();
602 this->setInitialBarycenter(
605 double const max_persistence = getMaxPersistence();
607 std::vector<double> min_diag_price(numberOfInputs_);
608 std::vector<double> min_price(numberOfInputs_);
609 for(
int i = 0; i < numberOfInputs_; i++) {
610 min_diag_price[i] = 0;
614 int const min_points_to_add = std::numeric_limits<int>::max();
615 this->enrichCurrentBidderDiagrams(2 * max_persistence, min_persistence,
616 min_diag_price, min_price,
617 min_points_to_add,
false);
619 bool converged =
false;
620 bool finished =
false;
626 std::pair<std::unique_ptr<KDT>, std::vector<KDT *>> pair;
627 bool use_kdt =
false;
630 if(!barycenter_goods_[0].empty()) {
631 pair = this->getKDTree();
635 std::vector<std::vector<MatchingType>> all_matchings(numberOfInputs_);
636 std::vector<int> sizes(numberOfInputs_);
637 for(
int i = 0; i < numberOfInputs_; i++) {
638 sizes[i] = current_bidder_diagrams_[i].size();
644 for(
size_t j = 0; j < barycenter_goods_[0].size(); j++) {
645 Good const &g = barycenter_goods_[0].at(j);
648 diagramType_,
true});
651 bool const actual_distance = (numberOfInputs_ == 2);
652 runMatchingAuction(&total_cost, sizes, *pair.first, pair.second,
653 &min_diag_price, &all_matchings, use_kdt,
664 updateBarycenter(all_matchings);
666 if(min_cost > total_cost) {
667 min_cost = total_cost;
668 last_min_cost_obtained = 0;
670 last_min_cost_obtained += 1;
673 converged = converged || last_min_cost_obtained > 1;
674 if(numberOfInputs_ == 2)
678 previous_matchings = std::move(all_matchings);
681 for(
size_t i = 0; i < barycenter_goods_.size(); ++i) {
682 for(
size_t j = 0; j < barycenter_goods_[i].size(); ++j) {
683 barycenter_goods_[i].at(j).setPrice(0);
686 for(
size_t i = 0; i < current_bidder_diagrams_.size(); ++i) {
687 for(
size_t j = 0; j < current_bidder_diagrams_[i].size(); ++j) {
688 current_bidder_diagrams_[i].at(j).setDiagonalPrice(0);
691 for(
int i = 0; i < numberOfInputs_; i++) {
692 min_diag_price[i] = 0;
696 barycenter.resize(0);
697 for(
size_t j = 0; j < barycenter_goods_[0].size(); j++) {
698 Good const &g = barycenter_goods_[0].at(j);
701 diagramType_,
true});
705 std::vector<std::vector<MatchingType>> corrected_matchings
706 = correctMatchings(previous_matchings);
707 return corrected_matchings;