參考SIR模型原理:
①https://blog.csdn.net/Feng512275/article/details/82859526/
②https://zhuanlan.zhihu.com/p/104072104?app=zhihulite
代碼參考:
https://github.com/pholme/sir
參考章節:
一.互聯網絡模型
構造了兩種互連子網。一個是通過隨機或優先連接兩個相同的子網絡形成的,包括scale-free-scale網絡和e-mail-e-mail網絡。這種互聯網絡可以用來表示現實世界中連接不同社區網絡所形成的網絡。互連密度是用參數γ來測量的,定義為γ=L/N。L表示互連的個數,N表示一個子網的大小。我們構造的另一種互連網絡是將一個網絡隨機分成兩個大小相同的互連子網。這種互聯網絡的結構表明,包含不同類型節點的網絡被划分為不同的互聯子網。每個子網包含相同類型的節點。我們使用的單一網絡包括友誼網絡和自動系統網絡。
二.傳染病傳播模型
傳染病傳播模型本文采用的傳染病傳播模型是最基本、研究最充分的傳染病傳播模型。網絡的元素可以分為三個部分包括易感者(感染者)、感染者(感染者)和康復者(康復者)。在每個時間步,如果易受感染的節點直接連接到一個受感染的節點,則易受感染的節點被感染的概率為λ。參數λ稱為擴展率。同時,受感染的節點會出現被移除的節點,其概率δ稱為恢復率。但是對於兩個互連網絡(A-B網絡),我們需要分別指定這些過程,表示λ^A(λ^B)網絡中節點之間的擴展速率,λ^AB(λ^BA)網絡中節點A(B)到節點B(A)的擴展速率.δ^A(δ^B)表示網絡中節點之間的恢復速率。在我們的研究中,在沒有失去一般性的情況下,我們讓δ^A=δ^B=1,但是我們需要使用相對小的數值作為傳播率(λ^A,λ^B,λ^AB,λ^BA)和連接密度γ。因此感染率仍然很小,如果傳播可以達到人口的很大一部分,那么單個節點的作用就不再重要,傳播將覆蓋幾乎所有的網絡,而與網絡的起源無關
三.編碼思路:
總共有S(0)、 I(1)、 R(2)三類節點。初始化所有節點為I即nodeState為1,然后SIR函數中處理I(1)節點,分別為感染其他節點和恢復
int SIR(ALGraph* G,int a[],int arrayLength,double inf,double rec) //傳入的分別為網絡,感染節點數組,數組長度,感染率,恢復率 { int infTatal=arrayLength;//感染節點總數 int Inf[G->vexnum]; int newInf[G->vexnum]; int i=0; //給感染數組賦初值 for(i=0;i<infTatal;i++) { Inf[i]=a[i]; G->adjlist[i].nodeState=0;//傳入的數組為感染態 } for(i=0;i<infTatal;i++) { newInf[i]=0; } double infection=inf;//感染概率 int count=infTatal;//當前網絡中的感染個數 srand((unsigned)time(NULL)); //設置種子,用於隨機數產生 while(count>0)//還能繼續感染 { int newInfLength=0;//表示新感染節點的個數 for(i=0;i<count;i++) { int vert=Inf[i];//當前的感染點 EdgeNode* p; p=G->adjlist[vert].firstedge; //用當前節點去感染其他節點 while(p!=NULL) { double test=rand()/(double)RAND_MAX;//rand()產生隨機數為[1,32767],RAND_MAX設置為32767,那么test范圍[0.1] int nodej=p->adjvex;//記錄當前連接的節點 //如果隨機數比感染概率小(能感染),且節點狀態為易感染,就感染該節點 if(test<infection&&G->adjlist[nodej].nodeState==1) { newInf[newInfLength]=nodej; G->adjlist[nodej].nodeState=0;//被感染 newInfLength++; } p=p->next; } } //感染節點恢復(不包括上一步新感染的) for(i=0;i<count;i++) { double recovRate=rec; double test_1=rand()/(double)RAND_MAX; //恢復分兩種情況:1.能恢復,改變nodeState為2;2.不能恢復,放入新感染數組 if(test_1<recovRate) { G->adjlist[Inf[i]].nodeState=2; } else { newInf[newInfLength]=Inf[i]; newInfLength++; } } //newInf數組中元素兩個來源:1.易感染節點被感染;2.感染節點未恢復 //再把新感染的數組newInf交給Inf進行下一次循環 for(i=0;i<newInfLength;i++) { Inf[i]=newInf[i]; } count=newInfLength;//記錄當前新感染的個數,作為繼續循環的依據 } int recnum=0;//統計傳播結束后,處於恢復狀態的節點個數 for(i=0;i<G->vexnum;i++) { if(G->adjlist[i].nodeState==2) { recnum++; } } return recnum; }
加上之前的全部代碼:(這里使用是一個感染節點數組作為傳播源)
#include<stdio.h> #include<stdlib.h> #include<string.h> #include<malloc.h> #include<math.h> #include<time.h> #define MaxVertexNum 90000 #define RAND_MAX 0x7fff //邊表節點 typedef struct node { int adjvex; struct node *next; int visit; }EdgeNode; //頂點表節點 typedef struct vnode { int vertex; int KS;//k-core float ec;//特征向量中心性 int is_infected; int nodeState; int is_visit; int layer;//層 int degree;//該節點的度 EdgeNode* firstedge;//指向邊表的指針 }VertexNode; typedef VertexNode AdjList[MaxVertexNum]; //圖的結構 typedef struct { AdjList adjlist;//頂點表(類似於數組的) int vexnum;//頂點個數 int arcnum;//邊個數 }ALGraph; //返回文件行數(網絡邊數),有換行符"\n"就為一行 int lines(char* str) { int c; FILE* fp; int lines=0; fp=fopen(str,"r"); if(fp) { while((c=fgetc(fp))!=EOF) if(c=='\n') lines++; fclose(fp); } return lines; } //返回文件最大數(網絡節點數) int max(char* str) { FILE* fp; char* p; int line=lines(str); int i=0; int a=0; int b=0; fp=fopen(str,"r"); char buf[256]; if((fp=fopen(str,"r"))==NULL) { perror("fail to read"); exit(1); } //把文件的內容給buf while(fgets(buf,line,fp)!=NULL) { //p就沒有 p=buf; //這一步沒有成功賦值 sscanf(p,"%d %d",&a,&b);//輸入源為p //i始終為最大的 if(a>i) i=a; if(b>i) i=b; } return i; } //創建圖 void createAlgraph(ALGraph* G,char* str) { FILE* fp; int line=lines(str); int node=max(str);//其中最大數 G->vexnum=node+1;//因為是從0開始,所以+1,多了一個0 G->arcnum=line; fp=fopen(str,"r"); char buf[1024]; int len; int m; int n; EdgeNode* s; char* p; int a=0; int b=0; int i=0; char StrLine[1024]; //每個節點的頂點表(vn1(i),vn2(i)) int vn1[line];//這里本來要用vn[line],如果是vc++就不能通過編譯,有多少行就有多少(i,j) int vn2[line]; //頂點錄入 for(int j=0;j<G->vexnum;j++) { G->adjlist[j].vertex=j; G->adjlist[j].firstedge=NULL; G->adjlist[j].nodeState=1; } if((fp=fopen(str,"r"))==NULL) { perror("faile to read"); exit(1); } while(!feof(fp))//因為行數等於邊數,則讀取行數個就可以把其他的節點的連接讀完 { fgets(StrLine,1024,fp); sscanf(StrLine,"%d%d",&a,&b); vn1[i]=a; vn2[i]=b; i++; } //邊節點放入鏈表 //一行就是一個坐標,有多少行就有多少坐標 for(int k=0;k<line;k++)//有多少行就有多少節點,每個節點對應一個邊鏈 { m=vn1[k],n=vn2[k]; int ii=0; EdgeNode* p; p=G->adjlist[m].firstedge; while(p!=NULL) { if(p->adjvex==n) { ii=1; break; } p=p->next; } if(ii!=1) { s=(EdgeNode*)malloc(sizeof(EdgeNode)); s->adjvex=n;//相連接的頂點 s->next=NULL; s->next=G->adjlist[m].firstedge;//類似於自己寫的鏈表 G->adjlist[m].firstedge=s; //無向圖 有來回 s=(EdgeNode*)malloc(sizeof(EdgeNode)); s->adjvex=m; s->next=NULL; s->next=G->adjlist[n].firstedge; G->adjlist[n].firstedge=s; } } //深度為每個節點后面連接的鏈長度 EdgeNode* q; for( i=0;i<G->vexnum;i++) { int k=0; q=G->adjlist[i].firstedge; while(q!=NULL) { k++; q=q->next; } G->adjlist[i].degree=k; } } //打印鄰接表 void printGraph(ALGraph* G) { EdgeNode* s; for(int i=0;i<G->vexnum;i++) { s=G->adjlist[i].firstedge;//s為一個帶adjvex,next指針的邊表節點 while(s) { printf("(%d,%d)",G->adjlist[i].vertex,s->adjvex); s=s->next; } printf("\n"); } } //所屬層插入 void insertLayer(ALGraph* G,int layer) { for(int i=0;i<G->vexnum;i++) { G->adjlist[i].layer=layer; } } //打印度中心性 void printDegreeCentrality(ALGraph* G) { for(int i=0;i<G->vexnum;i++) { printf("node %d dgree centrality is:%d\n",i,G->adjlist[i].degree); } } //計算特征向量中心性 void eigenvector_centrality(ALGraph *G) { float e[G->vexnum];//記錄上一次的指標(最終的特征向量中心性指標 ,因為會把最終的計算賦值給e);下面都用指標代表特征向量指標 float e1[G->vexnum];//記錄這一次的指標 float max = 0;//這一次的最大指標 float max1 = 0;//記錄上一次最大指標 int flag=0;//當flag=1時,代表找到各個指標 for(int i=0; i<G->vexnum; i++) { e[i]=1;//將每個點初始化為1 e1[i]=0; } EdgeNode *p; p=(EdgeNode*)malloc(sizeof(EdgeNode)); //循環開始 while(flag==0) { max1=max;//max1為上一次的最大值 max=0; for (int i=0; i<G->vexnum; i++) { p=G->adjlist[i].firstedge; while(p!=NULL) { e1[i]+=e[p->adjvex];//第一次的計算結果為他們各自的度 p=p->next; } if(e1[i]>max) max=e1[i];//記錄本次的最大指標 } for(int i=0; i<G->vexnum; i++) { if(e[i]!=e1[i]) break; if(i==G->vexnum-1) flag=1;//兩次計算結果相同結束循環 } if((1.0/max1-1.0/max)<0.01&&(1.0/max1-1.0/max)>-0.01) flag=1;//當差值較小時也可結束循環 //保留這次的結果到e中,並且將ei重置為0,方便下次計算 for(int i=0; i<G->vexnum; i++) { e[i]=e1[i]; e1[i]=0; } } for(int i=0; i<G->vexnum; i++) { e[i]=e[i]/max; G->adjlist[i].ec=e[i]; } } /* 1.SIR傳播模型,改變node.state(0為感染,1為易感染,2為恢復。所有節點初始化為1) 2.SIR函數主要操作感染點,感染點會做兩件事:①感染易感節點;②感染節點恢復 3.傳播完成的標志為不能再感染新的節點,並返回處於恢復節點的個數 */ int SIR(ALGraph* G,int a[],int arrayLength,double inf,double rec) //傳入的分別為網絡,感染節點數組,數組長度,感染率,恢復率 { int infTatal=arrayLength;//感染節點總數 int Inf[G->vexnum]; int newInf[G->vexnum]; int i=0; //給感染數組賦初值 for(i=0;i<infTatal;i++) { Inf[i]=a[i]; G->adjlist[i].nodeState=0;//傳入的數組為感染態 } for(i=0;i<infTatal;i++) { newInf[i]=0; } double infection=inf;//感染概率 int count=infTatal;//當前網絡中的感染個數 srand((unsigned)time(NULL)); //設置種子,用於隨機數產生 while(count>0)//還能繼續感染 { int newInfLength=0;//表示新感染節點的個數 for(i=0;i<count;i++) { int vert=Inf[i];//當前的感染點 EdgeNode* p; p=G->adjlist[vert].firstedge; //用當前節點去感染其他節點 while(p!=NULL) { double test=rand()/(double)RAND_MAX;//rand()產生隨機數為[1,32767],RAND_MAX設置為32767,那么test范圍[0.1] int nodej=p->adjvex;//記錄當前連接的節點 //如果隨機數比感染概率小(能感染),且節點狀態為易感染,就感染該節點 if(test<infection&&G->adjlist[nodej].nodeState==1) { newInf[newInfLength]=nodej; G->adjlist[nodej].nodeState=0;//被感染 newInfLength++; } p=p->next; } } //感染節點恢復(不包括上一步新感染的) for(i=0;i<count;i++) { double recovRate=rec; double test_1=rand()/(double)RAND_MAX; //恢復分兩種情況:1.能恢復,改變nodeState為2;2.不能恢復,放入新感染數組 if(test_1<recovRate) { G->adjlist[Inf[i]].nodeState=2; } else { newInf[newInfLength]=Inf[i]; newInfLength++; } } //newInf數組中元素兩個來源:1.易感染節點被感染;2.感染節點未恢復 //再把新感染的數組newInf交給Inf進行下一次循環 for(i=0;i<newInfLength;i++) { Inf[i]=newInf[i]; } count=newInfLength;//記錄當前新感染的個數,作為繼續循環的依據 } int recnum=0;//統計傳播結束后,處於恢復狀態的節點個數 for(i=0;i<G->vexnum;i++) { if(G->adjlist[i].nodeState==2) { recnum++; } } return recnum; } int main() { char* str1="E:\\data_set\\netsci1.txt"; char* str2="E:\\data_set\\netsci2.txt"; char* str3="E:\\data_set\\netsci3.txt"; //G1、G2為兩個網絡,G3為他們的連接 ALGraph* G1; ALGraph* G2; ALGraph* G3; G1=(ALGraph*)malloc(sizeof(ALGraph)); G2=(ALGraph*)malloc(sizeof(ALGraph)); G3=(ALGraph*)malloc(sizeof(ALGraph)); //創建三個表的信息 createAlgraph(G1,str1);//分別插入圖的地址,連接頂點信息 createAlgraph(G2,str2); createAlgraph(G3,str3); //插入層數 insertLayer(G1,1); insertLayer(G2,2); //計算特征向量中心性 eigenvector_centrality(G1); int a[1]={3}; int arrayLenth=sizeof(a)/sizeof(int); //SIR printf("recnum=%d",SIR(G1,a,arrayLenth,0.6,0.2)); return 0; }