Actual source code: network.c

slepc-3.19.2 2023-09-05
Report Typos and Errors
  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: #include "network.h"

 13: const char *GraphTypes[] = {"sym","asym","bip",NULL};
 14: const char *GraphWeights[] = {"unweighted","positive","posweighted","signed","multisigned","weighted","multiweighted","dynamic","multiposweighted",NULL};

 16: PetscErrorCode GraphCreate(MPI_Comm comm,Graph *graph)
 17: {
 18:   PetscFunctionBeginUser;
 19:   PetscCall(PetscNew(graph));
 20:   (*graph)->comm = comm;
 21:   PetscFunctionReturn(PETSC_SUCCESS);
 22: }

 24: PetscErrorCode GraphDestroy(Graph *graph)
 25: {
 26:   PetscFunctionBeginUser;
 27:   if (!*graph) PetscFunctionReturn(PETSC_SUCCESS);
 28:   PetscCall(MatDestroy(&((*graph)->adjacency)));
 29:   PetscCall(PetscFree(*graph));
 30:   PetscFunctionReturn(PETSC_SUCCESS);
 31: }

 33: PetscErrorCode GraphPreload(Graph graph,char *filename)
 34: {
 35:   PetscInt    i,nval,src,dst;
 36:   PetscBool   flg;
 37:   PetscMPIInt rank;
 38:   FILE        *file;
 39:   char        gtype[64],gweight[64],line[PETSC_MAX_PATH_LEN];

 41:   PetscFunctionBeginUser;
 42:   PetscCallMPI(MPI_Comm_rank(graph->comm,&rank));
 43:   if (rank==0) {
 44:     PetscCall(PetscFOpen(PETSC_COMM_SELF,filename,"r",&file));
 45:     /* process first line of the file */
 46:     nval = fscanf(file,"%%%s%s\n",gtype,gweight);
 47:     PetscCheck(nval==2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file");
 48:     for (i=0;i<(int)sizeof(GraphTypes);i++) {
 49:       PetscCheck(GraphTypes[i],PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Unknown graph type %s",gtype);
 50:       PetscCall(PetscStrcmp(gtype,GraphTypes[i],&flg));
 51:       if (flg) { graph->type = (GraphType)i; break; }
 52:     }
 53:     for (i=0;i<(int)sizeof(GraphWeights);i++) {
 54:       PetscCheck(GraphWeights[i],PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Unknown graph weight %s",gweight);
 55:       PetscCall(PetscStrcmp(gweight,GraphWeights[i],&flg));
 56:       if (flg) { graph->weight = (GraphWeight)i; break; }
 57:     }
 58:     /* skip second line of the file if it is a comment */
 59:     if (!fgets(line,PETSC_MAX_PATH_LEN,file)) line[0] = 0;
 60:     if (line[0]=='%') {
 61:       if (!fgets(line,PETSC_MAX_PATH_LEN,file)) line[0] = 0;
 62:     }
 63:     graph->nedges = 1;
 64:     nval = sscanf(line,"%" PetscInt_FMT "%" PetscInt_FMT,&src,&dst);
 65:     PetscCheck(nval==2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file");
 66:     graph->nvertices = PetscMax(src,dst);
 67:     /* read rest of file to count lines */
 68:     while (fgets(line,PETSC_MAX_PATH_LEN,file)) {
 69:       nval = sscanf(line,"%" PetscInt_FMT "%" PetscInt_FMT,&src,&dst);
 70:       PetscCheck(nval==2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file");
 71:       graph->nedges++;
 72:       graph->nvertices = PetscMax(graph->nvertices,PetscMax(src,dst));
 73:     }
 74:     PetscCall(PetscFClose(PETSC_COMM_SELF,file));
 75:   }
 76:   PetscCallMPI(MPI_Bcast(&graph->type,1,MPIU_INT,0,PETSC_COMM_WORLD));
 77:   PetscCallMPI(MPI_Bcast(&graph->weight,1,MPIU_INT,0,PETSC_COMM_WORLD));
 78:   PetscCallMPI(MPI_Bcast(&graph->nvertices,1,MPIU_INT,0,PETSC_COMM_WORLD));
 79:   PetscCallMPI(MPI_Bcast(&graph->nedges,1,MPIU_INT,0,PETSC_COMM_WORLD));
 80:   PetscFunctionReturn(PETSC_SUCCESS);
 81: }

 83: PetscErrorCode GraphPreallocate(Graph graph,char *filename)
 84: {
 85:   PetscInt i,nval,src,dst,Istart,Iend,*d_nnz,*o_nnz;
 86:   FILE     *file;
 87:   char     line[PETSC_MAX_PATH_LEN];

 89:   PetscFunctionBeginUser;
 90:   PetscCheck(graph->nvertices,PETSC_COMM_SELF,PETSC_ERR_ORDER,"Must call GraphPreload() first");
 91:   PetscCall(MatDestroy(&(graph->adjacency)));
 92:   PetscCall(MatCreate(graph->comm,&graph->adjacency));
 93:   PetscCall(MatSetSizes(graph->adjacency,PETSC_DECIDE,PETSC_DECIDE,graph->nvertices,graph->nvertices));
 94:   PetscCall(MatSetType(graph->adjacency,MATAIJ));
 95:   PetscCall(MatSetUp(graph->adjacency));   /* otherwise MatGetOwnershipRange() cannot be called */
 96:   PetscCall(MatGetOwnershipRange(graph->adjacency,&Istart,&Iend));
 97:   PetscCall(PetscCalloc2(Iend-Istart,&d_nnz,Iend-Istart,&o_nnz));

 99:   /* all process read the file */
100:   PetscCall(PetscFOpen(PETSC_COMM_SELF,filename,"r",&file));
101:   nval = fscanf(file,"%%%*s%*s\n");   /* first line of the file */
102:   if (!fgets(line,PETSC_MAX_PATH_LEN,file)) line[0] = 0;
103:   if (line[0]=='%') { /* skip second line of the file if it is a comment */
104:     if (!fgets(line,PETSC_MAX_PATH_LEN,file)) line[0] = 0;
105:   }
106:   nval = sscanf(line,"%" PetscInt_FMT "%" PetscInt_FMT,&src,&dst);
107:   PetscCheck(nval==2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file");
108:   src--; dst--;     /* adjust from 1-based to 0-based */
109:   if (src>=Istart && src<Iend) {
110:     if (dst>=Istart && dst<Iend) d_nnz[src-Istart]++;
111:     else o_nnz[src-Istart]++;
112:   }
113:   /* read rest of file */
114:   while (fgets(line,PETSC_MAX_PATH_LEN,file)) {
115:     nval = sscanf(line,"%" PetscInt_FMT "%" PetscInt_FMT,&src,&dst);
116:     PetscCheck(nval==2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file");
117:     src--; dst--;     /* adjust from 1-based to 0-based */
118:     if (src>=Istart && src<Iend) {
119:       if (dst>=Istart && dst<Iend) d_nnz[src-Istart]++;
120:       else o_nnz[src-Istart]++;
121:     }
122:   }
123:   PetscCall(PetscFClose(PETSC_COMM_SELF,file));

125:   for (i=Istart;i<Iend;i++) d_nnz[i-Istart]++;  /* diagonal entries */
126:   PetscCall(MatSeqAIJSetPreallocation(graph->adjacency,0,d_nnz));
127:   PetscCall(MatMPIAIJSetPreallocation(graph->adjacency,0,d_nnz,0,o_nnz));
128:   PetscCall(PetscFree2(d_nnz,o_nnz));
129:   PetscFunctionReturn(PETSC_SUCCESS);
130: }

132: PetscErrorCode GraphLoadUnweighted(Graph graph,char *filename)
133: {
134:   PetscInt i,nval,src,dst,Istart,Iend;
135:   FILE     *file;
136:   char     line[PETSC_MAX_PATH_LEN];

138:   PetscFunctionBeginUser;
139:   PetscCheck(graph->adjacency,PETSC_COMM_SELF,PETSC_ERR_ORDER,"Must call GraphPreallocate() first");
140:   PetscCall(MatGetOwnershipRange(graph->adjacency,&Istart,&Iend));
141:   /* all process read the file */
142:   PetscCall(PetscFOpen(PETSC_COMM_SELF,filename,"r",&file));
143:   nval = fscanf(file,"%%%*s%*s\n");   /* first line of the file */
144:   if (!fgets(line,PETSC_MAX_PATH_LEN,file)) line[0] = 0;
145:   if (line[0]=='%') { /* skip second line of the file if it is a comment */
146:     if (!fgets(line,PETSC_MAX_PATH_LEN,file)) line[0] = 0;
147:   }
148:   nval = sscanf(line,"%" PetscInt_FMT "%" PetscInt_FMT,&src,&dst);
149:   PetscCheck(nval==2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file");
150:   src--; dst--;     /* adjust from 1-based to 0-based */
151:   if (src>=Istart && src<Iend) PetscCall(MatSetValue(graph->adjacency,src,dst,1.0,INSERT_VALUES));
152:   /* read rest of file */
153:   while (fgets(line,PETSC_MAX_PATH_LEN,file)) {
154:     nval = sscanf(line,"%" PetscInt_FMT "%" PetscInt_FMT,&src,&dst);
155:     PetscCheck(nval==2,PETSC_COMM_SELF,PETSC_ERR_FILE_UNEXPECTED,"Badly formatted input file");
156:     src--; dst--;     /* adjust from 1-based to 0-based */
157:     if (src>=Istart && src<Iend) PetscCall(MatSetValue(graph->adjacency,src,dst,1.0,INSERT_VALUES));
158:   }
159:   PetscCall(PetscFClose(PETSC_COMM_SELF,file));
160:   for (i=Istart;i<Iend;i++) PetscCall(MatSetValue(graph->adjacency,i,i,0.0,INSERT_VALUES));  /* diagonal entries */
161:   PetscCall(MatAssemblyBegin(graph->adjacency,MAT_FINAL_ASSEMBLY));
162:   PetscCall(MatAssemblyEnd(graph->adjacency,MAT_FINAL_ASSEMBLY));
163:   PetscFunctionReturn(PETSC_SUCCESS);
164: }