@* Homology. These compute homology with $Z/pZ$ coefficients, for $p$ a prime. Since the coefficient is a field, this can just be done by computing the ranks of the boundary maps. @c #include "Morse.h" #include "globals.h" @ Compute the rank of a matrix with coefficients in $Z/2Z$ by Gaussian elimination (on the transpose). The |matrix| is $m\times n$ and the $i,j$ entry is at |matrix+m*j+i|. @c int gauss2(unsigned char *matrix,int m,int n) { unsigned char **perm; // column permutations int i,j,k,l; unsigned char *p,*q; perm = (unsigned char **)malloc(n*sizeof(unsigned char *)); for(i=0;i@; @@; @@; @i| to zero entries right of |i,j|@>@; } free(perm); return i; } @ @= for(;j= if(j==m) break; @ Note in fact |l| is ficticious and |j=l| already. @= if(k!=i) { p=perm[i]; perm[i]=perm[k]; perm[k]=p; } @ @i| to zero entries right of |i,j|@>= for(k=i+1;k@; @@; @@; @i| to zero entries right of |i,j|@>@; } free(perm); free(mult_tab); return i; } @ @d mod_sub(a,b) (((a)>=(b))?(a)-(b):(a)-(b)+prime) @i| to zero entries right of |i,j|@>= for(k=i+1;k@; j++; } j=0; list_read_init(crit[2]); while (!id_is_null(f=id_list_read(crit[2]))) // run through all critical faces { @@; @@; // mat[2] will be transposed j++; } @@; if(flag&4) @@; m = flag>>3; for(j=0;j<11;j++) { if((m>>j)&1) @@; } } @ @= { simplex_id v0, v1, w; v0 = FindGrad01(get_vertex(e,0),-100000); v1 = FindGrad01(get_vertex(e,1),-100000); p = mat[0]+n[0]*j; i=0; if(v0!=v1) { m=0; list_read_init(crit[0]); while(m<2) { w = id_list_read(crit[0]); if(id_is_null(w)) abort_message("missing critical vertex"); if(w==v0) {*p=1; m++;} else if(w==v1) {*p=-1; m++;} else *p=0; p++; i++; } } for(;i= { simplex_id te0,te1,te; int or0,or1; te0 = FindGrad23orientation(f,&or0); te1 = FindGrad23orientation(other_face_id(f),&or1); p = mat[2]+n[3]*j; i=0; m=0; k=2; if(id_is_null(te0)) k--; if(id_is_null(te1)) k--; list_read_init(crit[3]); while(m= { hlist *edges; SInt32 index; struct grad12_struct *q; simplex_id vlist[2]; hlist_initialize(&edges,sizeof(struct grad12_struct),113,2,0); find_all_grad12_paths(f,edges,8+4+2); p = mat[1]+n[1]*j; list_read_init(crit[1]); while(!id_is_null(e=id_list_read(crit[1]))) { get_edge_vertices(e,vlist); index = hlist_find(edges,vlist); if(index<0) *p=0; else { q = (struct grad12_struct *)hlist_entry(edges,index); *p = q->count; } p++; } hlist_abandon(&edges); } @ @= for(i=0;i= { for(i=0;i= { int prime = odd_primes[j]; for(i=0;i