25 namespace seqan3::detail
44 template <
typename alignment_configuration_t,
typename ...policies_t>
46 requires is_type_specialisation_of_v<alignment_configuration_t, configuration>
48 class pairwise_alignment_algorithm :
protected policies_t...
52 using traits_type = alignment_configuration_traits<alignment_configuration_t>;
54 using score_type =
typename traits_type::score_type;
56 using alignment_result_type =
typename traits_type::alignment_result_type;
58 static_assert(!std::same_as<alignment_result_type, empty_type>,
"Alignment result type was not configured.");
65 pairwise_alignment_algorithm() =
default;
66 pairwise_alignment_algorithm(pairwise_alignment_algorithm
const &) =
default;
67 pairwise_alignment_algorithm(pairwise_alignment_algorithm &&) =
default;
68 pairwise_alignment_algorithm & operator=(pairwise_alignment_algorithm
const &) =
default;
69 pairwise_alignment_algorithm & operator=(pairwise_alignment_algorithm &&) =
default;
70 ~pairwise_alignment_algorithm() =
default;
79 pairwise_alignment_algorithm(alignment_configuration_t
const & config) : policies_t(config)...
128 template <indexed_sequence_pair_range indexed_sequence_pairs_t,
typename callback_t>
130 requires std::invocable<callback_t, alignment_result_type>
132 void operator()(indexed_sequence_pairs_t && indexed_sequence_pairs, callback_t && callback)
136 for (
auto && [sequence_pair, idx] : indexed_sequence_pairs)
138 size_t const sequence1_size = std::ranges::distance(get<0>(sequence_pair));
139 size_t const sequence2_size = std::ranges::distance(get<1>(sequence_pair));
141 auto && [alignment_matrix, index_matrix] = this->acquire_matrices(sequence1_size, sequence2_size);
143 compute_matrix(get<0>(sequence_pair), get<1>(sequence_pair), alignment_matrix, index_matrix);
144 this->make_result_and_invoke(
std::forward<decltype(sequence_pair)>(sequence_pair),
147 this->optimal_coordinate,
154 template <indexed_sequence_pair_range indexed_sequence_pairs_t,
typename callback_t>
156 requires traits_type::is_vectorised && std::invocable<callback_t, alignment_result_type>
158 auto operator()(indexed_sequence_pairs_t && indexed_sequence_pairs, callback_t && callback)
160 using simd_collection_t =
std::vector<score_type, aligned_allocator<score_type,
alignof(score_type)>>;
161 using original_score_t =
typename traits_type::original_score_type;
164 auto seq1_collection = indexed_sequence_pairs | views::get<0> | views::get<0>;
165 auto seq2_collection = indexed_sequence_pairs | views::get<0> | views::get<1>;
167 this->initialise_tracker(seq1_collection, seq2_collection);
170 thread_local simd_collection_t simd_seq1_collection{};
171 thread_local simd_collection_t simd_seq2_collection{};
173 convert_batch_of_sequences_to_simd_vector(simd_seq1_collection, seq1_collection, traits_type::padding_symbol);
174 convert_batch_of_sequences_to_simd_vector(simd_seq2_collection, seq2_collection, traits_type::padding_symbol);
176 size_t const sequence1_size = std::ranges::distance(simd_seq1_collection);
177 size_t const sequence2_size = std::ranges::distance(simd_seq2_collection);
179 auto && [alignment_matrix, index_matrix] = this->acquire_matrices(sequence1_size, sequence2_size);
181 compute_matrix(simd_seq1_collection, simd_seq2_collection, alignment_matrix, index_matrix);
184 for (
auto && [sequence_pair, idx] : indexed_sequence_pairs)
186 original_score_t score = this->optimal_score[index] -
187 (this->padding_offsets[index] * this->scoring_scheme.padding_match_score());
188 matrix_coordinate coordinate{row_index_type{
size_t{this->optimal_coordinate.row[index]}},
189 column_index_type{
size_t{this->optimal_coordinate.col[index]}}};
190 this->make_result_and_invoke(
std::forward<decltype(sequence_pair)>(sequence_pair),
219 template <
typename simd_sequence_t,
220 std::ranges::forward_range sequence_collection_t,
223 requires std::ranges::output_range<simd_sequence_t, score_type>
225 void convert_batch_of_sequences_to_simd_vector(simd_sequence_t & simd_sequence,
226 sequence_collection_t & sequences,
227 padding_symbol_t
const & padding_symbol)
229 assert(
static_cast<size_t>(std::ranges::distance(sequences)) <= traits_type::alignments_per_vector);
231 simd_sequence.clear();
232 for (
auto && simd_vector_chunk : sequences | views::to_simd<score_type>(padding_symbol))
249 template <std::ranges::forward_range sequence1_t,
250 std::ranges::forward_range sequence2_t,
251 std::ranges::input_range alignment_matrix_t,
252 std::ranges::input_range index_matrix_t>
254 requires std::ranges::forward_range<std::ranges::range_reference_t<alignment_matrix_t>> &&
255 std::ranges::forward_range<std::ranges::range_reference_t<index_matrix_t>>
257 void compute_matrix(sequence1_t && sequence1,
258 sequence2_t && sequence2,
259 alignment_matrix_t && alignment_matrix,
260 index_matrix_t && index_matrix)
266 this->reset_optimum();
268 auto alignment_matrix_it = alignment_matrix.begin();
269 auto indexed_matrix_it = index_matrix.begin();
271 initialise_column(*alignment_matrix_it, *indexed_matrix_it, sequence2);
277 for (
auto sequence1_value : sequence1)
278 compute_column(*++alignment_matrix_it, *++indexed_matrix_it, sequence1_value, sequence2);
284 auto && alignment_column = *alignment_matrix_it;
285 auto && cell_index_column = *indexed_matrix_it;
287 auto alignment_column_it = alignment_column.begin();
288 auto cell_index_column_it = cell_index_column.begin();
290 this->track_last_column_cell(*alignment_column_it, *cell_index_column_it);
292 for ([[maybe_unused]]
auto && unused : sequence2)
293 this->track_last_column_cell(*++alignment_column_it, *++cell_index_column_it);
295 this->track_final_cell(*alignment_column_it, *cell_index_column_it);
316 template <std::ranges::input_range alignment_column_t,
317 std::ranges::input_range cell_index_column_t,
318 std::ranges::input_range sequence2_t>
319 void initialise_column(alignment_column_t && alignment_column,
320 cell_index_column_t && cell_index_column,
321 sequence2_t && sequence2)
327 auto first_column_it = alignment_column.begin();
328 auto cell_index_column_it = cell_index_column.begin();
329 *first_column_it = this->track_cell(this->initialise_origin_cell(), *cell_index_column_it);
335 for ([[maybe_unused]]
auto const & unused : sequence2)
338 *first_column_it = this->track_cell(this->initialise_first_column_cell(*first_column_it),
339 *++cell_index_column_it);
346 this->track_last_row_cell(*first_column_it, *cell_index_column_it);
367 template <std::ranges::input_range alignment_column_t,
368 std::ranges::input_range cell_index_column_t,
369 typename sequence1_value_t,
370 std::ranges::input_range sequence2_t>
374 void compute_column(alignment_column_t && alignment_column,
375 cell_index_column_t && cell_index_column,
376 sequence1_value_t
const & sequence1_value,
377 sequence2_t && sequence2)
379 using score_type =
typename traits_type::score_type;
385 auto alignment_column_it = alignment_column.begin();
386 auto cell_index_column_it = cell_index_column.begin();
388 auto cell = *alignment_column_it;
389 score_type diagonal = cell.best_score();
390 *alignment_column_it = this->track_cell(this->initialise_first_row_cell(cell), *cell_index_column_it);
396 for (
auto const & sequence2_value : sequence2)
398 auto cell = *++alignment_column_it;
399 score_type next_diagonal = cell.best_score();
400 *alignment_column_it = this->track_cell(
401 this->compute_inner_cell(diagonal, cell, this->scoring_scheme.score(sequence1_value, sequence2_value)),
402 *++cell_index_column_it);
403 diagonal = next_diagonal;
410 this->track_last_row_cell(*alignment_column_it, *cell_index_column_it);