Actual source code: embedgsvd.c
slepc-3.22.1 2024-10-28
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: static char help[] = "Computes a GSVD associated with a complex network analysis problem.\n"
12: "The command line options are:\n"
13: " -file <filename>, text file in TSV format containing undirected graph.\n"
14: " -beta <beta>, optional decay parameter.\n\n";
16: /*
17: Computes the Katz high-order proximity embedding of a graph via a partial GSVD.
19: [1] M. Ou et al. "Asymmetric transitivity preserving graph embedding" Proc. of
20: ACM SIGKDD 2016 - https://doi.org/10.1145/2939672.2939751
21: */
23: #include <slepcsvd.h>
24: #include "network.h"
26: PetscErrorCode SpectralRadius(Mat A,PetscReal *rho)
27: {
28: EPS eps;
29: PetscInt nconv;
30: PetscScalar kr,ki;
32: PetscFunctionBeginUser;
33: PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
34: PetscCall(EPSSetOperators(eps,A,NULL));
35: PetscCall(EPSSetProblemType(eps,EPS_NHEP));
36: PetscCall(EPSSetDimensions(eps,1,PETSC_DETERMINE,PETSC_DETERMINE));
37: PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE));
38: PetscCall(EPSSetFromOptions(eps));
39: PetscCall(EPSSolve(eps));
40: PetscCall(EPSGetConverged(eps,&nconv));
41: *rho = 0.0;
42: if (nconv) {
43: PetscCall(EPSGetEigenvalue(eps,0,&kr,&ki));
44: *rho = SlepcAbsEigenvalue(kr,ki);
45: }
46: PetscCall(EPSDestroy(&eps));
47: PetscFunctionReturn(PETSC_SUCCESS);
48: }
50: int main(int argc,char **argv)
51: {
52: Graph graph; /* directed graph */
53: Mat A,B; /* matrices */
54: SVD svd; /* singular value problem solver context */
55: PetscReal beta=0.0,rho=0.0;
56: char filename[PETSC_MAX_PATH_LEN];
57: PetscBool flg,terse;
59: PetscFunctionBeginUser;
60: PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
62: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
63: Load graph and build adjacency matrix
64: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
66: PetscCall(PetscOptionsGetString(NULL,NULL,"-graph",filename,sizeof(filename),&flg));
67: PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Must use -graph <name> to indicate a file containing the graph in TSV format");
69: PetscCall(GraphCreate(PETSC_COMM_WORLD,&graph));
70: PetscCall(GraphPreload(graph,filename));
71: PetscCheck(graph->type==GRAPH_DIRECTED,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"A directed graph is required");
72: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nGraph with %" PetscInt_FMT " vertices and %" PetscInt_FMT " edges.\n\n",graph->nvertices,graph->nedges));
73: if (graph->weight!=GRAPH_WEIGHT_UNWEIGHTED) PetscCall(PetscPrintf(PETSC_COMM_WORLD,"WARNING: ignoring weights in input graph\n"));
75: PetscCall(GraphPreallocate(graph,filename));
76: PetscCall(GraphLoadUnweighted(graph,filename));
78: PetscCall(MatViewFromOptions(graph->adjacency,NULL,"-adjacency"));
80: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
81: Build high-order proximity matrices (Katz): A=I-beta*G, B=beta*G
82: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
84: PetscCall(PetscOptionsGetReal(NULL,NULL,"-beta",&beta,NULL));
85: if (beta==0.0) { /* compute decay parameter beta */
86: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Computing spectral radius...\n"));
87: PetscCall(SpectralRadius(graph->adjacency,&rho));
88: if (rho) {
89: beta = 0.8*rho;
90: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Done, rho=%g, setting decay parameter beta=%g\n\n",(double)rho,(double)beta));
91: } else {
92: beta = 0.5;
93: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Failed, using default decay parameter beta=%g\n\n",(double)beta));
94: }
95: }
96: PetscCall(MatDuplicate(graph->adjacency,MAT_COPY_VALUES,&A));
97: PetscCall(MatScale(A,-beta));
98: PetscCall(MatShift(A,1.0));
99: PetscCall(MatDuplicate(graph->adjacency,MAT_COPY_VALUES,&B));
100: PetscCall(MatScale(B,beta));
102: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
103: Create the singular value solver and set various options
104: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
106: PetscCall(SVDCreate(PETSC_COMM_WORLD,&svd));
107: PetscCall(SVDSetOperators(svd,A,B));
108: PetscCall(SVDSetProblemType(svd,SVD_GENERALIZED));
109: PetscCall(SVDSetType(svd,SVDTRLANCZOS));
110: PetscCall(SVDSetFromOptions(svd));
112: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
113: Solve the problem and print solution
114: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
116: PetscCall(SVDSolve(svd));
118: /* show detailed info unless -terse option is given by user */
119: PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
120: if (terse) PetscCall(SVDErrorView(svd,SVD_ERROR_NORM,NULL));
121: else {
122: PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL));
123: PetscCall(SVDConvergedReasonView(svd,PETSC_VIEWER_STDOUT_WORLD));
124: PetscCall(SVDErrorView(svd,SVD_ERROR_NORM,PETSC_VIEWER_STDOUT_WORLD));
125: PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
126: }
127: PetscCall(SVDDestroy(&svd));
128: PetscCall(GraphDestroy(&graph));
129: PetscCall(MatDestroy(&A));
130: PetscCall(MatDestroy(&B));
131: PetscCall(SlepcFinalize());
132: return 0;
133: }
134: /*TEST
136: build:
137: depends: network.o
139: test:
140: args: -graph ${SLEPC_DIR}/share/slepc/datafiles/graphs/out.moreno_taro_taro -svd_nsv 4 -terse
141: filter: sed -e 's/4.38031/4.38032/' | sed -e 's/4.38033/4.38032/' | sed -e 's/3.75089/3.7509/' | sed -e 's/3.00071/3.00072/'
143: TEST*/