1 #include "BSprivate.h"
2
3 /* a data structure used only in factoring */
4 /* nodes in a list of i-nodes that need to be received/sent during */
5 /* factorization */
6 typedef struct __list_type {
7 int length; /* length of the columns in the i-node */
8 int num_cols; /* # of columns in the i-node */
9 int cur_len; /* don't know anymore */
10 FLOAT *d_buffer; /* don't know anymore */
11 int *i_buffer; /* don't know anymore */
12 int buf_size; /* don't know anymore */
13 int cl_ind; /* index of the clique, if local */
14 int in_ind; /* index of the i-node if local */
15 } list_type;
16
17 #include "BStree.h"
18
19 extern void BSupdt_matrix(BStree, int,int ,int *, FLOAT *, int, int *,
20 FLOAT *, int ,int ,BScl_2_inode *,BSinode *, BSprocinfo *);
21
22 /*+ BSfactorn - Compute the incomplete factor of a matrix
23
24 Input Parameters:
25 . A - The sparse matrix to be factored
26 . comm - the communication structure of the factoring
27 . procinfo - the usual processor info
28
29 Output Parameters:
30 . A - The factored sparse matrix
31
32 Returns:
33 0 if successful, otherwise a negative number whose absolute
34 value is the row number of the color (less one) where the
35 failure occured.
36
37 +*/
BSfactorn(BSpar_mat * A,BScomm * comm,BSprocinfo * procinfo)38 int BSfactorn(BSpar_mat *A, BScomm *comm, BSprocinfo *procinfo)
39 {
40 BMphase *to_phase, *from_phase;
41 BMmsg *msg;
42 int i, j, k;
43 int cl_ind, in_ind;
44 int count, size, ind;
45 int length, num_cols, gnum, offset;
46 int found;
47 int *intptr, *intptr2, *i_buffer;
48 FLOAT *d_buffer;
49 BScl_2_inode *clique2inode;
50 BSnumbering *color2clique;
51 BSinode *inodes;
52 int *data_ptr, data_len, msg_len;
53 FLOAT *msg_buf, *matrix;
54 FLOAT *nzptr;
55 char UP = 'U';
56 char DI = 'N';
57 char TR = 'N';
58 char SIDE = 'R';
59 FLOAT one = 1.0;
60 int info, fact_error;
61 list_type *list_data_ptr;
62 BStree_ptr tree_node_ptr;
63 BStree recv_tree;
64 BStree inode_tree;
65 int dummy = 0;
66 int *g_row_num, *iperm;
67
68 BMinit_comp_msg(comm->from_msg,procinfo); CHKERRN(0);
69 BMinit_comp_msg(comm->to_msg,procinfo); CHKERRN(0);
70 color2clique = A->color2clique;
71 clique2inode = A->clique2inode;
72 inodes = A->inodes->list;
73 g_row_num = A->global_row_num->numbers;
74 iperm = A->inv_perm->perm;
75
76 MY_INIT_TREE(inode_tree,(sizeof(int)*2));
77 for (i=0;i<color2clique->length-1;i++) {
78 for (cl_ind=color2clique->numbers[i];
79 cl_ind<color2clique->numbers[i+1];cl_ind++) {
80 for (in_ind=clique2inode->inode_index[cl_ind];
81 in_ind<clique2inode->inode_index[cl_ind+1];in_ind++) {
82 MY_INSERT_TREE_NODE(inode_tree,(inodes[in_ind].gcol_num),
83 found,tree_node_ptr,dummy);
84 data_ptr = (int *) MY_GET_TREE_DATA(tree_node_ptr);
85 data_ptr[0] = in_ind;
86 data_ptr[1] = cl_ind;
87 }
88 }
89 }
90
91 fact_error = FALSE;
92 /* now do this phase by phase */
93 for (i=0;i<color2clique->length-1;i++) {
94 to_phase = BMget_phase(comm->to_msg,i); CHKERRN(0);
95 from_phase = BMget_phase(comm->from_msg,i); CHKERRN(0);
96
97 /* send my part of each column in this color */
98 /* use the "from" part here */
99 msg = NULL;
100 while ((msg = BMnext_msg(from_phase,msg)) != NULL) {
101 CHKERRN(0);
102 /* allocate space for the message */
103 data_ptr = BMget_user(msg,&data_len); CHKERRN(0);
104 in_ind = data_ptr[1];
105 size = inodes[in_ind].length*inodes[in_ind].num_cols*sizeof(FLOAT)
106 + inodes[in_ind].length*sizeof(int);
107 MY_MALLOCN(msg_buf,(FLOAT *),size,1);
108 BMset_msg_ptr(msg,msg_buf); CHKERRN(0);
109 BMset_msg_size(msg,size); CHKERRN(0);
110 count = 0;
111 for (j=0;j<inodes[in_ind].num_cols;j++) {
112 nzptr = &(inodes[in_ind].nz[j*inodes[in_ind].length]);
113 for (k=0;k<inodes[in_ind].length;k++) {
114 msg_buf[count] = nzptr[k];
115 count++;
116 }
117 }
118 intptr = (int *) &(msg_buf[count]);
119 intptr2 = (int *) inodes[in_ind].row_num;
120 for (j=0;j<inodes[in_ind].length;j++) {
121 intptr[j] = g_row_num[iperm[intptr2[j]]];
122 }
123 BMsend_block_msg(from_phase,msg,procinfo); CHKERRN(0);
124 MY_FREE(msg_buf);
125 BMset_msg_ptr(msg,NULL); CHKERRN(0);
126 }
127 CHKERRN(0);
128 BMsend_block_msg(from_phase,NULL,procinfo); CHKERRN(0);
129
130
131 /* factor the cliques in this color */
132 for (cl_ind=color2clique->numbers[i];
133 cl_ind<color2clique->numbers[i+1];cl_ind++) {
134 if (procinfo->my_id == clique2inode->proc[cl_ind]) {
135 /* factor this clique */
136 /* stored in the upper diagonal */
137 size = clique2inode->d_mats[cl_ind].size;
138 matrix = clique2inode->d_mats[cl_ind].matrix;
139 DPOTRF(&UP,&size,matrix,&size,&info);
140 MLOG_flop((size*size*size)/6);
141 if (info != 0) fact_error = TRUE;
142 /* now invert the clique */
143 DTRTRI(&UP,&DI,&size,matrix,&size,&info);
144 MLOG_flop((size*size*size)/6);
145 if (info != 0) fact_error = TRUE;
146 }
147 }
148
149 /* figure out what I am going to receive and make room for it */
150 /* use the "to" part here */
151 MY_INIT_TREE(recv_tree,sizeof(list_type));
152 msg = NULL;
153 while ((msg = BMnext_msg(to_phase,msg)) != NULL) {
154 CHKERRN(0);
155 /* search list for this number */
156 data_ptr = BMget_user(msg,&data_len); CHKERRN(0);
157 gnum = data_ptr[0];
158 length = data_ptr[1];
159 num_cols = data_ptr[2];
160 MY_INSERT_TREE_NODE(recv_tree,gnum,found,tree_node_ptr,dummy);
161 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
162 if (found) {
163 list_data_ptr->length += length;
164 } else {
165 list_data_ptr->length = length;
166 list_data_ptr->num_cols = num_cols;
167 list_data_ptr->cur_len = 0;
168 list_data_ptr->d_buffer = NULL;
169 list_data_ptr->i_buffer = NULL;
170 list_data_ptr->cl_ind = -1;
171 list_data_ptr->in_ind = -1;
172 }
173 }
174 CHKERRN(0);
175 /* add any of my parts to the list, if necessary */
176 for (cl_ind=color2clique->numbers[i];
177 cl_ind<color2clique->numbers[i+1];cl_ind++) {
178 /* only work on columns that I own */
179 if (clique2inode->proc[cl_ind] != procinfo->my_id) continue;
180 for (in_ind=clique2inode->inode_index[cl_ind];
181 in_ind<clique2inode->inode_index[cl_ind+1];in_ind++) {
182 MY_INSERT_TREE_NODE(recv_tree,(inodes[in_ind].gcol_num),
183 found,tree_node_ptr,dummy);
184 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
185 if (found) {
186 list_data_ptr->cl_ind = cl_ind;
187 list_data_ptr->in_ind = in_ind;
188 list_data_ptr->length += inodes[in_ind].length;
189 } else {
190 list_data_ptr->length = inodes[in_ind].length;
191 list_data_ptr->num_cols = inodes[in_ind].num_cols;
192 list_data_ptr->cur_len = 0;
193 list_data_ptr->d_buffer = NULL;
194 list_data_ptr->i_buffer = NULL;
195 list_data_ptr->cl_ind = cl_ind;
196 list_data_ptr->in_ind = in_ind;
197 }
198 }
199 }
200 /* make sure that no one has been received out of order */
201 if (procinfo->error_check) {
202 MY_FIRST_IN_TREE(recv_tree,tree_node_ptr);
203 while (! IS_TREE_PTR_NULL(tree_node_ptr)) {
204 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
205 MY_NEXT_IN_TREE(tree_node_ptr);
206 if (list_data_ptr->cl_ind == -1) {
207 MY_SETERRCN(FACTOR_ERROR,"Received column out of order\n");
208 }
209 }
210 }
211
212 /* allocate space for each whole column */
213 MY_FIRST_IN_TREE(recv_tree,tree_node_ptr);
214 while (! IS_TREE_PTR_NULL(tree_node_ptr)) {
215 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
216 MY_NEXT_IN_TREE(tree_node_ptr);
217 list_data_ptr->buf_size =
218 (list_data_ptr->num_cols*list_data_ptr->length*sizeof(FLOAT)) +
219 (list_data_ptr->length*sizeof(int));
220 MY_MALLOCN(list_data_ptr->d_buffer,(FLOAT *),
221 list_data_ptr->buf_size,2);
222 list_data_ptr->i_buffer = (int *)
223 &(list_data_ptr->d_buffer[list_data_ptr->num_cols*
224 list_data_ptr->length]);
225 }
226
227 /* put my part of the column into the collected column */
228 MY_FIRST_IN_TREE(recv_tree,tree_node_ptr);
229 while (! IS_TREE_PTR_NULL(tree_node_ptr)) {
230 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
231 MY_NEXT_IN_TREE(tree_node_ptr);
232 in_ind = list_data_ptr->in_ind;
233 length = inodes[in_ind].length;
234 num_cols = inodes[in_ind].num_cols;
235 nzptr = inodes[in_ind].nz;
236 d_buffer = list_data_ptr->d_buffer;
237 for (j=0;j<num_cols;j++) {
238 for (k=0;k<length;k++) {
239 d_buffer[k] = nzptr[(j*length)+k];
240 }
241 d_buffer += list_data_ptr->length;
242 }
243 i_buffer = list_data_ptr->i_buffer;
244 intptr = inodes[in_ind].row_num;
245 /* while copying, convert i_buffer from new local to new global */
246 for (j=0;j<length;j++) {
247 i_buffer[j] = g_row_num[iperm[intptr[j]]];
248 }
249 list_data_ptr->cur_len += length;
250 }
251
252 /* collect all the parts of "my" columns */
253 while ((msg = BMrecv_block_msg(to_phase,procinfo)) != NULL) {
254 CHKERRN(0);
255 BMcheck_on_async_block(from_phase); CHKERRN(0);
256 msg_buf = (FLOAT *) BMget_msg_ptr(msg); CHKERRN(0);
257 BMset_msg_ptr(msg,NULL); CHKERRN(0);
258 data_ptr = BMget_user(msg,&msg_len); CHKERRN(0);
259 gnum = data_ptr[0];
260 length = data_ptr[1];
261 num_cols = data_ptr[2];
262 MY_SEARCH_TREE(recv_tree,gnum,found,tree_node_ptr);
263 if (found) {
264 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
265 d_buffer = &(list_data_ptr->d_buffer[list_data_ptr->cur_len]);
266 count = 0;
267 for (j=0;j<num_cols;j++) {
268 for (k=0;k<length;k++) {
269 d_buffer[k] = msg_buf[count];
270 count++;
271 }
272 d_buffer += list_data_ptr->length;
273 }
274 i_buffer = &(list_data_ptr->i_buffer[list_data_ptr->cur_len]);
275 intptr = (int *) &(msg_buf[num_cols*length]);
276 for (j=0;j<length;j++) {
277 i_buffer[j] = intptr[j];
278 }
279 list_data_ptr->cur_len += length;
280 } else {
281 MY_SETERRCN(FACTOR_ERROR,"Can't find matching column\n");
282 }
283 }
284 CHKERRN(0);
285 BMfree_block_msg(to_phase); CHKERRN(0);
286 BMfinish_on_async_block(from_phase); CHKERRN(0);
287
288 /* see if the factorizations went okay, if not, then return */
289 /* we have to wait until here to make sure all messages are received */
290 /* make anything that needs free'ing is free'ed */
291 /* everyone agrees on success/error */
292 GISUM(&fact_error,1,&info,procinfo->procset);
293 if (fact_error != 0) {
294 /* free up space and return with the negative phase number */
295 MY_FIRST_IN_TREE(recv_tree,tree_node_ptr);
296 while (! IS_TREE_PTR_NULL(tree_node_ptr)) {
297 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
298 MY_NEXT_IN_TREE(tree_node_ptr);
299 MY_FREE(list_data_ptr->d_buffer);
300 }
301 MY_FREE_TREE(recv_tree);
302 return(-(i+1));
303 }
304
305 /* factor these columns */
306 MY_FIRST_IN_TREE(recv_tree,tree_node_ptr);
307 while (! IS_TREE_PTR_NULL(tree_node_ptr)) {
308 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
309 MY_NEXT_IN_TREE(tree_node_ptr);
310 if (list_data_ptr->length <= 0) continue;
311 cl_ind = list_data_ptr->cl_ind;
312 in_ind = list_data_ptr->in_ind;
313 j = inodes[in_ind].num_cols;
314 offset = inodes[in_ind].gcol_num - clique2inode->g_offset[cl_ind];
315 size = clique2inode->d_mats[cl_ind].size;
316 matrix = clique2inode->d_mats[cl_ind].matrix;
317 matrix += ((offset*size)+offset);
318 DTRMM(&SIDE,&UP,&TR,&DI,&(list_data_ptr->length),&j,&one,
319 matrix,&size,list_data_ptr->d_buffer,&(list_data_ptr->length));
320 MLOG_flop(2*((j*(j+1))/2)*list_data_ptr->length*j);
321 }
322
323 /* send these columns out */
324 msg = NULL;
325 while ((msg = BMnext_msg(to_phase,msg)) != NULL) {
326 CHKERRN(0);
327 data_ptr = BMget_user(msg,&data_len); CHKERRN(0);
328 gnum = data_ptr[0];
329 MY_SEARCH_TREE(recv_tree,gnum,found,tree_node_ptr);
330 if (found) {
331 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
332 } else {
333 MY_SETERRCN(FACTOR_ERROR,"Can't find matching column\n");
334 }
335 BMset_msg_ptr(msg,list_data_ptr->d_buffer); CHKERRN(0);
336 BMset_msg_size(msg,list_data_ptr->buf_size); CHKERRN(0);
337 BMsend_block_msg(to_phase,msg,procinfo); CHKERRN(0);
338 BMset_msg_ptr(msg,NULL); CHKERRN(0);
339 }
340 CHKERRN(0);
341 BMsend_block_msg(to_phase,NULL,procinfo); CHKERRN(0);
342
343 /* retrieve my part of the column and put it in the data structure */
344 MY_FIRST_IN_TREE(recv_tree,tree_node_ptr);
345 while (! IS_TREE_PTR_NULL(tree_node_ptr)) {
346 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
347 MY_NEXT_IN_TREE(tree_node_ptr);
348 if (list_data_ptr->length <= 0) continue;
349 in_ind = list_data_ptr->in_ind;
350 length = inodes[in_ind].length;
351 num_cols = inodes[in_ind].num_cols;
352 nzptr = inodes[in_ind].nz;
353 d_buffer = list_data_ptr->d_buffer;
354 for (j=0;j<num_cols;j++) {
355 for (k=0;k<length;k++) {
356 nzptr[(j*length)+k] = d_buffer[k];
357 }
358 d_buffer += list_data_ptr->length;
359 }
360 }
361
362 /* now update the remaining matrix using the local columns */
363 MY_FIRST_IN_TREE(recv_tree,tree_node_ptr);
364 while (! IS_TREE_PTR_NULL(tree_node_ptr)) {
365 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
366 MY_NEXT_IN_TREE(tree_node_ptr);
367 if (list_data_ptr->length <= 0) continue;
368 in_ind = list_data_ptr->in_ind;
369 BSupdt_matrix(inode_tree,list_data_ptr->length,
370 list_data_ptr->num_cols,
371 list_data_ptr->i_buffer,list_data_ptr->d_buffer,
372 inodes[in_ind].length,
373 inodes[in_ind].row_num,inodes[in_ind].nz,
374 color2clique->numbers[i+1],
375 color2clique->numbers[color2clique->length-1],clique2inode,
376 inodes,procinfo); CHKERRN(0);
377 }
378
379 /* receive the incoming columns and update the remaining matrix */
380 /* strip out the "local" part and put it in the matrix */
381 while ((msg = BMrecv_block_msg(from_phase,procinfo)) != NULL) {
382 CHKERRN(0);
383 BMcheck_on_async_block(to_phase); CHKERRN(0);
384 msg_buf = (FLOAT *) BMget_msg_ptr(msg); CHKERRN(0);
385 BMset_msg_ptr(msg,NULL); CHKERRN(0);
386 data_ptr = BMget_user(msg,&msg_len); CHKERRN(0);
387 cl_ind = data_ptr[0];
388 in_ind = data_ptr[1];
389 /* let gnum be first global row location of my part */
390 gnum = g_row_num[iperm[inodes[in_ind].row_num[0]]];
391 num_cols = inodes[in_ind].num_cols;
392 length = inodes[in_ind].length;
393 msg_len = BMget_msg_size(msg); CHKERRN(0);
394 msg_len /= ((sizeof(FLOAT)*num_cols)+sizeof(int));
395
396 /* find my local part */
397 i_buffer = (int *) &(msg_buf[num_cols*msg_len]);
398 for (ind=0;ind<msg_len;ind++) {
399 if (gnum == i_buffer[ind]) break;
400 }
401 nzptr = inodes[in_ind].nz;
402 d_buffer = &(msg_buf[ind]);
403 for (j=0;j<num_cols;j++) {
404 for (k=0;k<length;k++) {
405 nzptr[(j*length)+k] = d_buffer[k];
406 }
407 d_buffer += msg_len;
408 }
409
410 /* update the remaining matrix */
411 BSupdt_matrix(inode_tree,msg_len,num_cols,i_buffer,msg_buf,length,
412 inodes[in_ind].row_num,nzptr,color2clique->numbers[i+1],
413 color2clique->numbers[color2clique->length-1],clique2inode,
414 inodes,procinfo); CHKERRN(0);
415 }
416 CHKERRN(0);
417 BMfree_block_msg(from_phase); CHKERRN(0);
418 BMfinish_on_async_block(to_phase); CHKERRN(0);
419
420 /* now free up the list */
421 MY_FIRST_IN_TREE(recv_tree,tree_node_ptr);
422 while (! IS_TREE_PTR_NULL(tree_node_ptr)) {
423 list_data_ptr = (list_type *) MY_GET_TREE_DATA(tree_node_ptr);
424 MY_NEXT_IN_TREE(tree_node_ptr);
425 MY_FREE(list_data_ptr->d_buffer);
426 }
427 MY_FREE_TREE(recv_tree);
428 }
429
430 MY_FREE_TREE(inode_tree);
431 BMfinish_comp_msg(comm->from_msg,procinfo); CHKERRN(0);
432 BMfinish_comp_msg(comm->to_msg,procinfo); CHKERRN(0);
433 return(0);
434 }
435
436 /*+ BSupdt_matrix - Update the matrix using a bunch of i-nodes
437
438 Input Parameters:
439 I don't really know what these things do now.
440
441 Output Parameters:
442 I don't really know what these things do now.
443
444 Returns:
445 void
446
447 Notes:
448 This code is very complicated and could be done much better.
449
450 +*/
BSupdt_matrix(BStree inode_tree,int len,int num_cols,int * updt_index,FLOAT * updt_inode,int len2,int * my_inode_ind,FLOAT * my_inode,int start,int end,BScl_2_inode * clique2inode,BSinode * inodes,BSprocinfo * procinfo)451 void BSupdt_matrix(BStree inode_tree,int len,int num_cols,int *updt_index,
452 FLOAT *updt_inode,int len2,int *my_inode_ind,FLOAT *my_inode,
453 int start,int end,BScl_2_inode *clique2inode,BSinode *inodes,
454 BSprocinfo *procinfo)
455 {
456 int i, j;
457 int updt_count, offset, toffset;
458 int cl_ind, in_ind;
459 BSpermutation *iperm;
460 FLOAT *work, *matrix;
461 int *intptr;
462 FLOAT *nzptr;
463 int tcount, count2, gnum, size;
464 int found;
465 BStree_ptr tree_node_ptr;
466 int *tree_data_ptr;
467
468 /* sort the updt_inode */
469 iperm = BSalloc_permutation(len); CHKERR(0);
470 for (i=0;i<len;i++) iperm->perm[i] = i;
471 BSheap_sort1(len,updt_index,iperm->perm); CHKERR(0);
472 MY_MALLOC(work,(FLOAT *),sizeof(FLOAT)*len,1);
473 for (i=0;i<num_cols;i++) {
474 BSiperm_dvec(&(updt_inode[i*len]),work,iperm); CHKERR(0);
475 for (j=0;j<len;j++) updt_inode[(i*len)+j] = work[j];
476 }
477 MY_FREE(work);
478 BSfree_permutation(iperm); CHKERR(0);
479
480 MY_MALLOC(work,(FLOAT *),sizeof(FLOAT)*num_cols,1);
481 updt_count = 0;
482 while (updt_count < len) {
483 MY_SEARCH_TREE(inode_tree,updt_index[updt_count],found,tree_node_ptr);
484 if (found) {
485 tree_data_ptr = (int *) MY_GET_TREE_DATA(tree_node_ptr);
486 in_ind = tree_data_ptr[0];
487 cl_ind = tree_data_ptr[1];
488 /********** update the inode */
489 tcount = 0;
490 count2 = 0;
491 intptr = inodes[in_ind].row_num;
492 while ((tcount<len2) && (count2<inodes[in_ind].length)) {
493 if (my_inode_ind[tcount] == intptr[count2]) {
494 for (j=0;j<num_cols;j++) {
495 work[j] = my_inode[(j*len2)+tcount];
496 }
497 nzptr = &(inodes[in_ind].nz[count2]);
498 for (i=0;i<inodes[in_ind].num_cols;i++) {
499 for (j=0;j<num_cols;j++) {
500 (*nzptr) -= work[j]*
501 updt_inode[(j*len)+updt_count+i];
502 }
503 nzptr += inodes[in_ind].length;
504 }
505 MLOG_flop(2*num_cols*inodes[in_ind].num_cols);
506 tcount++;
507 count2++;
508 } else if (my_inode_ind[tcount] < intptr[count2]) {
509 tcount++;
510 } else {
511 count2++;
512 }
513 }
514 /* update a protion of the clique if it is local */
515 if (clique2inode->proc[cl_ind] == procinfo->my_id) {
516 /***** update the clique (stored in UPPER triangle) */
517 gnum = clique2inode->g_offset[cl_ind];
518 offset = updt_index[updt_count]-gnum;
519 matrix = clique2inode->d_mats[cl_ind].matrix;
520 size = clique2inode->d_mats[cl_ind].size;
521 for (i=0;i<inodes[in_ind].num_cols;i++) {
522 tcount = updt_count + i;
523 for (j=0;j<num_cols;j++) {
524 work[j] = updt_inode[(j*len)+tcount];
525 }
526 while (tcount < len) {
527 toffset = updt_index[tcount]-gnum;
528 if (toffset < size) {
529 for (j=0;j<num_cols;j++) {
530 matrix[(toffset*size)+offset+i] -= work[j]*
531 updt_inode[(j*len)+tcount];
532 }
533 MLOG_flop(2*num_cols);
534 } else {
535 break;
536 }
537 tcount++;
538 }
539 }
540 }
541 /* end of clique update */
542 updt_count += inodes[in_ind].num_cols;
543 } else {
544 updt_count++;
545 }
546 }
547 MY_FREE(work);
548 }
549