SNAP Library 3.0, User Reference  2016-07-20 17:56:49
SNAP, a general purpose, high performance system for analysis and manipulation of large networks
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
motifcluster.h File Reference
#include "Snap.h"

Go to the source code of this file.

Classes

class  TSweepCut
 
class  MotifCluster
 
class  ChibaNishizekiWeighter
 

Typedefs

typedef TVec< THash< TInt, TInt > > WeightVH
 

Enumerations

enum  MotifType {
  M1, M2, M3, M4,
  M5, M6, M7, M8,
  M9, M10, M11, M12,
  M13, bifan, semiclique, triangle,
  clique3, clique4, clique5, clique6,
  clique7, clique8, clique9, edge
}
 

Functions

void SymeigsSmallest (const TSparseColMatrix &A, int nev, TFltV &evals, TFullColMatrix &evecs, double tol=kDefaultTol, int maxiter=kMaxIter)
 

Variables

const double kDefaultTol = 1e-12
 
const int kMaxIter = 300
 

Typedef Documentation

typedef TVec< THash<TInt, TInt> > WeightVH

Definition at line 8 of file motifcluster.h.

Enumeration Type Documentation

enum MotifType
Enumerator
M1 
M2 
M3 
M4 
M5 
M6 
M7 
M8 
M9 
M10 
M11 
M12 
M13 
bifan 
semiclique 
triangle 
clique3 
clique4 
clique5 
clique6 
clique7 
clique8 
clique9 
edge 

Definition at line 10 of file motifcluster.h.

10  {
11  M1, // u --> v, v --> w, w --> u
12  M2, // u <--> v, v --> w, w --> u
13  M3, // u <--> v, v <--> w, w --> u
14  M4, // u <--> v, v <--> w, w <--> u
15  M5, // u --> v, v --> w, u --> w
16  M6, // u <--> v, w --> u, w --> v
17  M7, // u <--> v, u --> w, v --> w
18  M8, // u --> v, u --> w
19  M9, // u --> v, w --> u
20  M10, // v --> u, w --> u
21  M11, // u <--> v, u --> w
22  M12, // u <--> v, w --> u
23  M13, // u <--> v, u <--> w
24  bifan, // u --> w, u --> x, v --> w, v --> x
25  semiclique, // u -- v, u -- w, u -- x, v -- w, v -- x
26  triangle, // undirected cliques
27  clique3, // |
28  clique4, // |
29  clique5, // |
30  clique6, // |
31  clique7, // |
32  clique8, // |
33  clique9, // |
34  edge, // (undirected) edges
35 };

Function Documentation

void SymeigsSmallest ( const TSparseColMatrix A,
int  nev,
TFltV evals,
TFullColMatrix evecs,
double  tol = kDefaultTol,
int  maxiter = kMaxIter 
)

Definition at line 955 of file motifcluster.cpp.

956  {
957  // type of problem
958  int mode = 1;
959  // communication variable
960  int ido = 0;
961  // standard eigenvalue problem
962  char bmat = 'I';
963  // dimension of problem
964  int n = A.GetCols();
965  // type of eigenvector (smallest algebraic)
966  char which[] = "SA";
967  // Space for residual
968  double *resid = new double[n];
969  // Number of lanczos vectors
970  int ncv = MIN(MAX(2 * nev, 20), n - 1);
971  // Space for Lanczos basis vectors
972  double *V = new double[ncv * n];
973 
974  // dimension of basis vectors
975  int ldv = n;
976  // Various input / output data
977  int *iparam = new int[11];
978  // exact shifts
979  iparam[0] = 1;
980  // maximum number of iterations
981  iparam[2] = maxiter;
982  // mode
983  iparam[6] = mode;
984  // Output data
985  int *ipntr = new int[11];
986  // work array
987  double *workd = new double[3 * n];
988  // output workspace
989  int lworkl = ncv * (ncv + 8);
990  double *workl = new double[lworkl];
991  // Communication information
992  int info = 0;
993 
994  TFltV Ax(n);
995  // Communication loop. We keep applying A * x until ARPACK tells us to stop.
996  while (true) {
997  F77_NAME(dsaupd)(&ido, &bmat, &n, &which[0], &nev, &tol, &resid[0], &ncv,
998  &V[0], &ldv, &iparam[0], &ipntr[0], &workd[0], &workl[0],
999  &lworkl, &info);
1000  TFltV result(n);
1001  result.PutAll(0.0);
1002  double *load = &workd[ipntr[0] - 1];
1003  for (int i = 0; i < n; i++) {
1004  result[i] = load[i];
1005  }
1006 
1007  if (ido == 1) {
1008  // Need another matrix-vector product.
1009  A.Multiply(result, Ax);
1010  double *store = &workd[ipntr[1] - 1];
1011  for (int i = 0; i < n; i++) {
1012  store[i] = Ax[i];
1013  }
1014  } else if (ido == 99) {
1015  // We have finished.
1016  break;
1017  } else {
1018  // There was an error.
1019  TExcept::Throw("ARPACK error");
1020  }
1021  }
1022 
1023  // Extract the eigenvectors and eigenvalues
1024 
1025  // compute Ritz vectors
1026  int rvec(true);
1027  // Form of Ritz vectors
1028  char howmny = 'A';
1029  // select doesn't matter with howmny set to A
1030  int *select = new int[ncv];
1031  // Store Ritz values
1032  double *d = new double[nev];
1033  // Shift
1034  double sigmar = 0;
1035 
1036  F77_NAME(dseupd)(&rvec, &howmny, select, &d[0], &V[0], &ldv, &sigmar, &bmat,
1037  &n, &which[0], &nev, &tol, &resid[0], &ncv, &V[0], &ldv,
1038  &iparam[0], &ipntr[0], &workd[0], &workl[0], &lworkl, &info);
1039  if (info != 0) {
1040  TExcept::Throw("ARPACK error");
1041  }
1042 
1043  // Get eigenvalues and eigenvalues in sorted order
1044  TVec< TKeyDat<TFlt, TInt> > d_inds(nev);
1045  for (int i = 0; i < nev; i++) {
1046  d_inds[i] = TKeyDat<TFlt, TInt>(d[i], i);
1047  }
1048  d_inds.Sort();
1049 
1050  evals = TFltV(nev);
1051  evecs.ColV = TVec<TFltV>(nev);
1052  for (int i = 0; i < nev; i++) {
1053  double eval = d_inds[i].Key;
1054  int ind = d_inds[i].Dat;
1055  evals[ind] = eval;
1056  TFltV& col = evecs.ColV[ind];
1057  for (int j = ind * n; j < ind * n + n; j++) {
1058  col.Add(V[j]);
1059  }
1060  }
1061  evecs.RowN = n;
1062  evecs.ColN = nev;
1063 
1064  // Clean up
1065  delete[] resid;
1066  delete[] V;
1067  delete[] iparam;
1068  delete[] ipntr;
1069  delete[] workd;
1070  delete[] workl;
1071  delete[] select;
1072  delete[] d;
1073 }
Definition: ds.h:345
void F77_NAME() dsaupd(int *ido, char *bmat, int *n, char *which, int *nev, double *tol, double *resid, int *ncv, double *V, int *ldv, int *iparam, int *ipntr, double *workd, double *workl, int *lworkl, int *info)
static void Throw(const TStr &MsgStr)
Definition: ut.h:187
#define MIN(a, b)
Definition: bd.h:346
int GetCols() const
Definition: linalg.h:47
void Multiply(const TFltVV &B, int ColId, TFltV &Result) const
Definition: linalg.h:24
TVec< TFltV > ColV
Definition: linalg.h:131
TVec< TFlt > TFltV
Definition: ds.h:1531
void F77_NAME() dseupd(int *rvec, char *HowMny, int *select, double *d, double *Z, int *ldz, double *sigma, char *bmat, int *n, char *which, int *nev, double *tol, double *resid, int *ncv, double *V, int *ldv, int *iparam, int *ipntr, double *workd, double *workl, int *lworkl, int *info)
#define F77_NAME(name)
Definition: motifcluster.cpp:7
TSizeTy Add()
Adds a new element at the end of the vector, after its current last element.
Definition: ds.h:574
#define MAX(a, b)
Definition: bd.h:350

Variable Documentation

const double kDefaultTol = 1e-12

Definition at line 5 of file motifcluster.h.

const int kMaxIter = 300

Definition at line 6 of file motifcluster.h.