+
+
+/*
+ * Create the cross_link reference in order to have a byblock list.
+ */
+void create_cross_link(struct constraintmatrix *myconstraints, int k){
+
+ int i, j;
+ int blk;
+ struct sparseblock *p;
+ struct sparseblock *q;
+
+ struct sparseblock *prev;
+
+ /*
+ * Now, cross link.
+ */
+ prev=NULL;
+ for (i=1; i<=k; i++){
+ p=myconstraints[i].blocks;
+ while (p != NULL){
+ if (p->nextbyblock == NULL){
+ blk=p->blocknum;
+
+ /*
+ * link in the remaining blocks.
+ */
+ for (j=i+1; j<=k; j++){
+ q=myconstraints[j].blocks;
+
+ while (q != NULL){
+ if (q->blocknum == p->blocknum){
+ if (p->nextbyblock == NULL){
+ p->nextbyblock=q;
+ q->nextbyblock=NULL;
+ prev=q;
+ }
+ else{
+ prev->nextbyblock=q;
+ q->nextbyblock=NULL;
+ prev=q;
+ }
+ break;
+ }
+ q=q->next;
+ }
+ }
+ }
+ p=p->next;
+ }
+ }
+}
+
+
+
+
+void addentry(struct constraintmatrix *constraints,
+ struct blockmatrix *C,
+ int matno,
+ int blkno,
+ int indexi,
+ int indexj,
+ double ent,
+ int blocksize)
+{
+ struct sparseblock *p;
+ struct sparseblock *p_sav;
+
+ p=constraints[matno].blocks;
+
+ if (matno != 0.0) {
+ if (p == NULL){
+ /*
+ * We haven't yet allocated any blocks.
+ */
+ p=(struct sparseblock *)calloc(1, sizeof(struct sparseblock));
+
+ //two entries because this library ignores indices starting in zerox
+ p->constraintnum=matno;
+ p->blocknum=blkno;
+ p->numentries=1;
+ p->next=NULL;
+
+ p->entries=calloc(p->numentries+1, sizeof(double));
+ p->iindices=calloc(p->numentries+1, sizeof(int));
+ p->jindices=calloc(p->numentries+1, sizeof(int));
+
+ p->entries[p->numentries]=ent;
+ p->iindices[p->numentries]=indexi;
+ p->jindices[p->numentries]=indexj;
+
+ p->blocksize=blocksize;
+
+ constraints[matno].blocks=p;
+ } else {
+ /*
+ * We have some existing blocks. See whether this block is already
+ * in the chain.
+ */
+ while (p != NULL){
+ if (p->blocknum == blkno){
+ /*
+ * Found the right block.
+ */
+ p->constraintnum=matno;
+ p->blocknum=blkno;
+ p->numentries=p->numentries+1;
+
+ p->entries = realloc(p->entries, (p->numentries+1) * sizeof(double) );
+ p->iindices = realloc(p->iindices, (p->numentries+1) * sizeof(int) );
+ p->jindices = realloc(p->jindices, (p->numentries+1) * sizeof(int) );
+
+ p->entries[p->numentries]=ent;
+ p->iindices[p->numentries]=indexi;
+ p->jindices[p->numentries]=indexj;
+
+ return;
+ }
+ p_sav=p;
+ p=p->next;
+ }
+
+ /*
+ * If we get here, we have a non-empty structure but not the right block
+ * inside hence create a new structure.
+ */
+
+ p=(struct sparseblock *)calloc(1, sizeof(struct sparseblock));
+
+ //two entries because this library ignores indices starting in zerox
+ p->constraintnum=matno;
+ p->blocknum=blkno;
+ p->numentries=1;
+ p->next=NULL;
+
+ p->entries=calloc(p->numentries+1, sizeof(double));
+ p->iindices=calloc(p->numentries+1, sizeof(int));
+ p->jindices=calloc(p->numentries+1, sizeof(int));
+
+ p->entries[p->numentries]=ent;
+ p->iindices[p->numentries]=indexi;
+ p->jindices[p->numentries]=indexj;
+
+ p->blocksize=blocksize;
+
+ p_sav->next=p;
+ }
+ } else {
+ if (ent != 0.0){
+ int blksz=C->blocks[blkno].blocksize;
+ if (C->blocks[blkno].blockcategory == DIAG){
+ C->blocks[blkno].data.vec[indexi]=ent;
+ }else{
+ C->blocks[blkno].data.mat[ijtok(indexi,indexj,blksz)]=ent;
+ C->blocks[blkno].data.mat[ijtok(indexj,indexi,blksz)]=ent;
+ };
+ };
+
+ }
+}
+
