Graphviz  2.41.20170921.2350
matrix_ops.c
Go to the documentation of this file.
1 /* $Id$ $Revision$ */
2 /* vim:set shiftwidth=4 ts=8: */
3 
4 /*************************************************************************
5  * Copyright (c) 2011 AT&T Intellectual Property
6  * All rights reserved. This program and the accompanying materials
7  * are made available under the terms of the Eclipse Public License v1.0
8  * which accompanies this distribution, and is available at
9  * http://www.eclipse.org/legal/epl-v10.html
10  *
11  * Contributors: See CVS logs. Details at http://www.graphviz.org/
12  *************************************************************************/
13 
14 
15 #include "matrix_ops.h"
16 #include "memory.h"
17 #include <stdlib.h>
18 #include <stdio.h>
19 #include <math.h>
20 
21 static double p_iteration_threshold = 1e-3;
22 
23 int
24 power_iteration(double **square_mat, int n, int neigs, double **eigs,
25  double *evals, int initialize)
26 {
27  /* compute the 'neigs' top eigenvectors of 'square_mat' using power iteration */
28 
29  int i, j;
30  double *tmp_vec = N_GNEW(n, double);
31  double *last_vec = N_GNEW(n, double);
32  double *curr_vector;
33  double len;
34  double angle;
35  double alpha;
36  int iteration = 0;
37  int largest_index;
38  double largest_eval;
39  int Max_iterations = 30 * n;
40 
41  double tol = 1 - p_iteration_threshold;
42 
43  if (neigs >= n) {
44  neigs = n;
45  }
46 
47  for (i = 0; i < neigs; i++) {
48  curr_vector = eigs[i];
49  /* guess the i-th eigen vector */
50  choose:
51  if (initialize)
52  for (j = 0; j < n; j++)
53  curr_vector[j] = rand() % 100;
54  /* orthogonalize against higher eigenvectors */
55  for (j = 0; j < i; j++) {
56  alpha = -dot(eigs[j], 0, n - 1, curr_vector);
57  scadd(curr_vector, 0, n - 1, alpha, eigs[j]);
58  }
59  len = norm(curr_vector, 0, n - 1);
60  if (len < 1e-10) {
61  /* We have chosen a vector colinear with prvious ones */
62  goto choose;
63  }
64  vecscale(curr_vector, 0, n - 1, 1.0 / len, curr_vector);
65  iteration = 0;
66  do {
67  iteration++;
68  cpvec(last_vec, 0, n - 1, curr_vector);
69 
70  right_mult_with_vector_d(square_mat, n, n, curr_vector,
71  tmp_vec);
72  cpvec(curr_vector, 0, n - 1, tmp_vec);
73 
74  /* orthogonalize against higher eigenvectors */
75  for (j = 0; j < i; j++) {
76  alpha = -dot(eigs[j], 0, n - 1, curr_vector);
77  scadd(curr_vector, 0, n - 1, alpha, eigs[j]);
78  }
79  len = norm(curr_vector, 0, n - 1);
80  if (len < 1e-10 || iteration > Max_iterations) {
81  /* We have reached the null space (e.vec. associated with e.val. 0) */
82  goto exit;
83  }
84 
85  vecscale(curr_vector, 0, n - 1, 1.0 / len, curr_vector);
86  angle = dot(curr_vector, 0, n - 1, last_vec);
87  } while (fabs(angle) < tol);
88  evals[i] = angle * len; /* this is the Rayleigh quotient (up to errors due to orthogonalization):
89  u*(A*u)/||A*u||)*||A*u||, where u=last_vec, and ||u||=1
90  */
91  }
92  exit:
93  for (; i < neigs; i++) {
94  /* compute the smallest eigenvector, which are */
95  /* probably associated with eigenvalue 0 and for */
96  /* which power-iteration is dangerous */
97  curr_vector = eigs[i];
98  /* guess the i-th eigen vector */
99  for (j = 0; j < n; j++)
100  curr_vector[j] = rand() % 100;
101  /* orthogonalize against higher eigenvectors */
102  for (j = 0; j < i; j++) {
103  alpha = -dot(eigs[j], 0, n - 1, curr_vector);
104  scadd(curr_vector, 0, n - 1, alpha, eigs[j]);
105  }
106  len = norm(curr_vector, 0, n - 1);
107  vecscale(curr_vector, 0, n - 1, 1.0 / len, curr_vector);
108  evals[i] = 0;
109 
110  }
111 
112 
113  /* sort vectors by their evals, for overcoming possible mis-convergence: */
114  for (i = 0; i < neigs - 1; i++) {
115  largest_index = i;
116  largest_eval = evals[largest_index];
117  for (j = i + 1; j < neigs; j++) {
118  if (largest_eval < evals[j]) {
119  largest_index = j;
120  largest_eval = evals[largest_index];
121  }
122  }
123  if (largest_index != i) { /* exchange eigenvectors: */
124  cpvec(tmp_vec, 0, n - 1, eigs[i]);
125  cpvec(eigs[i], 0, n - 1, eigs[largest_index]);
126  cpvec(eigs[largest_index], 0, n - 1, tmp_vec);
127 
128  evals[largest_index] = evals[i];
129  evals[i] = largest_eval;
130  }
131  }
132 
133  free(tmp_vec);
134  free(last_vec);
135 
136  return (iteration <= Max_iterations);
137 }
138 
139 
140 
141 void
142 mult_dense_mat(double **A, float **B, int dim1, int dim2, int dim3,
143  float ***CC)
144 {
145 /*
146  A is dim1 x dim2, B is dim2 x dim3, C = A x B
147 */
148 
149  double sum;
150  int i, j, k;
151  float *storage;
152  float **C = *CC;
153  if (C != NULL) {
154  storage = (float *) realloc(C[0], dim1 * dim3 * sizeof(A[0]));
155  *CC = C = (float **) realloc(C, dim1 * sizeof(A));
156  } else {
157  storage = (float *) malloc(dim1 * dim3 * sizeof(A[0]));
158  *CC = C = (float **) malloc(dim1 * sizeof(A));
159  }
160 
161  for (i = 0; i < dim1; i++) {
162  C[i] = storage;
163  storage += dim3;
164  }
165 
166  for (i = 0; i < dim1; i++) {
167  for (j = 0; j < dim3; j++) {
168  sum = 0;
169  for (k = 0; k < dim2; k++) {
170  sum += A[i][k] * B[k][j];
171  }
172  C[i][j] = (float) (sum);
173  }
174  }
175 }
176 
177 void
178 mult_dense_mat_d(double **A, float **B, int dim1, int dim2, int dim3,
179  double ***CC)
180 {
181 /*
182  A is dim1 x dim2, B is dim2 x dim3, C = A x B
183 */
184  double **C = *CC;
185  double *storage;
186  int i, j, k;
187  double sum;
188 
189  if (C != NULL) {
190  storage = (double *) realloc(C[0], dim1 * dim3 * sizeof(double));
191  *CC = C = (double **) realloc(C, dim1 * sizeof(double *));
192  } else {
193  storage = (double *) malloc(dim1 * dim3 * sizeof(double));
194  *CC = C = (double **) malloc(dim1 * sizeof(double *));
195  }
196 
197  for (i = 0; i < dim1; i++) {
198  C[i] = storage;
199  storage += dim3;
200  }
201 
202  for (i = 0; i < dim1; i++) {
203  for (j = 0; j < dim3; j++) {
204  sum = 0;
205  for (k = 0; k < dim2; k++) {
206  sum += A[i][k] * B[k][j];
207  }
208  C[i][j] = sum;
209  }
210  }
211 }
212 
213 void
214 mult_sparse_dense_mat_transpose(vtx_data * A, double **B, int dim1,
215  int dim2, float ***CC)
216 {
217 /*
218  A is dim1 x dim1 and sparse, B is dim2 x dim1, C = A x B
219 */
220 
221  float *storage;
222  int i, j, k;
223  double sum;
224  float *ewgts;
225  int *edges;
226  int nedges;
227  float **C = *CC;
228  if (C != NULL) {
229  storage = (float *) realloc(C[0], dim1 * dim2 * sizeof(A[0]));
230  *CC = C = (float **) realloc(C, dim1 * sizeof(A));
231  } else {
232  storage = (float *) malloc(dim1 * dim2 * sizeof(A[0]));
233  *CC = C = (float **) malloc(dim1 * sizeof(A));
234  }
235 
236  for (i = 0; i < dim1; i++) {
237  C[i] = storage;
238  storage += dim2;
239  }
240 
241  for (i = 0; i < dim1; i++) {
242  edges = A[i].edges;
243  ewgts = A[i].ewgts;
244  nedges = A[i].nedges;
245  for (j = 0; j < dim2; j++) {
246  sum = 0;
247  for (k = 0; k < nedges; k++) {
248  sum += ewgts[k] * B[j][edges[k]];
249  }
250  C[i][j] = (float) (sum);
251  }
252  }
253 }
254 
255 
256 
257 /* Copy a range of a double vector to a double vector */
258 void cpvec(double *copy, int beg, int end, double *vec)
259 {
260  int i;
261 
262  copy = copy + beg;
263  vec = vec + beg;
264  for (i = end - beg + 1; i; i--) {
265  *copy++ = *vec++;
266  }
267 }
268 
269 /* Returns scalar product of two double n-vectors. */
270 double dot(double *vec1, int beg, int end, double *vec2)
271 {
272  int i;
273  double sum;
274 
275  sum = 0.0;
276  vec1 = vec1 + beg;
277  vec2 = vec2 + beg;
278  for (i = end - beg + 1; i; i--) {
279  sum += (*vec1++) * (*vec2++);
280  }
281  return (sum);
282 }
283 
284 
285 /* Scaled add - fills double vec1 with vec1 + alpha*vec2 over range*/
286 void scadd(double *vec1, int beg, int end, double fac, double *vec2)
287 {
288  int i;
289 
290  vec1 = vec1 + beg;
291  vec2 = vec2 + beg;
292  for (i = end - beg + 1; i; i--) {
293  (*vec1++) += fac * (*vec2++);
294  }
295 }
296 
297 /* Scale - fills vec1 with alpha*vec2 over range, double version */
298 void vecscale(double *vec1, int beg, int end, double alpha, double *vec2)
299 {
300  int i;
301 
302  vec1 += beg;
303  vec2 += beg;
304  for (i = end - beg + 1; i; i--) {
305  (*vec1++) = alpha * (*vec2++);
306  }
307 }
308 
309 /* Returns 2-norm of a double n-vector over range. */
310 double norm(double *vec, int beg, int end)
311 {
312  return (sqrt(dot(vec, beg, end, vec)));
313 }
314 
315 
316 #ifndef __cplusplus
317 
318 /* inline */
319 void orthog1(int n, double *vec /* vector to be orthogonalized against 1 */
320  )
321 {
322  int i;
323  double *pntr;
324  double sum;
325 
326  sum = 0.0;
327  pntr = vec;
328  for (i = n; i; i--) {
329  sum += *pntr++;
330  }
331  sum /= n;
332  pntr = vec;
333  for (i = n; i; i--) {
334  *pntr++ -= sum;
335  }
336 }
337 
338 #define RANGE 500
339 
340 /* inline */
341 void init_vec_orth1(int n, double *vec)
342 {
343  /* randomly generate a vector orthogonal to 1 (i.e., with mean 0) */
344  int i;
345 
346  for (i = 0; i < n; i++)
347  vec[i] = rand() % RANGE;
348 
349  orthog1(n, vec);
350 }
351 
352 /* inline */
353 void
354 right_mult_with_vector(vtx_data * matrix, int n, double *vector,
355  double *result)
356 {
357  int i, j;
358 
359  double res;
360  for (i = 0; i < n; i++) {
361  res = 0;
362  for (j = 0; j < matrix[i].nedges; j++)
363  res += matrix[i].ewgts[j] * vector[matrix[i].edges[j]];
364  result[i] = res;
365  }
366  /* orthog1(n,vector); */
367 }
368 
369 /* inline */
370 void
371 right_mult_with_vector_f(float **matrix, int n, double *vector,
372  double *result)
373 {
374  int i, j;
375 
376  double res;
377  for (i = 0; i < n; i++) {
378  res = 0;
379  for (j = 0; j < n; j++)
380  res += matrix[i][j] * vector[j];
381  result[i] = res;
382  }
383  /* orthog1(n,vector); */
384 }
385 
386 /* inline */
387 void
388 vectors_subtraction(int n, double *vector1, double *vector2,
389  double *result)
390 {
391  int i;
392  for (i = 0; i < n; i++) {
393  result[i] = vector1[i] - vector2[i];
394  }
395 }
396 
397 /* inline */
398 void
399 vectors_addition(int n, double *vector1, double *vector2, double *result)
400 {
401  int i;
402  for (i = 0; i < n; i++) {
403  result[i] = vector1[i] + vector2[i];
404  }
405 }
406 
407 #ifdef UNUSED
408 /* inline */
409 void
410 vectors_mult_addition(int n, double *vector1, double alpha,
411  double *vector2)
412 {
413  int i;
414  for (i = 0; i < n; i++) {
415  vector1[i] = vector1[i] + alpha * vector2[i];
416  }
417 }
418 #endif
419 
420 /* inline */
421 void
422 vectors_scalar_mult(int n, double *vector, double alpha, double *result)
423 {
424  int i;
425  for (i = 0; i < n; i++) {
426  result[i] = vector[i] * alpha;
427  }
428 }
429 
430 /* inline */
431 void copy_vector(int n, double *source, double *dest)
432 {
433  int i;
434  for (i = 0; i < n; i++)
435  dest[i] = source[i];
436 }
437 
438 /* inline */
439 double vectors_inner_product(int n, double *vector1, double *vector2)
440 {
441  int i;
442  double result = 0;
443  for (i = 0; i < n; i++) {
444  result += vector1[i] * vector2[i];
445  }
446 
447  return result;
448 }
449 
450 /* inline */
451 double max_abs(int n, double *vector)
452 {
453  double max_val = -1e50;
454  int i;
455  for (i = 0; i < n; i++)
456  if (fabs(vector[i]) > max_val)
457  max_val = fabs(vector[i]);
458 
459  return max_val;
460 }
461 
462 #ifdef UNUSED
463 /* inline */
464 void orthogvec(int n, double *vec1, /* vector to be orthogonalized */
465  double *vec2 /* normalized vector to be orthogonalized against */
466  )
467 {
468  double alpha;
469  if (vec2 == NULL) {
470  return;
471  }
472 
473  alpha = -vectors_inner_product(n, vec1, vec2);
474 
475  vectors_mult_addition(n, vec1, alpha, vec2);
476 }
477 
478  /* sparse matrix extensions: */
479 
480 /* inline */
481 void mat_mult_vec(vtx_data * L, int n, double *vec, double *result)
482 {
483  /* compute result= -L*vec */
484  int i, j;
485  double sum;
486  int *edges;
487  float *ewgts;
488 
489  for (i = 0; i < n; i++) {
490  sum = 0;
491  edges = L[i].edges;
492  ewgts = L[i].ewgts;
493  for (j = 0; j < L[i].nedges; j++) {
494  sum -= ewgts[j] * vec[edges[j]];
495  }
496  result[i] = sum;
497  }
498 }
499 #endif
500 
501 /* inline */
502 void
504  int dim1, int dim2,
505  double *vector, double *result)
506 {
507  /* matrix is dim2 x dim1, vector has dim2 components, result=matrix^T x vector */
508  int i, j;
509 
510  double res;
511  for (i = 0; i < dim1; i++) {
512  res = 0;
513  for (j = 0; j < dim2; j++)
514  res += matrix[j][i] * vector[j];
515  result[i] = res;
516  }
517 }
518 
519 /* inline */
520 void
521 right_mult_with_vector_d(double **matrix,
522  int dim1, int dim2,
523  double *vector, double *result)
524 {
525  /* matrix is dim1 x dim2, vector has dim2 components, result=matrix x vector */
526  int i, j;
527 
528  double res;
529  for (i = 0; i < dim1; i++) {
530  res = 0;
531  for (j = 0; j < dim2; j++)
532  res += matrix[i][j] * vector[j];
533  result[i] = res;
534  }
535 }
536 
537 
538 /*****************************
539 ** Single precision (float) **
540 ** version **
541 *****************************/
542 
543 /* inline */
544 void orthog1f(int n, float *vec)
545 {
546  int i;
547  float *pntr;
548  float sum;
549 
550  sum = 0.0;
551  pntr = vec;
552  for (i = n; i; i--) {
553  sum += *pntr++;
554  }
555  sum /= n;
556  pntr = vec;
557  for (i = n; i; i--) {
558  *pntr++ -= sum;
559  }
560 }
561 
562 #ifdef UNUSED
563 /* inline */
564 void right_mult_with_vectorf
565  (vtx_data * matrix, int n, float *vector, float *result) {
566  int i, j;
567 
568  float res;
569  for (i = 0; i < n; i++) {
570  res = 0;
571  for (j = 0; j < matrix[i].nedges; j++)
572  res += matrix[i].ewgts[j] * vector[matrix[i].edges[j]];
573  result[i] = res;
574  }
575 }
576 
577 /* inline */
578 void right_mult_with_vector_fd
579  (float **matrix, int n, float *vector, double *result) {
580  int i, j;
581 
582  float res;
583  for (i = 0; i < n; i++) {
584  res = 0;
585  for (j = 0; j < n; j++)
586  res += matrix[i][j] * vector[j];
587  result[i] = res;
588  }
589 }
590 #endif
591 
592 /* inline */
594  (float *packed_matrix, int n, float *vector, float *result) {
595  /* packed matrix is the upper-triangular part of a symmetric matrix arranged in a vector row-wise */
596  int i, j, index;
597  float vector_i;
598 
599  float res;
600  for (i = 0; i < n; i++) {
601  result[i] = 0;
602  }
603  for (index = 0, i = 0; i < n; i++) {
604  res = 0;
605  vector_i = vector[i];
606  /* deal with main diag */
607  res += packed_matrix[index++] * vector_i;
608  /* deal with off diag */
609  for (j = i + 1; j < n; j++, index++) {
610  res += packed_matrix[index] * vector[j];
611  result[j] += packed_matrix[index] * vector_i;
612  }
613  result[i] += res;
614  }
615 }
616 
617 /* inline */
618 void
619 vectors_substractionf(int n, float *vector1, float *vector2, float *result)
620 {
621  int i;
622  for (i = 0; i < n; i++) {
623  result[i] = vector1[i] - vector2[i];
624  }
625 }
626 
627 /* inline */
628 void
629 vectors_additionf(int n, float *vector1, float *vector2, float *result)
630 {
631  int i;
632  for (i = 0; i < n; i++) {
633  result[i] = vector1[i] + vector2[i];
634  }
635 }
636 
637 /* inline */
638 void
639 vectors_mult_additionf(int n, float *vector1, float alpha, float *vector2)
640 {
641  int i;
642  for (i = 0; i < n; i++) {
643  vector1[i] = vector1[i] + alpha * vector2[i];
644  }
645 }
646 
647 /* inline */
648 void vectors_scalar_multf(int n, float *vector, float alpha, float *result)
649 {
650  int i;
651  for (i = 0; i < n; i++) {
652  result[i] = (float) vector[i] * alpha;
653  }
654 }
655 
656 /* inline */
657 void copy_vectorf(int n, float *source, float *dest)
658 {
659  int i;
660  for (i = 0; i < n; i++)
661  dest[i] = source[i];
662 }
663 
664 /* inline */
665 double vectors_inner_productf(int n, float *vector1, float *vector2)
666 {
667  int i;
668  double result = 0;
669  for (i = 0; i < n; i++) {
670  result += vector1[i] * vector2[i];
671  }
672 
673  return result;
674 }
675 
676 /* inline */
677 void set_vector_val(int n, double val, double *result)
678 {
679  int i;
680  for (i = 0; i < n; i++)
681  result[i] = val;
682 }
683 
684 /* inline */
685 void set_vector_valf(int n, float val, float* result)
686 {
687  int i;
688  for (i = 0; i < n; i++)
689  result[i] = val;
690 }
691 
692 /* inline */
693 double max_absf(int n, float *vector)
694 {
695  int i;
696  float max_val = -1e30f;
697  for (i = 0; i < n; i++)
698  if (fabs(vector[i]) > max_val)
699  max_val = (float) (fabs(vector[i]));
700 
701  return max_val;
702 }
703 
704 /* inline */
705 void square_vec(int n, float *vec)
706 {
707  int i;
708  for (i = 0; i < n; i++) {
709  vec[i] *= vec[i];
710  }
711 }
712 
713 /* inline */
714 void invert_vec(int n, float *vec)
715 {
716  int i;
717  float v;
718  for (i = 0; i < n; i++) {
719  if ((v = vec[i]) != 0.0)
720  vec[i] = 1.0f / v;
721  }
722 }
723 
724 /* inline */
725 void sqrt_vec(int n, float *vec)
726 {
727  int i;
728  double d;
729  for (i = 0; i < n; i++) {
730  /* do this in two steps to avoid a bug in gcc-4.00 on AIX */
731  d = sqrt(vec[i]);
732  vec[i] = (float) d;
733  }
734 }
735 
736 /* inline */
737 void sqrt_vecf(int n, float *source, float *target)
738 {
739  int i;
740  double d;
741  float v;
742  for (i = 0; i < n; i++) {
743  if ((v = source[i]) >= 0.0) {
744  /* do this in two steps to avoid a bug in gcc-4.00 on AIX */
745  d = sqrt(v);
746  target[i] = (float) d;
747  }
748  }
749 }
750 
751 /* inline */
752 void invert_sqrt_vec(int n, float *vec)
753 {
754  int i;
755  double d;
756  float v;
757  for (i = 0; i < n; i++) {
758  if ((v = vec[i]) > 0.0) {
759  /* do this in two steps to avoid a bug in gcc-4.00 on AIX */
760  d = 1. / sqrt(v);
761  vec[i] = (float) d;
762  }
763  }
764 }
765 
766 #ifdef UNUSED
767 /* inline */
768 void init_vec_orth1f(int n, float *vec)
769 {
770  /* randomly generate a vector orthogonal to 1 (i.e., with mean 0) */
771  int i;
772 
773  for (i = 0; i < n; i++)
774  vec[i] = (float) (rand() % RANGE);
775 
776  orthog1f(n, vec);
777 }
778 
779 
780  /* sparse matrix extensions: */
781 
782 /* inline */
783 void mat_mult_vecf(vtx_data * L, int n, float *vec, float *result)
784 {
785  /* compute result= -L*vec */
786  int i, j;
787  float sum;
788  int *edges;
789  float *ewgts;
790 
791  for (i = 0; i < n; i++) {
792  sum = 0;
793  edges = L[i].edges;
794  ewgts = L[i].ewgts;
795  for (j = 0; j < L[i].nedges; j++) {
796  sum -= ewgts[j] * vec[edges[j]];
797  }
798  result[i] = sum;
799  }
800 }
801 #endif
802 
803 #endif
void vectors_additionf(int n, float *vector1, float *vector2, float *result)
Definition: matrix_ops.c:629
void invert_sqrt_vec(int n, float *vec)
Definition: matrix_ops.c:752
void sqrt_vec(int n, float *vec)
Definition: matrix_ops.c:725
void orthog1(int n, double *vec)
Definition: matrix_ops.c:319
void vecscale(double *vec1, int beg, int end, double alpha, double *vec2)
Definition: matrix_ops.c:298
void invert_vec(int n, float *vec)
Definition: matrix_ops.c:714
void set_vector_val(int n, double val, double *result)
Definition: matrix_ops.c:677
int power_iteration(double **square_mat, int n, int neigs, double **eigs, double *evals, int initialize)
Definition: matrix_ops.c:24
void copy_vectorf(int n, float *source, float *dest)
Definition: matrix_ops.c:657
void cpvec(double *copy, int beg, int end, double *vec)
Definition: matrix_ops.c:258
void set_vector_valf(int n, float val, float *result)
Definition: matrix_ops.c:685
#define C
Definition: pack.c:29
void vectors_addition(int n, double *vector1, double *vector2, double *result)
Definition: matrix_ops.c:399
void copy_vector(int n, double *source, double *dest)
Definition: matrix_ops.c:431
void mult_dense_mat(double **A, float **B, int dim1, int dim2, int dim3, float ***CC)
Definition: matrix_ops.c:142
void right_mult_with_vector_transpose(double **matrix, int dim1, int dim2, double *vector, double *result)
Definition: matrix_ops.c:503
void square_vec(int n, float *vec)
Definition: matrix_ops.c:705
double max_absf(int n, float *vector)
Definition: matrix_ops.c:693
void vectors_mult_additionf(int n, float *vector1, float alpha, float *vector2)
Definition: matrix_ops.c:639
void right_mult_with_vector_f(float **matrix, int n, double *vector, double *result)
Definition: matrix_ops.c:371
void vectors_scalar_multf(int n, float *vector, float alpha, float *result)
Definition: matrix_ops.c:648
#define dot(v, w)
Definition: geom.c:400
double vectors_inner_productf(int n, float *vector1, float *vector2)
Definition: matrix_ops.c:665
int nedges
Definition: sparsegraph.h:80
void mult_sparse_dense_mat_transpose(vtx_data *A, double **B, int dim1, int dim2, float ***CC)
Definition: matrix_ops.c:214
void sqrt_vecf(int n, float *source, float *target)
Definition: matrix_ops.c:737
int * edges
Definition: sparsegraph.h:81
void mult_dense_mat_d(double **A, float **B, int dim1, int dim2, int dim3, double ***CC)
Definition: matrix_ops.c:178
double vectors_inner_product(int n, double *vector1, double *vector2)
Definition: matrix_ops.c:439
#define NULL
Definition: logic.h:39
void init_vec_orth1(int n, double *vec)
Definition: matrix_ops.c:341
double max_abs(int n, double *vector)
Definition: matrix_ops.c:451
float * ewgts
Definition: sparsegraph.h:82
#define alpha
Definition: shapes.c:3902
#define RANGE
Definition: matrix_ops.c:338
void orthog1f(int n, float *vec)
Definition: matrix_ops.c:544
void vectors_substractionf(int n, float *vector1, float *vector2, float *result)
Definition: matrix_ops.c:619
double norm(double *vec, int beg, int end)
Definition: matrix_ops.c:310
void vectors_subtraction(int n, double *vector1, double *vector2, double *result)
Definition: matrix_ops.c:388
void right_mult_with_vector_ff(float *packed_matrix, int n, float *vector, float *result)
Definition: matrix_ops.c:594
void scadd(double *vec1, int beg, int end, double fac, double *vec2)
Definition: matrix_ops.c:286
void vectors_scalar_mult(int n, double *vector, double alpha, double *result)
Definition: matrix_ops.c:422
#define N_GNEW(n, t)
Definition: agxbuf.c:20
void right_mult_with_vector_d(double **matrix, int dim1, int dim2, double *vector, double *result)
Definition: matrix_ops.c:521
void right_mult_with_vector(vtx_data *matrix, int n, double *vector, double *result)
Definition: matrix_ops.c:354