1use std::{
2 io::{Write, stderr, stdout},
3 sync::Arc,
4};
5
6use algebra::module::{
7 Module,
8 homomorphism::{FreeModuleHomomorphism, ModuleHomomorphism},
9};
10use ext::{
11 chain_complex::{AugmentedChainComplex, BoundedChainComplex, ChainComplex, FreeChainComplex},
12 utils,
13 yoneda::yoneda_representative_element,
14};
15use fp::{matrix::Matrix, vector::FpVector};
16use itertools::Itertools;
17use sseq::coordinates::{Bidegree, BidegreeElement};
18use tensor_product_chain_complex::TensorChainComplex;
19
20fn main() -> anyhow::Result<()> {
21 ext::utils::init_logging()?;
22
23 let resolution = Arc::new(utils::query_module_only("Module", None, false)?);
24 let module = resolution.target().module(0);
25 let p = resolution.prime();
26
27 if resolution.target().max_s() != 1 || !module.is_unit() || p != 2 {
28 panic!("Can only run Steenrod on the sphere");
29 }
30
31 let b = Bidegree::n_s(
32 query::raw("n of Ext class", str::parse),
33 query::raw("s of Ext class", str::parse),
34 );
35
36 resolution.compute_through_bidegree(b + b);
37
38 let class: Vec<u32> =
39 query::vector("Input Ext class", resolution.number_of_gens_in_bidegree(b));
40
41 let yoneda = Arc::new(yoneda_representative_element(
42 Arc::clone(&resolution),
43 b,
44 &class,
45 ));
46
47 print!("Dimensions of Yoneda representative: 1");
48 for s in 0..=b.s() {
49 print!(" {}", yoneda.module(s).total_dimension());
50 }
51 println!();
52
53 let square = Arc::new(TensorChainComplex::new(
54 Arc::clone(&yoneda),
55 Arc::clone(&yoneda),
56 ));
57 let doubled_b = b + b;
58
59 tracing::info_span!("Computing quasi-inverses").in_scope(|| {
60 square.compute_through_bidegree(doubled_b);
61 for s in 0..=doubled_b.s() {
62 square
63 .differential(s)
64 .compute_auxiliary_data_through_degree(doubled_b.t());
65 }
66 });
67
68 eprintln!("Computing Steenrod operations: ");
69
70 let mut delta = Vec::with_capacity(b.s() as usize);
71
72 for i in 0..=b.s() {
73 let mut maps: Vec<Arc<FreeModuleHomomorphism<_>>> =
74 Vec::with_capacity(doubled_b.s() as usize - 1);
75
76 for s in 0..=doubled_b.s() - i {
77 let source = resolution.module(s);
78 let target = square.module(s + i);
79
80 let map = FreeModuleHomomorphism::new(Arc::clone(&source), Arc::clone(&target), 0);
81 maps.push(Arc::new(map));
82 }
83 delta.push(maps);
84 }
85
86 let tracing_span = tracing::info_span!("Computing Steenrod operations");
87 let _tracing_guard = tracing_span.enter();
88
89 for i in 0..=b.s() {
91 let shift_s = Bidegree::s_t(i, 0);
92 let start = std::time::Instant::now();
95
96 for s in 0..=(doubled_b - shift_s).s() {
97 if i == 0 && s == 0 {
98 let map = &delta[0][0];
99 map.add_generators_from_matrix_rows(
100 0,
101 Matrix::from_vec(p, &[vec![1]]).as_slice_mut(),
102 );
103 map.extend_by_zero(doubled_b.t());
104 continue;
105 }
106
107 let square = Arc::clone(&square);
108
109 let source = resolution.module(s);
110 let target = square.module(s + i);
111
112 let dtarget_module = square.module(s + i - 1);
113
114 let d_res = resolution.differential(s);
115 let d_target = square.differential(s + i);
116
117 let map = Arc::clone(&delta[i as usize][s as usize]);
118 let prev_map = match s {
119 0 => None,
120 _ => Some(Arc::clone(&delta[i as usize][s as usize - 1])),
121 };
122
123 let prev_delta = match i {
124 0 => None,
125 _ => Some(Arc::clone(&delta[i as usize - 1][s as usize])),
126 };
127
128 for t in 0..=doubled_b.t() {
129 let b = Bidegree::s_t(s, t);
130
131 let num_gens = source.number_of_gens_in_degree(t);
132
133 let fx_dim = target.dimension(t);
134 let fdx_dim = dtarget_module.dimension(t);
135
136 if fx_dim == 0 || fdx_dim == 0 || num_gens == 0 {
137 map.extend_by_zero(t);
138 continue;
139 }
140
141 let mut output_matrix = Matrix::new(p, num_gens, fx_dim);
142 let mut result = FpVector::new(p, fdx_dim);
143 for j in 0..num_gens {
144 if let Some(m) = &prev_delta {
145 let prevd = m.output(t, j);
147
148 square.swap(&mut result, prevd, b + shift_s - Bidegree::s_t(1, 0));
150 result += prevd;
151 }
152
153 if let Some(m) = &prev_map {
154 let dx = d_res.output(t, j);
155 m.apply(result.as_slice_mut(), 1, t, dx.as_slice());
156 }
157 assert!(d_target.apply_quasi_inverse(
158 output_matrix.row_mut(j),
159 t,
160 result.as_slice(),
161 ));
162
163 result.set_to_zero();
164 }
165 map.add_generators_from_matrix_rows(t, output_matrix.as_slice_mut());
166 }
167 }
168
169 let final_map = &delta[i as usize][(doubled_b - shift_s).s() as usize];
170 let num_gens = resolution.number_of_gens_in_bidegree(doubled_b - shift_s);
171 print!(
172 "Sq^{} {basis_string} = [{result}]",
173 (b - shift_s).s(),
174 basis_string =
175 BidegreeElement::new(b, FpVector::from_slice(p, &class)).to_basis_string(),
176 result = (0..num_gens)
177 .map(|k| format!("{}", final_map.output(doubled_b.t(), k).entry(0)))
178 .format(", "),
179 );
180 stdout().flush().unwrap();
181 eprint!(" ({:?})", start.elapsed());
182 stderr().flush().unwrap();
183 println!();
184 }
185
186 Ok(())
187}
188
189mod sum_module {
190 use std::sync::Arc;
191
192 use algebra::module::{
193 Module, ZeroModule,
194 block_structure::{BlockStructure, GeneratorBasisEltPair},
195 };
196 use bivec::BiVec;
197 use fp::vector::FpSliceMut;
198 use once::OnceBiVec;
199
200 pub struct SumModule<M: Module> {
201 algebra: Arc<M::Algebra>,
203 min_degree: i32,
204 pub modules: Vec<Arc<M>>,
205 pub block_structures: OnceBiVec<BlockStructure>,
206 }
207
208 impl<M: Module> std::fmt::Display for SumModule<M> {
209 fn fmt(&self, f: &mut std::fmt::Formatter) -> std::fmt::Result {
210 if self.modules.is_empty() {
211 write!(f, "0")
212 } else {
213 write!(f, "{}", self.modules[0])?;
214 for m in self.modules.iter().skip(1) {
215 write!(f, " (+) {m}")?;
216 }
217 Ok(())
218 }
219 }
220 }
221
222 impl<M: Module> SumModule<M> {
223 pub fn new(algebra: Arc<M::Algebra>, modules: Vec<Arc<M>>, min_degree: i32) -> Self {
224 Self {
225 algebra,
226 modules,
227 min_degree,
228 block_structures: OnceBiVec::new(min_degree),
229 }
230 }
231
232 pub fn get_module_num(&self, degree: i32, index: usize) -> usize {
233 self.block_structures[degree]
234 .index_to_generator_basis_elt(index)
235 .generator_index
236 }
237
238 pub fn offset(&self, degree: i32, module_num: usize) -> usize {
239 self.block_structures[degree]
240 .generator_to_block(degree, module_num)
241 .start
242 }
243 }
244
245 impl<M: Module> Module for SumModule<M> {
246 type Algebra = M::Algebra;
247
248 fn algebra(&self) -> Arc<Self::Algebra> {
249 Arc::clone(&self.algebra)
250 }
251
252 fn min_degree(&self) -> i32 {
253 self.min_degree
254 }
255
256 fn compute_basis(&self, degree: i32) {
257 for module in &self.modules {
258 module.compute_basis(degree);
259 }
260 for i in self.block_structures.len()..=degree {
261 let mut block_sizes = BiVec::new(i);
262 block_sizes.push(self.modules.iter().map(|m| m.dimension(i)).collect());
263 self.block_structures
264 .push(BlockStructure::new(&block_sizes));
265 }
266 }
267
268 fn max_computed_degree(&self) -> i32 {
269 self.block_structures.len()
270 }
271
272 fn dimension(&self, degree: i32) -> usize {
273 self.block_structures
274 .get(degree)
275 .map(BlockStructure::total_dimension)
276 .unwrap_or(0)
277 }
278
279 fn act_on_basis(
280 &self,
281 mut result: FpSliceMut,
282 coeff: u32,
283 op_degree: i32,
284 op_index: usize,
285 mod_degree: i32,
286 mod_index: usize,
287 ) {
288 let target_degree = mod_degree + op_degree;
289 let GeneratorBasisEltPair {
290 generator_index: module_num,
291 basis_index,
292 ..
293 } = self.block_structures[mod_degree].index_to_generator_basis_elt(mod_index);
294 let range =
295 self.block_structures[target_degree].generator_to_block(target_degree, *module_num);
296 let module = &self.modules[*module_num];
297
298 module.act_on_basis(
299 result.slice_mut(range.start, range.end),
300 coeff,
301 op_degree,
302 op_index,
303 mod_degree,
304 *basis_index,
305 );
306 }
307
308 fn basis_element_to_string(&self, degree: i32, index: usize) -> String {
309 let GeneratorBasisEltPair {
310 generator_index: module_num,
311 basis_index,
312 ..
313 } = self.block_structures[degree].index_to_generator_basis_elt(index);
314 self.modules[*module_num].basis_element_to_string(degree, *basis_index)
315 }
316
317 fn max_degree(&self) -> Option<i32> {
318 self.modules
319 .iter()
320 .map(|m| m.max_degree())
321 .max()
322 .unwrap_or(Some(self.min_degree))
323 }
324 }
325
326 impl<M: Module> ZeroModule for SumModule<M> {
327 fn zero_module(algebra: Arc<M::Algebra>, min_degree: i32) -> Self {
328 Self::new(algebra, vec![], min_degree)
329 }
330 }
331
332 #[cfg(test)]
333 mod tests {
334 use algebra::{AdemAlgebra, module::FDModule};
335
336 use super::*;
337
338 #[test]
339 fn test_sum_modules() {
340 let k = r#"{"type" : "finite dimensional module", "p": 2, "gens": {"x0": 0}, "actions": []}"#;
341 let k2 = r#"{"type" : "finite dimensional module", "p": 2, "gens": {"x0": 0, "y0":0}, "actions": []}"#;
342 let zero =
343 r#"{"type" : "finite dimensional module", "p": 2, "gens": {}, "actions": []}"#;
344 let c2 = r#"{"type" : "finite dimensional module", "p": 2, "gens": {"x0": 0, "x1": 1}, "actions": ["Sq1 x0 = x1"]}"#;
345 let ceta = r#"{"type" : "finite dimensional module", "p": 2, "gens": {"x0": 0, "x2": 2}, "actions": ["Sq2 x0 = x2"]}"#;
346 let c2sumceta = r#"{"type" : "finite dimensional module", "p": 2, "gens": {"x0": 0, "x1": 1,"y0": 0, "y2": 2}, "actions": ["Sq1 x0 = x1", "Sq2 y0 = y2"]}"#;
347
348 test_sum_module(vec![], zero);
349 test_sum_module(vec![k, k], k2);
350 test_sum_module(vec![c2, ceta], c2sumceta);
351 }
352
353 fn test_sum_module(modules: Vec<&str>, desc: &str) {
354 let p = fp::prime::ValidPrime::new(2);
355 let algebra = Arc::new(AdemAlgebra::new(p, false));
356
357 let modules: Vec<Arc<FDModule<AdemAlgebra>>> = modules
358 .into_iter()
359 .map(|s| {
360 let m = serde_json::from_str(s).unwrap();
361 Arc::new(FDModule::from_json(Arc::clone(&algebra), &m).unwrap())
362 })
363 .collect::<Vec<_>>();
364
365 let sum_1 = FDModule::from(&SumModule::new(Arc::clone(&algebra), modules, 0));
366
367 let desc = serde_json::from_str(desc).unwrap();
368 let sum_2 = FDModule::from_json(Arc::clone(&algebra), &desc).unwrap();
369
370 if let Err(msg) = sum_1.test_equal(&sum_2) {
371 panic!("Test case failed. {msg}");
372 }
373 }
374 }
375}
376
377mod tensor_product_chain_complex {
378 use std::sync::Arc;
379
380 use algebra::{
381 Algebra, Bialgebra,
382 module::{Module, TensorModule, ZeroModule, homomorphism::ModuleHomomorphism},
383 };
384 use ext::chain_complex::ChainComplex;
385 use fp::{
386 matrix::AugmentedMatrix,
387 vector::{FpSlice, FpSliceMut, FpVector},
388 };
389 use once::OnceBiVec;
390 use sseq::coordinates::Bidegree;
391
392 use super::sum_module::SumModule;
393
394 pub type Stm<M, N> = SumModule<TensorModule<M, N>>;
395
396 pub struct TensorChainComplex<A, CC1, CC2>
397 where
398 A: Algebra + Bialgebra,
399 CC1: ChainComplex<Algebra = A>,
400 CC2: ChainComplex<Algebra = A>,
401 {
402 left_cc: Arc<CC1>,
403 right_cc: Arc<CC2>,
404 modules: OnceBiVec<Arc<Stm<CC1::Module, CC2::Module>>>,
405 zero_module: Arc<Stm<CC1::Module, CC2::Module>>,
406 differentials: OnceBiVec<Arc<TensorChainMap<A, CC1, CC2>>>,
407 }
408
409 impl<A, CC1, CC2> TensorChainComplex<A, CC1, CC2>
410 where
411 A: Algebra + Bialgebra,
412 CC1: ChainComplex<Algebra = A>,
413 CC2: ChainComplex<Algebra = A>,
414 {
415 pub fn new(left_cc: Arc<CC1>, right_cc: Arc<CC2>) -> Self {
416 Self {
417 modules: OnceBiVec::new(0),
418 differentials: OnceBiVec::new(0),
419 zero_module: Arc::new(SumModule::zero_module(
420 left_cc.algebra(),
421 left_cc.min_degree() + right_cc.min_degree(),
422 )),
423 left_cc,
424 right_cc,
425 }
426 }
427
428 fn left_cc(&self) -> Arc<CC1> {
429 Arc::clone(&self.left_cc)
430 }
431
432 fn right_cc(&self) -> Arc<CC2> {
433 Arc::clone(&self.right_cc)
434 }
435
436 fn left_min_shift(&self) -> Bidegree {
437 Bidegree::s_t(0, self.left_cc.min_degree())
438 }
439
440 fn right_min_shift(&self) -> Bidegree {
441 Bidegree::s_t(0, self.right_cc.min_degree())
442 }
443 }
444
445 impl<A, CC> TensorChainComplex<A, CC, CC>
446 where
447 A: Algebra + Bialgebra,
448 CC: ChainComplex<Algebra = A>,
449 {
450 pub fn swap(&self, result: &mut FpVector, vec: &FpVector, b: Bidegree) {
453 let s = b.s();
454 for left_s in 0..=s {
455 let right_s = s - left_s;
456 let module = &self.modules[s];
457
458 let source_offset = module.offset(b.t(), left_s as usize);
459 let target_offset = module.offset(b.t(), right_s as usize);
460
461 for left_t in 0..=b.t() {
462 let right_t = b.t() - left_t;
463
464 let left_dim = module.modules[left_s as usize].left.dimension(left_t);
465 let right_dim = module.modules[left_s as usize].right.dimension(right_t);
466
467 if left_dim == 0 || right_dim == 0 {
468 continue;
469 }
470
471 let source_inner_offset = module.modules[left_s as usize].offset(b.t(), left_t);
472 let target_inner_offset =
473 module.modules[right_s as usize].offset(b.t(), right_t);
474
475 for i in 0..left_dim {
476 for j in 0..right_dim {
477 let value =
478 vec.entry(source_offset + source_inner_offset + i * right_dim + j);
479 if value != 0 {
480 result.add_basis_element(
481 target_offset + target_inner_offset + j * left_dim + i,
482 value,
483 );
484 }
485 }
486 }
487 }
488 }
489 }
490 }
491
492 impl<A, CC1, CC2> ChainComplex for TensorChainComplex<A, CC1, CC2>
493 where
494 A: Algebra + Bialgebra,
495 CC1: ChainComplex<Algebra = A>,
496 CC2: ChainComplex<Algebra = A>,
497 {
498 type Algebra = A;
499 type Homomorphism = TensorChainMap<A, CC1, CC2>;
500 type Module = Stm<CC1::Module, CC2::Module>;
501
502 fn algebra(&self) -> Arc<A> {
503 self.left_cc.algebra()
504 }
505
506 fn min_degree(&self) -> i32 {
507 self.left_cc.min_degree() + self.right_cc.min_degree()
508 }
509
510 fn zero_module(&self) -> Arc<Self::Module> {
511 Arc::clone(&self.zero_module)
512 }
513
514 fn has_computed_bidegree(&self, b: Bidegree) -> bool {
515 self.left_cc
516 .has_computed_bidegree(b - self.left_min_shift())
517 && self
518 .right_cc
519 .has_computed_bidegree(b - self.left_min_shift())
520 && self.differentials.len() > b.s()
521 }
522
523 fn module(&self, s: i32) -> Arc<Self::Module> {
524 Arc::clone(&self.modules[s])
525 }
526
527 fn differential(&self, s: i32) -> Arc<Self::Homomorphism> {
528 Arc::clone(&self.differentials[s])
529 }
530
531 fn compute_through_bidegree(&self, b: Bidegree) {
532 self.left_cc
533 .compute_through_bidegree(b - self.right_min_shift());
534 self.right_cc
535 .compute_through_bidegree(b - self.left_min_shift());
536
537 self.modules.extend(b.s(), |i| {
538 let new_module_list: Vec<Arc<TensorModule<CC1::Module, CC2::Module>>> = (0..=i)
539 .map(|j| {
540 Arc::new(TensorModule::new(
541 self.left_cc.module(j),
542 self.right_cc.module(i - j),
543 ))
544 })
545 .collect::<Vec<_>>();
546 Arc::new(SumModule::new(
547 self.algebra(),
548 new_module_list,
549 self.min_degree(),
550 ))
551 });
552
553 for module in self.modules.values() {
554 module.compute_basis(b.t());
555 }
556
557 self.differentials.extend(b.s(), |s| {
558 if s == 0 {
559 Arc::new(TensorChainMap {
560 left_cc: self.left_cc(),
561 right_cc: self.right_cc(),
562 source_s: 0,
563 source: self.module(0),
564 target: self.zero_module(),
565 quasi_inverses: OnceBiVec::new(self.min_degree()),
566 })
567 } else {
568 Arc::new(TensorChainMap {
569 left_cc: self.left_cc(),
570 right_cc: self.right_cc(),
571 source_s: s,
572 source: self.module(s),
573 target: self.module(s - 1),
574 quasi_inverses: OnceBiVec::new(self.min_degree()),
575 })
576 }
577 });
578 }
579
580 fn next_homological_degree(&self) -> i32 {
581 self.modules.len()
582 }
583 }
584
585 pub struct TensorChainMap<A, CC1, CC2>
586 where
587 A: Algebra + Bialgebra,
588 CC1: ChainComplex<Algebra = A>,
589 CC2: ChainComplex<Algebra = A>,
590 {
591 left_cc: Arc<CC1>,
592 right_cc: Arc<CC2>,
593 source_s: i32,
594 source: Arc<Stm<CC1::Module, CC2::Module>>,
595 target: Arc<Stm<CC1::Module, CC2::Module>>,
596 quasi_inverses: OnceBiVec<Vec<Option<Vec<(usize, usize, FpVector)>>>>,
597 }
598
599 impl<A, CC1, CC2> ModuleHomomorphism for TensorChainMap<A, CC1, CC2>
600 where
601 A: Algebra + Bialgebra,
602 CC1: ChainComplex<Algebra = A>,
603 CC2: ChainComplex<Algebra = A>,
604 {
605 type Source = Stm<CC1::Module, CC2::Module>;
606 type Target = Stm<CC1::Module, CC2::Module>;
607
608 fn source(&self) -> Arc<Self::Source> {
609 Arc::clone(&self.source)
610 }
611
612 fn target(&self) -> Arc<Self::Target> {
613 Arc::clone(&self.target)
614 }
615
616 fn degree_shift(&self) -> i32 {
617 0
618 }
619
620 fn apply_to_basis_element(
622 &self,
623 mut result: FpSliceMut,
624 coeff: u32,
625 degree: i32,
626 input_idx: usize,
627 ) {
628 let left_s = self.source.get_module_num(degree, input_idx);
631 let right_s = self.source_s as usize - left_s;
632
633 let source_module = &self.source.modules[left_s];
634
635 let first_offset = self.source.offset(degree, left_s);
636 let inner_index = input_idx - first_offset;
637
638 let left_t = source_module.seek_module_num(degree, inner_index);
641 let right_t = degree - left_t;
642
643 let inner_index = inner_index - source_module.offset(degree, left_t);
644
645 let source_right_dim = source_module.right.dimension(right_t);
646 let right_index = inner_index % source_right_dim;
647 let left_index = inner_index / source_right_dim;
648
649 if right_s > 0 {
651 let target_module = &self.target.modules[left_s];
652 let target_offset = self.target.offset(degree, left_s)
653 + self.target.modules[left_s].offset(degree, left_t);
654 let target_right_dim = target_module.right.dimension(right_t);
655
656 let result = result.slice_mut(
657 target_offset + left_index * target_right_dim,
658 target_offset + (left_index + 1) * target_right_dim,
659 );
660 self.right_cc
661 .differential(right_s as i32)
662 .apply_to_basis_element(result, coeff, right_t, right_index);
663 }
664
665 if left_s > 0 {
667 let target_module = &self.target.modules[left_s - 1];
668 let target_offset = self.target.offset(degree, left_s - 1)
669 + self.target.modules[left_s - 1].offset(degree, left_t);
670 let target_right_dim = target_module.right.dimension(right_t);
671
672 let mut dl = FpVector::new(self.prime(), target_module.left.dimension(left_t));
673 self.left_cc
674 .differential(left_s as i32)
675 .apply_to_basis_element(dl.as_slice_mut(), coeff, left_t, left_index);
676 for i in 0..dl.len() {
677 result.add_basis_element(
678 target_offset + i * target_right_dim + right_index,
679 dl.entry(i),
680 );
681 }
682 }
683 }
684
685 fn compute_auxiliary_data_through_degree(&self, degree: i32) {
686 self.quasi_inverses
687 .extend(degree, |i| self.calculate_quasi_inverse(i));
688 }
689
690 fn apply_quasi_inverse(&self, mut result: FpSliceMut, degree: i32, input: FpSlice) -> bool {
691 let qis = &self.quasi_inverses[degree];
692 assert_eq!(input.len(), qis.len());
693
694 for (i, x) in input.iter_nonzero() {
695 if let Some(qi) = &qis[i] {
696 for (offset_start, offset_end, data) in qi.iter() {
697 result
698 .slice_mut(*offset_start, *offset_end)
699 .add(data.as_slice(), x);
700 }
701 }
702 }
703 true
704 }
705 }
706
707 impl<A, CC1, CC2> TensorChainMap<A, CC1, CC2>
708 where
709 A: Algebra + Bialgebra,
710 CC1: ChainComplex<Algebra = A>,
711 CC2: ChainComplex<Algebra = A>,
712 {
713 #[allow(clippy::range_minus_one)]
714 fn calculate_quasi_inverse(
715 &self,
716 degree: i32,
717 ) -> Vec<Option<Vec<(usize, usize, FpVector)>>> {
718 let p = self.prime();
719 let mut quasi_inverse_list: Vec<Option<Vec<(usize, usize, FpVector)>>> =
721 vec![None; self.target.dimension(degree)];
722
723 for left_t in self.left_cc.min_degree()..=degree - self.right_cc.min_degree() {
724 let right_t = degree - left_t;
725
726 let source_dim = self
727 .source
728 .modules
729 .iter()
730 .map(|m| m.left.dimension(left_t) * m.right.dimension(right_t))
731 .sum();
732 let target_dim = self
733 .target
734 .modules
735 .iter()
736 .map(|m| m.left.dimension(left_t) * m.right.dimension(right_t))
737 .sum();
738
739 if source_dim == 0 || target_dim == 0 {
740 continue;
741 }
742
743 let mut matrix = AugmentedMatrix::new(p, source_dim, [target_dim, source_dim]);
744
745 let mut target_offset = 0;
747 let mut row_count = 0;
748 for s in 0..=self.source_s - 1 {
749 let source_module = &self.source.modules[s as usize]; let target_module = &self.target.modules[s as usize]; let source_right_dim = source_module.right.dimension(right_t);
753 let source_left_dim = source_module.left.dimension(left_t);
754 let target_right_dim = target_module.right.dimension(right_t);
755 let target_left_dim = target_module.left.dimension(left_t);
756 assert_eq!(target_left_dim, source_left_dim);
757
758 let mut result = FpVector::new(p, target_right_dim);
759 for ri in 0..source_right_dim {
760 self.right_cc
761 .differential(self.source_s - s)
762 .apply_to_basis_element(result.as_slice_mut(), 1, right_t, ri);
763 for li in 0..source_left_dim {
764 let mut row = matrix.row_mut(row_count + li * source_right_dim + ri);
765 row.slice_mut(
766 target_offset + li * target_right_dim,
767 target_offset + (li + 1) * target_right_dim,
768 )
769 .assign(result.as_slice());
770 }
771 result.set_to_zero();
772 }
773 target_offset += target_right_dim * target_left_dim;
774 row_count += source_right_dim * source_left_dim;
775 }
776
777 let mut target_offset = 0;
779 let mut row_count = {
780 let m = &self.source.modules[0usize];
781 m.left.dimension(left_t) * m.right.dimension(right_t)
782 };
783 for s in 1..=self.source_s {
784 let source_module = &self.source.modules[s as usize]; let target_module = &self.target.modules[s as usize - 1]; let source_right_dim = source_module.right.dimension(right_t);
788 let source_left_dim = source_module.left.dimension(left_t);
789 let target_right_dim = target_module.right.dimension(right_t);
790 let target_left_dim = target_module.left.dimension(left_t);
791 assert_eq!(target_right_dim, source_right_dim);
792
793 let mut result = FpVector::new(p, target_left_dim);
794 for li in 0..source_left_dim {
795 self.left_cc.differential(s).apply_to_basis_element(
796 result.as_slice_mut(),
797 1,
798 left_t,
799 li,
800 );
801 for ri in 0..source_right_dim {
802 let mut row = matrix.row_mut(row_count);
803 for (i, x) in result.iter_nonzero() {
804 row.add_basis_element(target_offset + i * target_right_dim + ri, x);
805 }
806 row_count += 1;
807 }
808 result.set_to_zero();
809 }
810 target_offset += target_right_dim * target_left_dim;
811 }
812
813 matrix.segment(1, 1).add_identity();
814 matrix.row_reduce();
815
816 let mut index = 0;
817 let mut row = 0;
818 for s in 0..self.source_s as usize {
819 let target_module = &self.target.modules[s]; let target_right_dim = target_module.right.dimension(right_t);
822 let target_left_dim = target_module.left.dimension(left_t);
823
824 for li in 0..target_left_dim {
825 for ri in 0..target_right_dim {
826 if matrix.pivots()[index] >= 0 {
827 let true_index = self.target.offset(degree, s)
828 + self.target.modules[s].offset(degree, left_t)
829 + li * target_right_dim
830 + ri;
831 let mut entries = Vec::new();
832 let mut offset = 0;
833 for s_ in 0..=self.source_s as usize {
834 let dim = {
835 let m = &self.source.modules[s_];
836 m.left.dimension(left_t) * m.right.dimension(right_t)
837 };
838 if dim == 0 {
839 continue;
840 }
841
842 let mut entry = FpVector::new(p, dim);
843 entry.as_slice_mut().assign(
844 matrix
845 .row_segment(row, 1, 1)
846 .restrict(offset, offset + dim),
847 );
848
849 if !entry.is_zero() {
850 let true_slice_start = self.source.offset(degree, s_)
851 + self.source.modules[s_].offset(degree, left_t);
852 let true_slice_end = true_slice_start + dim;
853 entries.push((true_slice_start, true_slice_end, entry));
854 }
855
856 offset += dim;
857 }
858 assert!(quasi_inverse_list[true_index].is_none());
859 assert!(!entries.is_empty());
860 quasi_inverse_list[true_index] = Some(entries);
861 row += 1;
862 }
863 index += 1;
864 }
865 }
866 }
867 }
868 quasi_inverse_list
869 }
870 }
871
872 #[cfg(test)]
873 mod tests {
874 use ext::{
875 resolution_homomorphism::ResolutionHomomorphism, utils::construct,
876 yoneda::yoneda_representative_element,
877 };
878 use rstest::rstest;
879
880 use super::*;
881
882 #[rstest]
883 #[trace]
884 #[case(Bidegree::n_s(0, 1), &[1], &[1])]
885 #[case(Bidegree::n_s(0, 2), &[1], &[1])]
886 #[case(Bidegree::n_s(1, 1), &[1], &[1])]
887 #[case(Bidegree::n_s(3, 1), &[1], &[1])]
888 #[case(Bidegree::n_s(14, 4), &[1], &[1])]
889 fn test_square_cc(#[case] b: Bidegree, #[case] class: &[u32], #[case] output: &[u32]) {
890 let doubled_b = b + b;
891 let resolution = Arc::new(construct("S_2", None).unwrap());
892 let p = resolution.prime();
893 resolution.compute_through_bidegree(doubled_b);
894
895 let yoneda = Arc::new(yoneda_representative_element(
896 Arc::clone(&resolution),
897 b,
898 class,
899 ));
900
901 let square = Arc::new(TensorChainComplex::new(
902 Arc::clone(&yoneda),
903 Arc::clone(&yoneda),
904 ));
905 square.compute_through_bidegree(doubled_b);
906
907 let f = ResolutionHomomorphism::new(
908 "".to_string(),
909 Arc::clone(&resolution),
910 square,
911 Bidegree::zero(),
912 );
913 f.extend_step_raw(Bidegree::zero(), Some(vec![FpVector::from_slice(p, &[1])]));
914
915 f.extend(doubled_b);
916 let final_map = f.get_map(doubled_b.s());
917
918 for (i, &v) in output.iter().enumerate() {
919 assert_eq!(final_map.output(doubled_b.t(), i).len(), 1);
920 assert_eq!(final_map.output(doubled_b.t(), i).entry(0), v);
921 }
922 }
923 }
924}