1use std::{collections::BTreeSet, iter};
2
3use bitvec::vec::BitVec;
4use permutation::Permutation;
5use rustc_hash::FxHashSet;
6use tagged_vec::TaggedVec;
7
8use crate::{
9 index::{
10 DirectedEdgeIndex, DirectedNodeIndex, EdgeIndex, GraphIndexInteger, NodeIndex,
11 OptionalEdgeIndex,
12 },
13 io::gfa1::PlainGfaEdgeData,
14};
15
16static DEBUG_REMOVE_MULTIEDGES: bool = false;
17
18#[cfg(test)]
19mod tests;
20
21#[derive(Debug, Clone, PartialEq, Eq, PartialOrd, Ord, Hash)]
22pub struct BidirectedAdjacencyArray<IndexType: GraphIndexInteger, NodeData, EdgeData> {
23 pub(crate) node_array: TaggedVec<DirectedNodeIndex<IndexType>, DirectedEdgeIndex<IndexType>>,
30
31 pub(crate) edge_array: TaggedVec<DirectedEdgeIndex<IndexType>, DirectedNodeIndex<IndexType>>,
36
37 pub(crate) node_data: TaggedVec<NodeIndex<IndexType>, NodeData>,
43
44 pub(crate) edge_data_keys: TaggedVec<DirectedEdgeIndex<IndexType>, EdgeDataKey<IndexType>>,
51
52 pub(crate) edge_data: TaggedVec<EdgeIndex<IndexType>, BidirectedEdgeData<IndexType, EdgeData>>,
56}
57
58#[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash)]
59pub(crate) struct EdgeDataKey<IndexType: GraphIndexInteger> {
60 inverse: DirectedEdgeIndex<IndexType>,
61 data_index: OptionalEdgeIndex<IndexType>,
62}
63
64#[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Hash)]
65pub(crate) struct BidirectedEdgeData<IndexType, EdgeData> {
66 forward: DirectedEdgeIndex<IndexType>,
67 reverse: DirectedEdgeIndex<IndexType>,
68 data: EdgeData,
69}
70
71pub struct DirectedEdge<IndexType> {
72 from: DirectedNodeIndex<IndexType>,
73 to: DirectedNodeIndex<IndexType>,
74 index: DirectedEdgeIndex<IndexType>,
75}
76
77pub struct DirectedEdgeDataView<'a, IndexType, EdgeData> {
78 is_forward: bool,
79 edge: EdgeIndex<IndexType>,
80 data: &'a EdgeData,
81}
82
83pub struct EdgeView<'a, IndexType, EdgeData> {
84 from: DirectedNodeIndex<IndexType>,
85 to: DirectedNodeIndex<IndexType>,
86 forward: DirectedEdgeIndex<IndexType>,
87 reverse: DirectedEdgeIndex<IndexType>,
88 data: &'a EdgeData,
89}
90
91#[derive(Debug, Clone, Hash, PartialEq, Eq)]
92pub struct BidirectedEdge<IndexType, EdgeData> {
93 pub from: NodeIndex<IndexType>,
94 pub from_forward: bool,
96 pub to: NodeIndex<IndexType>,
97 pub to_forward: bool,
99 pub data: EdgeData,
100}
101
102impl<IndexType: GraphIndexInteger, NodeData, EdgeData>
103 BidirectedAdjacencyArray<IndexType, NodeData, EdgeData>
104{
105 pub fn new(
123 nodes: TaggedVec<NodeIndex<IndexType>, NodeData>,
124 edges: TaggedVec<EdgeIndex<IndexType>, BidirectedEdge<IndexType, EdgeData>>,
125 ) -> Self {
126 let mut node_array = TaggedVec::from_iter(iter::repeat_n(
127 DirectedEdgeIndex::from_usize(0),
128 nodes.len() * 2 + 1,
129 ));
130
131 for edge in edges.iter_values() {
133 let from_directed_forward =
134 DirectedNodeIndex::from_bidirected(edge.from, edge.from_forward);
135 node_array[from_directed_forward].increment();
136 let from_directed_reverse =
137 DirectedNodeIndex::from_bidirected(edge.to, edge.to_forward).invert();
138 node_array[from_directed_reverse].increment();
139 }
140
141 let directed_edge_count =
143 node_array
144 .iter_values_mut()
145 .fold(DirectedEdgeIndex::zero(), |sum, element| {
146 let sum = sum.add(*element);
147 *element = sum;
148 sum
149 });
150 assert_eq!(
151 directed_edge_count,
152 node_array.iter_values().last().copied().unwrap(),
153 );
154
155 let mut edge_array = TaggedVec::from_iter(iter::repeat_n(
157 DirectedNodeIndex::from_usize(0),
158 directed_edge_count.into_usize(),
159 ));
160 let mut edge_data_keys = TaggedVec::from_iter(iter::repeat_n(
161 EdgeDataKey {
162 inverse: DirectedEdgeIndex::zero(),
163 data_index: OptionalEdgeIndex::new_none(),
164 },
165 directed_edge_count.into_usize(),
166 ));
167 let mut edge_data = TaggedVec::new();
168
169 for (edge_index, edge) in edges.into_iter(..) {
172 let from_directed_forward =
173 DirectedNodeIndex::from_bidirected(edge.from, edge.from_forward);
174 let to_directed_forward = DirectedNodeIndex::from_bidirected(edge.to, edge.to_forward);
175 let edge_index_forward = {
176 node_array[from_directed_forward].decrement();
177 node_array[from_directed_forward]
178 };
179
180 let from_directed_reverse = to_directed_forward.invert();
181 let to_directed_reverse = from_directed_forward.invert();
182 let edge_index_reverse = {
183 node_array[from_directed_reverse].decrement();
184 node_array[from_directed_reverse]
185 };
186
187 edge_array[edge_index_forward] = to_directed_forward;
188 edge_array[edge_index_reverse] = to_directed_reverse;
189
190 edge_data_keys[edge_index_forward] = EdgeDataKey {
191 inverse: edge_index_reverse,
192 data_index: edge_index.into(),
193 };
194 edge_data_keys[edge_index_reverse] = EdgeDataKey {
195 inverse: edge_index_forward,
196 data_index: OptionalEdgeIndex::new_none(),
197 };
198
199 let data_index = edge_data.push(BidirectedEdgeData {
200 forward: edge_index_forward,
201 reverse: edge_index_reverse,
202 data: edge.data,
203 });
204 assert_eq!(edge_index, data_index);
205 }
206
207 Self {
208 node_array,
209 edge_array,
210 node_data: nodes,
211 edge_data_keys,
212 edge_data,
213 }
214 }
215
216 pub fn reorder_edges(
220 &mut self,
221 node: DirectedNodeIndex<IndexType>,
222 comparator: impl FnMut(
223 DirectedNodeIndex<IndexType>,
224 DirectedNodeIndex<IndexType>,
225 ) -> std::cmp::Ordering,
226 ) {
227 self.reorder_edges_with_buffers(
228 node,
229 comparator,
230 &mut Permutation::one(0),
231 &mut BitVec::new(),
232 );
233 }
234
235 fn reorder_edges_with_buffers(
236 &mut self,
237 node: DirectedNodeIndex<IndexType>,
238 mut comparator: impl FnMut(
239 DirectedNodeIndex<IndexType>,
240 DirectedNodeIndex<IndexType>,
241 ) -> std::cmp::Ordering,
242 permutation: &mut Permutation,
243 bitvec: &mut BitVec,
244 ) {
245 let offset = self.node_array[node].into_usize();
246 let limit = self.node_array[node.add(DirectedNodeIndex::from_usize(1))].into_usize();
247
248 permutation.assign_from_sort_by(
250 &self.edge_array.as_untagged_slice()[offset..limit],
251 |a, b| comparator(*a, *b),
252 );
253
254 permutation
256 .apply_slice_in_place(&mut self.edge_array.as_untagged_mut_slice()[offset..limit]);
257
258 bitvec.clear();
260 bitvec.resize(limit - offset, false);
261 for edge_index in offset..limit {
262 if bitvec[edge_index - offset] {
263 continue;
264 }
265
266 let permuted_edge_index =
267 DirectedEdgeIndex::from_usize(permutation.apply_idx(edge_index - offset) + offset);
268 let edge_index = DirectedEdgeIndex::from_usize(edge_index);
269 let inverse_edge_index = self.edge_data_keys[edge_index].inverse;
270
271 self.edge_data_keys[inverse_edge_index].inverse = permuted_edge_index;
272
273 if (offset..limit).contains(&inverse_edge_index.into_usize()) {
274 bitvec.set(inverse_edge_index.into_usize() - offset, true);
277
278 let permuted_inverse_edge_index = DirectedEdgeIndex::from_usize(
279 permutation.apply_idx(inverse_edge_index.into_usize() - offset) + offset,
280 );
281 self.edge_data_keys[edge_index].inverse = permuted_inverse_edge_index;
282 }
283 }
284
285 permutation
286 .apply_slice_in_place(&mut self.edge_data_keys.as_untagged_mut_slice()[offset..limit]);
287
288 bitvec.clear();
290 bitvec.resize(limit - offset, false);
291 for edge_index in offset..limit {
292 if bitvec[edge_index - offset] {
293 continue;
294 }
295
296 let edge_index = DirectedEdgeIndex::from_usize(edge_index);
297 let inverse_edge_index = self.edge_data_keys[edge_index].inverse;
298
299 let edge_data_key = self.edge_data_keys[edge_index]
300 .data_index
301 .into_option()
302 .xor(
303 self.edge_data_keys[inverse_edge_index]
304 .data_index
305 .into_option(),
306 )
307 .unwrap();
308 let edge_data = &mut self.edge_data[edge_data_key];
309
310 if (offset..limit).contains(&edge_data.forward.into_usize()) {
311 edge_data.forward = DirectedEdgeIndex::from_usize(
312 permutation.apply_idx(edge_data.forward.into_usize() - offset) + offset,
313 );
314 }
315
316 if (offset..limit).contains(&edge_data.reverse.into_usize()) {
317 edge_data.reverse = DirectedEdgeIndex::from_usize(
318 permutation.apply_idx(edge_data.reverse.into_usize() - offset) + offset,
319 );
320 }
321
322 if (offset..limit).contains(&inverse_edge_index.into_usize()) {
323 bitvec.set(inverse_edge_index.into_usize() - offset, true);
326 }
327 }
328 }
329
330 pub fn remove_multiedges(&mut self) -> BTreeSet<EdgeIndex<IndexType>> {
338 let mut deleted_edges = BTreeSet::default();
339 let mut existing_edges = FxHashSet::default();
340 let mut directed_edge_decrement = DirectedEdgeIndex::zero();
341
342 for from_node in self
343 .node_array
344 .iter_indices(..DirectedNodeIndex::from_usize(self.node_array.len() - 1))
345 {
346 if DEBUG_REMOVE_MULTIEDGES {
347 println!("from_node: {from_node}");
348 }
349
350 let edge_offset = self.node_array[from_node];
351 let edge_limit = self.node_array[from_node.add(1.into())];
352
353 self.node_array[from_node] = self.node_array[from_node].sub(directed_edge_decrement);
355
356 existing_edges.clear();
357 for directed_edge in self.edge_array.iter_indices(edge_offset..edge_limit) {
358 if DEBUG_REMOVE_MULTIEDGES {
359 println!(" directed_edge: {directed_edge}");
360 }
361
362 let to_node = self.edge_array[directed_edge];
363 if DEBUG_REMOVE_MULTIEDGES {
364 println!(" to_node: {to_node}");
365 }
366 let decremented_directed_edge = directed_edge.sub(directed_edge_decrement);
367 if DEBUG_REMOVE_MULTIEDGES {
368 println!(" decremented_directed_edge: {decremented_directed_edge}");
369 }
370 let inverse_directed_edge = self.edge_data_keys[directed_edge].inverse;
371 if DEBUG_REMOVE_MULTIEDGES {
372 println!(" inverse_directed_edge: {inverse_directed_edge}");
373 }
374
375 if existing_edges.insert(to_node) {
376 self.edge_array[decremented_directed_edge] = self.edge_array[directed_edge];
380
381 self.edge_data_keys[decremented_directed_edge] =
382 self.edge_data_keys[directed_edge];
383 self.edge_data_keys[inverse_directed_edge].inverse = decremented_directed_edge;
384
385 let edge = self.edge_data_keys[directed_edge]
386 .data_index
387 .into_option()
388 .or(self.edge_data_keys[inverse_directed_edge]
389 .data_index
390 .into_option())
391 .unwrap();
392 if DEBUG_REMOVE_MULTIEDGES {
393 println!(" edge: {edge}");
394 }
395
396 if self.edge_data[edge].forward == directed_edge {
397 self.edge_data[edge].forward = decremented_directed_edge;
398 } else if self.edge_data[edge].reverse == directed_edge {
399 self.edge_data[edge].reverse = decremented_directed_edge;
400 }
401 } else {
402 if DEBUG_REMOVE_MULTIEDGES {
404 println!(" DELETE");
405 }
406 if let Some(edge) = self.edge_data_keys[directed_edge].data_index.into_option()
407 {
408 if DEBUG_REMOVE_MULTIEDGES {
409 println!(" edge: {edge}");
410 }
411 deleted_edges.insert(edge);
412 }
413 directed_edge_decrement.increment();
414 }
415 }
416 }
417
418 let deleted_edges = deleted_edges;
419
420 self.edge_array.splice(
422 DirectedEdgeIndex::from_usize(self.edge_array.len()).sub(directed_edge_decrement)..,
423 iter::empty(),
424 );
425 self.edge_data_keys.splice(
426 DirectedEdgeIndex::from_usize(self.edge_data_keys.len()).sub(directed_edge_decrement)..,
427 iter::empty(),
428 );
429
430 let mut deleted_edges_iter = deleted_edges.iter().copied().peekable();
431 let mut edge_decrement = EdgeIndex::from_usize(0);
432 for edge in self.edge_data.iter_indices(..) {
433 if Some(&edge) == deleted_edges_iter.peek() {
434 deleted_edges_iter.next();
436 edge_decrement.increment();
437 } else {
438 let decremented_edge = edge.sub(edge_decrement);
440
441 if self.edge_data_keys[self.edge_data[edge].forward]
442 .data_index
443 .is_some()
444 {
445 self.edge_data_keys[self.edge_data[edge].forward].data_index =
446 decremented_edge.into();
447 } else if self.edge_data_keys[self.edge_data[edge].reverse]
448 .data_index
449 .is_some()
450 {
451 self.edge_data_keys[self.edge_data[edge].reverse].data_index =
452 decremented_edge.into();
453 }
454 }
455 }
456
457 let mut deleted_edges_iter = deleted_edges.iter().copied().peekable();
459 let mut current_index = EdgeIndex::from_usize(0);
460 self.edge_data.retain(|_| {
461 if Some(¤t_index) == deleted_edges_iter.peek() {
462 deleted_edges_iter.next();
464 current_index.increment();
465 false
466 } else {
467 current_index.increment();
469 true
470 }
471 });
472
473 deleted_edges
474 }
475
476 pub fn node_count(&self) -> usize {
477 self.node_data.len()
478 }
479
480 pub fn edge_count(&self) -> usize {
481 self.edge_data.len()
482 }
483
484 pub fn iter_nodes(&self) -> impl Iterator<Item = NodeIndex<IndexType>> {
485 self.node_data.iter_indices(..)
486 }
487
488 pub fn iter_edges(&self) -> impl Iterator<Item = EdgeIndex<IndexType>> {
489 self.edge_data.iter_indices(..)
490 }
491
492 pub fn iter_outgoing_edges(
493 &self,
494 node: DirectedNodeIndex<IndexType>,
495 ) -> impl Iterator<Item = DirectedEdge<IndexType>> {
496 let start = self.node_array[node];
497 let end = self.node_array[node.add(DirectedNodeIndex::from_usize(1))];
498 self.edge_array
499 .iter(start..end)
500 .map(move |(edge_index, &to_node)| DirectedEdge {
501 from: node,
502 to: to_node,
503 index: edge_index,
504 })
505 }
506
507 pub fn iter_incident_edges(
509 &self,
510 node: NodeIndex<IndexType>,
511 ) -> impl Iterator<Item = EdgeIndex<IndexType>> {
512 let forward_node = DirectedNodeIndex::from_bidirected(node, true);
513 let reverse_node = DirectedNodeIndex::from_bidirected(node, false);
514 self.iter_outgoing_edges(forward_node)
515 .chain(self.iter_outgoing_edges(reverse_node))
516 .filter_map(|directed_edge| {
517 let directed_edge_data = self.directed_edge_data(directed_edge.index());
518 if directed_edge.from() == directed_edge.to()
519 || directed_edge.from() == directed_edge.to().invert()
520 {
521 directed_edge_data
522 .is_forward()
523 .then_some(directed_edge_data.edge())
524 } else {
525 Some(directed_edge_data.edge())
526 }
527 })
528 }
529
530 pub fn node_data(&self, node: NodeIndex<IndexType>) -> &NodeData {
531 &self.node_data[node]
532 }
533
534 pub fn edge(&self, edge: EdgeIndex<IndexType>) -> EdgeView<'_, IndexType, EdgeData> {
535 let bidirected_edge_data = &self.edge_data[edge];
536
537 let forward_to = self.edge_array[bidirected_edge_data.forward];
538 let reverse_to = self.edge_array[bidirected_edge_data.reverse];
539
540 if forward_to == reverse_to {
541 let from = forward_to.invert();
543 let to = forward_to;
544 EdgeView {
545 from,
546 to,
547 forward: bidirected_edge_data.forward,
548 reverse: bidirected_edge_data.reverse,
549 data: &bidirected_edge_data.data,
550 }
551 } else if forward_to.invert() == reverse_to {
552 let from = forward_to;
554 let to = forward_to;
555 EdgeView {
556 from,
557 to,
558 forward: bidirected_edge_data.forward,
559 reverse: bidirected_edge_data.reverse,
560 data: &bidirected_edge_data.data,
561 }
562 } else {
563 let from = reverse_to.invert();
565 let to = forward_to;
566 EdgeView {
567 from,
568 to,
569 forward: bidirected_edge_data.forward,
570 reverse: bidirected_edge_data.reverse,
571 data: &bidirected_edge_data.data,
572 }
573 }
574 }
575
576 pub fn directed_edge_data<'this>(
577 &'this self,
578 directed_edge: DirectedEdgeIndex<IndexType>,
579 ) -> DirectedEdgeDataView<'this, IndexType, EdgeData> {
580 let key = &self.edge_data_keys[directed_edge];
581 if let Some(edge) = key.data_index.into_option() {
582 DirectedEdgeDataView {
583 is_forward: true,
584 edge,
585 data: &self.edge_data[edge].data,
586 }
587 } else {
588 let inverse_key = &self.edge_data_keys[key.inverse];
589 let Some(edge) = inverse_key.data_index.into_option() else {
590 panic!(
591 "Edge data for edge {:?} and its inverse {:?} are both missing",
592 directed_edge, key.inverse
593 );
594 };
595 DirectedEdgeDataView {
596 is_forward: false,
597 edge,
598 data: &self.edge_data[edge].data,
599 }
600 }
601 }
602
603 pub fn directed_edge_into_bidirected(
604 &self,
605 directed_edge: DirectedEdgeIndex<IndexType>,
606 ) -> EdgeIndex<IndexType> {
607 let key = &self.edge_data_keys[directed_edge];
608 if let Some(edge) = key.data_index.into_option() {
609 edge
610 } else {
611 let inverse_key = &self.edge_data_keys[key.inverse];
612 inverse_key
613 .data_index
614 .expect("Edge data for directed edge and its inverse are both missing")
615 }
616 }
617}
618
619impl<IndexType> DirectedEdge<IndexType> {
620 pub fn from(&self) -> DirectedNodeIndex<IndexType>
621 where
622 IndexType: Copy,
623 {
624 self.from
625 }
626
627 pub fn to(&self) -> DirectedNodeIndex<IndexType>
628 where
629 IndexType: Copy,
630 {
631 self.to
632 }
633
634 pub fn index(&self) -> DirectedEdgeIndex<IndexType>
635 where
636 IndexType: Copy,
637 {
638 self.index
639 }
640}
641
642impl<'a, IndexType, EdgeData> DirectedEdgeDataView<'a, IndexType, EdgeData> {
643 pub fn is_forward(&self) -> bool {
644 self.is_forward
645 }
646
647 pub fn edge(&self) -> EdgeIndex<IndexType>
648 where
649 IndexType: Copy,
650 {
651 self.edge
652 }
653
654 pub fn data(&self) -> &EdgeData {
655 self.data
656 }
657}
658
659impl<'a, IndexType, EdgeData> EdgeView<'a, IndexType, EdgeData> {
660 pub fn from(&self) -> DirectedNodeIndex<IndexType>
661 where
662 IndexType: Copy,
663 {
664 self.from
665 }
666
667 pub fn to(&self) -> DirectedNodeIndex<IndexType>
668 where
669 IndexType: Copy,
670 {
671 self.to
672 }
673
674 pub fn forward(&self) -> DirectedEdgeIndex<IndexType>
675 where
676 IndexType: Copy,
677 {
678 self.forward
679 }
680
681 pub fn reverse(&self) -> DirectedEdgeIndex<IndexType>
682 where
683 IndexType: Copy,
684 {
685 self.reverse
686 }
687
688 pub fn data(&self) -> &EdgeData {
689 self.data
690 }
691}
692
693impl<IndexType: GraphIndexInteger, EdgeData> BidirectedEdge<IndexType, EdgeData> {
694 pub fn new(
695 from: DirectedNodeIndex<IndexType>,
696 to: DirectedNodeIndex<IndexType>,
697 data: EdgeData,
698 ) -> Self {
699 Self {
700 from: from.into_bidirected(),
701 from_forward: from.is_forward(),
702 to: to.into_bidirected(),
703 to_forward: to.is_forward(),
704 data,
705 }
706 }
707}
708
709impl<IndexType: GraphIndexInteger> BidirectedEdge<IndexType, PlainGfaEdgeData> {
710 pub fn new_gfa(
711 from: DirectedNodeIndex<IndexType>,
712 to: DirectedNodeIndex<IndexType>,
713 overlap: u16,
714 ) -> Self {
715 Self {
716 from: from.into_bidirected(),
717 from_forward: from.is_forward(),
718 to: to.into_bidirected(),
719 to_forward: to.is_forward(),
720 data: PlainGfaEdgeData::new(overlap),
721 }
722 }
723}