penna_dead.cpp

PEnna生物模型中的死亡和繁殖规律

#include<stdlib.h>
#include<time.h>
#include<stdio.h>
#define N0  4000      //不代表群体的初始值  群体的初始值需要用第一个函数统计
#define L 32          //基因串的长度
#define T 2           //坏基因的最大个数
#define R 8           //最小繁殖年龄
#define A  16         //病毒串的长度
#define K  60
#define pppp  0.999      //得传染病的概率
#define pp    0.001      //传染病的治愈率
#define ppp	  0.999      //传染病的死亡率
struct Animal 
{    
	    int s[L];   //个体基因串
	    int age;    //个体年龄
        int ill;    //个体的致病基因个数
		int  H;     // Hamming距离
		double p;   //得传染病的概率
	    int  z;               //表示个体的状态，值为1时表示健康，值为-1时表示传染病，值为0时表示死亡
};
struct Animal w[K][K];

int v[A];

int creat(int v[A],struct Animal w[K][K])  //个体初始化
{   
	int m1,m2,m5;
	//m3=100;
	m5=0;
    srand((unsigned)time(NULL));
	for(int i=0;i<A;i++)
	{
        v[i]=rand()%32;
	}
 
	for(i=0;i<N0;i++)
	{    
		 //if(m3==0)m3=100;
         m1=rand()%K;
		 m2=rand()%K;
		 w[m1][m2].age=rand()%L;
		 w[m1][m2].z=1;
		 //m3--;
		 //printf("w[%d][%d].age=%d\n",m1,m1,w[m1][m2].age);
			for(int n1=0;n1<L;n1++)
			{
	            w[m1][m2].s[n1]=rand()%2;
		//	printf("w[%d][%d].s[%d]=%d\n",m1,m2,n1,w[m1][m2].s[n1]);
			}
	}
	for(i=0;i<K;i++)
	{
        for(int j=0;j<K;j++)
		{
			if(w[i][j].z==1)m5++;
		}
	}
	//printf("m5=%d\n",m5);
	return(m5);
}

void Hamming(struct Animal w[K][K])
{   
    srand((unsigned)time(NULL));
	int n0,n1;
    for(int i=0;i<K;i++)
    {  
	   for(int j=0;j<K;j++)
	   {
        
         if(w[i][j].z==1)
		 {  
		    n0=rand()%L;
	        if(n0<=L-A)
			{
		          for(n1=n0;n1<n0+A;n1++)
				  {
			  
			             if(w[i][j].s[n1]==v[n1-n0])
						 {
				            w[i][j].H++;
						 }
						  
				  }
				 
		          w[i][j].p=w[i][j].H/(double)A;
                  //printf("w[%d][%d].p=%f\n",i,j,w[i][j].p);
		          if(w[i][j].p>=(1.0-pppp))
				  {
					  w[i][j].z=-1;
					 // N1++;
				  }
		          else w[i][j].z=1;
			}
	        else
			{
		            for( n1=L-A;n1<L;n1++)
					{
			              if(w[i][j].s[n1]==v[n1-n0])
						  {
				                w[i][j].H++;
						  }
					}
					//printf("w[%d][%d].H=%d\n",i,j,w[i][j].H);
		            w[i][j].p=w[i][j].H/(double)A;
				//	printf("w[%d][%d].p=%f\n",i,j,w[i][j].p);
		            if(w[i][j].p>=(1.0-pppp))
					{
					  w[i][j].z=-1;
					  //N1++;
					}
		            else w[i][j].z=1;
			}
		 }
	   }
	   
   }
    for(i=0;i<K;i++)
	{
	    for(int j=0;j<K;j++)
		{
		   if(w[i][j].z!=1&&w[i][j].z!=-1)w[i][j].z=0;
		}
	}
}


void computer(struct Animal w[K][K])
{
//用邻居规则对个体进行分析
	 int m4=creat(v,w);
     srand((unsigned)time(NULL));

	 for(int t=0;t<2;t++)     // 一个时刻到来，统计每个个体的致死基因个数
	 {     
		   printf("m4=%d\n",m4);
		   for(int i=0;i<K;i++)
		   {      
			      for(int j=0;j<K;j++)
				  {
                        
						if(w[i][j].z==1)
						{ 
							for(int n2=0;n2<w[i][j].age;n2++)if(w[i][j].s[n2]==1)w[i][j].ill++;
						//	printf("w[%d][%d].ill=%d\n",i,j,w[i][j].ill);
						    if(w[i][j].ill>=T||w[i][j].age>31)
							{
								w[i][j].z=0;
								m4--;
								//printf("m4=%d\n",m4);
							}
							else if((w[(i+1+K)%K][j].z==-1||w[(i-1+K)%K][j].z==-1||w[i][(j-1+K)%K].z==-1||w[i][(j+1+K)%K].z==-1)&&(rand()/32767.0>(1-pppp)))w[i][j].z=-1;
							else w[i][j].z=1;
						}
						else if(w[i][j].z==-1)
						{   
							for(int n2=0;n2<w[i][j].age;n2++)if(w[i][j].s[n2]==1)w[i][j].ill++;
						    if(w[i][j].ill>=T||w[i][j].age>31)
							{
								w[i][j].z=0;
								m4--;
								//printf("m4=%d\n",m4);
							}
							else if(rand()/32767.0>=(1.0-pp))w[i][j].z=1;
							else 
							{
                                  if(rand()/32767.0>=(1.0-ppp))
								  {
									  w[i][j].z=0;
									  m4--;
								  }
								  else w[i][j].z=-1;
							}
						}
						else 
						{
							if((w[(i+1+K)%K][j].age>=R&&w[(i+1+K)%K][j].z!=0)||(w[(i-1+K)%K][j].age>=R&&w[(i-1+K)%K][j].z!=0)||(w[i][(j-1+K)%K].age>=R&&w[i][(j-1+K)%K].z!=0)||(w[i][(j+1+K)%K].age>=R&&w[i][(j+1+K)%K].z!=0))
							{
								
							
								      w[i][j].age=0;
						              m4++;
									  do
									  {
                                          double p0=rand()/32767.0;
										  if(p0>0&&p0<=0.25)
											  {
											         if(w[(i+1+K)%K][j].z!=0){for(int n3=0;n3<31;n3++)w[i][j].s[n3]=w[(i+1+K)%K][j].s[n3];w[i][j].z=w[(i+1+K)%K][j].z;}
                                                     else continue;//w[i][j].z=w[(i+1+K)%K][j].z;
											  }
										  else if(p0>0.25&&p0<=0.5)
											  {
											          if(w[i][(j-1+K)%K].z!=0){for(int n3=0;n3<31;n3++)w[i][j].s[n3]=w[i][(j-1+K)%K].s[n3];w[i][j].z=w[i][(j-1+K)%K].z;}
											          else continue;//w[i][j].z=w[i][(j-1+K)%K].z;
											  }
										  else if(p0>0.5&&p0<=0.75)
											  {
											         if(w[(i-1+K)%K][j].z!=0){for(int n3=0;n3<31;n3++)w[i][j].s[n3]=w[(i-1+K)%K][j].s[n3];w[i][j].z=w[(i-1+K)%K][j].z;}
											         else continue;//w[i][j].z=w[(i-1+K)%K][j].z;
											  }
										  else
											  {
											        if(w[i][(j+1+K)%K].z!=0){for(int n3=0;n3<31;n3++)w[i][j].s[n3]=w[i][(j+1+K)%K].s[n3];w[i][j].z=w[i][(j+1+K)%K].z;}
											        else continue;//w[i][j].z=w[i][(j+1+K)%K].z;
											  }
									  }while(w[i][j].z!=0);

                                   w[i][j].z=1;
							}
						}
						
					 w[i][j].age++;	
				  }
				  
		   }
               
		  //printf("m4=%d\n",m4);
	 }
	 return;
}
void main()
{
	creat(v,w);
	//printf("m3=%d\n",m3);
	Hamming(w);
	computer(w);
}