📄 replik.c
字号:
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <string.h>
/******* RANDOM NUMBER GENERATOR BY ZIFF **********/
#define A 471
#define B 1586
#define C 6988
#define D 9689
#define M 16383
#define RIMAX 2147483648.0 /* = 2^31 */
#define RandomInteger (++nd, ra[nd & M] = ra[(nd-A) & M] ^ ra[(nd-B) & M] ^ ra[(nd-C) & M] ^ ra[(nd-D) & M])
void seed(long seed);
static long ra[M+1], nd;
/**************************************************/
#define STEPNUMBER 600 /* number of Monte-Carlo steps */
#define SIDE 400
#define SIZE (SIDE*SIDE) /* lattice size */
#define a 0.02 /* addition mutation probability */#define b 0.02 /* deletion mutation probability */
#define fgvT_A 0.1 /* template function parameters */
#define fgvT_B 200.0001
#define fgvT_H 3.000001
#define fgvR_A 0.1 /* replicase function parameters */
#define fgvR_B 200.0001
#define fgvR_H 3.000001
#define fgvF_A 0.9 /* fidelity function parameters */#define fgvF_B 5.0001
#define fgvF_H 2.0001#define NUKLPAR 6.0000001 /* monomerdensity dependence parameter */#define DEATH 0.001 /* probability of decay */
#define DIFFRATE 0.0 /* diffusion rate */
#define INIC_N 5
#define CATCH 5
#define SAVEL 60000 /**************************/
#define RESULTL1 60000 /* intervals of data-saving */
#define RESULTL2 600000 /* */
#define RESULTL3 6000000 /**************************/
#define MAXL_N 70 /* in the program to give a maximum possible replicator length, */
/* before defining arrays */
typedef char matrix1[SIDE+1][SIDE+1]; /* array definitions */
typedef char matrix2[SIDE+1][SIDE+1][MAXL_N+1];
typedef long int vector1[SIZE+1];
typedef long int vector2[SIDE+1];
typedef long int vector3[9];
typedef long int vector4[5];
typedef double vector5[MAXL_N+1];typedef long int vector6[MAXL_N+1];
typedef double restomb1[MAXL_N+1]; typedef double restomb2[MAXL_N+1][9];
void seed(long seed); /* random number initialisation */
void trq_values(void); /* evaluates the template, replicase and fidelity values of replicators on each cell */
/* at initialisation */
void trq_update(void); /* overwrites template, replicase, fidelity values */
void inic(void); /* inoculation of the square lattice with short replicators */
void mc_stepi(void);
void xy_choose(void); /* picks a site at random */
void neighbourhood(void); /* counts the number and position of neighbouring empty and occupied sites among the*/
/* four nearest neighbours */
long randl(long); /* random number between 0 and parameter */
double randd(void); /* random double between 0 and 1 */
void toff_marg(void); /* Toffoli-Margolus algorithm */
void mc_step(void); /* performs one Monte-Carlo step, by picking each site of the grid in a random sequence,*/
/* and applying the update algorithm on them */
void save(void); /***************************/void result1(void); /* data recording procedures */void result2(void); /* */void result3(void); /***************************/matrix1 t_n, r_n, q_n, sq_length;matrix2 sq_seq;
vector5 fgv_template, fgv_replicase, fgv_fidelity;
double t, r, q; /* template, replicase activity, fidelity */
double n_denz; /* density of inbuilt (not free!) monomers */
long int length, length_o, length_n; /*variables, containing temporarily the length of replicators */
long int x, y, xs, ys; /* coordinates */
long int xs1, ys1, xs2, ys2; /* coordinates of replicase and empty site */
long int max_A, max_B, max_C, max_D; /* maximal contiguous sequnce of A,B,C,D monomers respectively */
long int i, j, k, l, m; /* cycle variables */
double rd;
long int rl;
double n_factor;
vector6 data;char prev_n, nukl;vector4 max_N, len_N;
vector1 availabl;
vector2 left, right, up, down;
long int available;
long int s;
int life, empty;
vector3 lifes, empties;
long int position;
int s1, s2;
FILE * tf;FILE * tf2;/***********************************************************
* Random generator initialization *
* by a simple Congruential generator *
***********************************************************/
void seed(long seed)
{
int i;
if(seed<0) { puts("SEED error."); exit(1); }
ra[0]= (long) fmod(16807.0*(double)seed, 2147483647.0);
for(i=1; i<=M; i++)
{
ra[i] = (long)fmod( 16807.0 * (double) ra[i-1], 2147483647.0);
}
}
long randl(long num) /* random integer number between 0 and num-1 */
{
return(RandomInteger % num);
}
double randd(void) /* random real number between 0 and 1 */
{
return((double) RandomInteger / RIMAX);
}
void trq_values(void) /* evaluates the template, replicase and fidelity values of replicators on each cell */ { /* at initialisation */ for(i=1; i<=SIDE; i++) { for(j=1; j<=SIDE; j++) { length=sq_length[i][j];
if(!length) { t_n[i][j]=0; r_n[i][j]=0; q_n[i][j]=0; continue; } max_N[1]=0; max_N[2]=0; max_N[3]=0; max_N[4]=0;
len_N[0]=0; len_N[1]=0; len_N[2]=0; len_N[3]=0; len_N[4]=0;
prev_n=0;
for(k=1; k<=length; k++)
{
if(sq_seq[i][j][k]==prev_n)
{
len_N[prev_n]++;
}
else
{
if(len_N[prev_n]>max_N[prev_n])
{
max_N[prev_n]=len_N[prev_n];
}
len_N[prev_n]=0; prev_n=sq_seq[i][j][k];
len_N[prev_n]++;
}
}
t_n[i][j]=max_N[1];
r_n[i][j]=max_N[2];
q_n[i][j]=max_N[3];
} } } void trq_update(void) /* overwrites template, replicase, fidelity values */ { length=sq_length[xs2][ys2];
max_N[1]=0; max_N[2]=0; max_N[3]=0; max_N[4]=0;
len_N[0]=0; len_N[1]=0; len_N[2]=0; len_N[3]=0; len_N[4]=0;
prev_n=0;
for(k=1; k<=length; k++)
{
if(sq_seq[xs2][ys2][k]==prev_n)
{
len_N[prev_n]++;
}
else
{
if(len_N[prev_n]>max_N[prev_n])
{
max_N[prev_n]=len_N[prev_n];
}
len_N[prev_n]=0; prev_n=sq_seq[xs2][ys2][k];
len_N[prev_n]++;
}
}
t_n[xs2][ys2]=max_N[1];
r_n[xs2][ys2]=max_N[2];
q_n[xs2][ys2]=max_N[3];
} void inic(void)
{
long int i;
seed(555);
position=0; for(i=1; i<=SIDE; i++)
{
left[i]=i-1;
right[i]=i+1;
up[i]=i-1;
down[i]=i+1;
}
left[1]=SIDE;
up[1]=SIDE;
right[SIDE]=1;
down[SIDE]=1;
/* inoculation of the square lattice with short replicators */
n_denz=0; for(i=1; i<=SIDE; i++)
{
for(j=1; j<=SIDE; j++)
{
sq_length[i][j]=0;
for(k=1; k<=MAXL_N; k++)
{
sq_seq[i][j][k]=0;
}
if(randd()<0.5)
{
sq_length[i][j]=INIC_N;
length=sq_length[i][j];
n_denz+=length; for(k=1; k<=length; k++)
{
sq_seq[i][j][k]=randl(4)+1;
}
}
}
}
n_denz=n_denz/SIZE; trq_values(); for(i=0; i<=MAXL_N; i++) /* sets the values of template, replicase, fidelity functions */ { fgv_template[i]=fgvT_A+(1-fgvT_A)*pow(i,fgvT_H)/(fgvT_B+pow(i,fgvT_H));
fgv_replicase[i]=fgvR_A+(1-fgvR_A)*pow(i,fgvR_H)/(fgvR_B+pow(i,fgvR_H));
fgv_fidelity[i]=fgvF_A+(1-fgvF_A)*pow(i,fgvF_H)/(fgvF_B+pow(i,fgvF_H));
}
fgv_template[0]=0; fgv_replicase[0]=0; }
void mc_stepi(void)
{
for (j=1; j<=SIZE; j++)
{
availabl[j]=j;
}
available=SIZE;
}
void xy_choose(void)
{
rl=randl(available)+1;
s=availabl[rl];
x=fmod((s-1),SIDE)+1;
y=(s-1)/SIDE+1;
availabl[rl]=availabl[available];
available--;
}
void neighbourhood(void)
{
life=0;
empty=0;
for(k=1; k<=4; k++)
{
if (k==1) { xs=left[x]; ys=y; }
else if (k==2) { xs=right[x]; ys=y; }
else if (k==3) { xs=x; ys=up[y]; }
else { xs=x; ys=down[y]; }
length_n=sq_length[xs][ys];
if(length_n>=5)
{
life++;
lifes[life]=k;
}
else if(length_n==0) {
empty++;
empties[empty]=k;
}
}
}
void toff_marg(void)
{
long int xd1,xd2,xd3,xd4,xd5,xd6,xd7;
long int yd1,yd2,yd3,yd4,yd5,yd6,yd7;
long int temp_l;
double temp_t, temp_r, temp_q;
char temp_seq[MAXL_N+1];
rl=randl(8)+1;
switch(rl)
{
case 1: case 2:
xd1=randl(SIDE)+1; yd1=randl(SIDE)+1;
xd2=xd1; yd2=down[yd1];
xd3=left[xd1]; yd3=down[yd1];
xd4=left[xd1]; yd4=yd1;
xd5=xd1; yd5=up[yd1];
xd6=right[xd1]; yd6=up[yd1];
xd7=right[xd1]; yd7=yd1;
break;
case 3: case 4:
xd1=randl(SIDE)+1; yd1=randl(SIDE)+1;
xd2=left[xd1]; yd2=yd1;
xd3=left[xd1]; yd3=up[yd1];
xd4=xd1; yd4=up[yd1];
xd5=right[xd1]; yd5=yd1;
xd6=right[xd1]; yd6=down[yd1];
xd7=xd1; yd7=down[yd1];
break;
case 5: case 6:
xd1=randl(SIDE)+1; yd1=randl(SIDE)+1;
xd2=xd1; yd2=down[yd1];
xd3=left[xd1]; yd3=down[yd1];
xd4=left[xd1]; yd4=yd1;
xd5=xd1; yd5=up[yd1];
xd6=right[xd1]; yd6=up[yd1];
xd7=right[xd1]; yd7=yd1;
break;
case 7: case 8:
xd1=randl(SIDE)+1; yd1=randl(SIDE)+1;
xd2=left[xd1]; yd2=yd1;
xd3=left[xd1]; yd3=up[yd1];
xd4=xd1; yd4=up[yd1];
xd5=right[xd1]; yd5=yd1;
xd6=right[xd1]; yd6=down[yd1];
xd7=xd1; yd7=down[yd1];
break;
}
switch(rl)
{
case 2: case 4: case 6: case 8: /* rotate left */
length=sq_length[xd1][yd1];
temp_l=length;
temp_t=t_n[xd1][yd1]; temp_r=r_n[xd1][yd1]; temp_q=q_n[xd1][yd1];
for(k=1; k<=length; k++)
{
temp_seq[k]=sq_seq[xd1][yd1][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
temp_seq[k]=0;
}
length=sq_length[xd2][yd2]; /* 1. exchange */
sq_length[xd1][yd1]=sq_length[xd2][yd2];
t_n[xd1][yd1]=t_n[xd2][yd2]; r_n[xd1][yd1]=r_n[xd2][yd2]; q_n[xd1][yd1]=q_n[xd2][yd2];
for(k=1; k<=length; k++)
{
sq_seq[xd1][yd1][k]=sq_seq[xd2][yd2][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd1][yd1][k]=0;
}
length=sq_length[xd3][yd3]; /* 2. exchange */
sq_length[xd2][yd2]=sq_length[xd3][yd3];
t_n[xd2][yd2]=t_n[xd3][yd3]; r_n[xd2][yd2]=r_n[xd3][yd3]; q_n[xd2][yd2]=q_n[xd3][yd3];
for(k=1; k<=length; k++)
{
sq_seq[xd2][yd2][k]=sq_seq[xd3][yd3][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd2][yd2][k]=0;
}
length=sq_length[xd4][yd4]; /* 3. exchange */
sq_length[xd3][yd3]=sq_length[xd4][yd4];
t_n[xd3][yd3]=t_n[xd4][yd4]; r_n[xd3][yd3]=r_n[xd4][yd4]; q_n[xd3][yd3]=q_n[xd4][yd4];
for(k=1; k<=length; k++)
{
sq_seq[xd3][yd3][k]=sq_seq[xd4][yd4][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd3][yd3][k]=0;
}
length=temp_l; /* 4. exchange */
sq_length[xd4][yd4]=temp_l;
t_n[xd4][yd4]=temp_t; r_n[xd4][yd4]=temp_r; q_n[xd4][yd4]=temp_q;
for(k=1; k<=length; k++)
{
sq_seq[xd4][yd4][k]=temp_seq[k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd4][yd4][k]=0;
}
break;
case 1: case 3: case 5: case 7: /* rotate right */
length=sq_length[xd1][yd1];
temp_l=length;
temp_t=t_n[xd1][yd1]; temp_r=r_n[xd1][yd1]; temp_q=q_n[xd1][yd1];
for(k=1; k<=length; k++)
{
temp_seq[k]=sq_seq[xd1][yd1][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
temp_seq[k]=0;
}
length=sq_length[xd4][yd4]; /* 1. exchange */
sq_length[xd1][yd1]=sq_length[xd4][yd4];
t_n[xd1][yd1]=t_n[xd4][yd4]; r_n[xd1][yd1]=r_n[xd4][yd4]; q_n[xd1][yd1]=q_n[xd4][yd4];
for(k=1; k<=length; k++)
{
sq_seq[xd1][yd1][k]=sq_seq[xd4][yd4][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd1][yd1][k]=0;
}
length=sq_length[xd3][yd3]; /* 2. exchange */
sq_length[xd4][yd4]=sq_length[xd3][yd3];
t_n[xd4][yd4]=t_n[xd3][yd3]; r_n[xd4][yd4]=r_n[xd3][yd3]; q_n[xd4][yd4]=q_n[xd3][yd3];
for(k=1; k<=length; k++)
{
sq_seq[xd4][yd4][k]=sq_seq[xd3][yd3][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd4][yd4][k]=0;
}
length=sq_length[xd2][yd2]; /* 3. exchange */
sq_length[xd3][yd3]=sq_length[xd2][yd2];
t_n[xd3][yd3]=t_n[xd2][yd2]; r_n[xd3][yd3]=r_n[xd2][yd2]; q_n[xd3][yd3]=q_n[xd2][yd2];
for(k=1; k<=length; k++)
{
sq_seq[xd3][yd3][k]=sq_seq[xd2][yd2][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd3][yd3][k]=0;
}
length=temp_l; /* 4. exchange */
sq_length[xd2][yd2]=temp_l;
t_n[xd2][yd2]=temp_t; r_n[xd2][yd2]=temp_r; q_n[xd2][yd2]=temp_q;
for(k=1; k<=length; k++)
{
sq_seq[xd2][yd2][k]=temp_seq[k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd2][yd2][k]=0;
}
break;
}
length=sq_length[xd2][yd2];
temp_l=length;
temp_t=t_n[xd2][yd2]; temp_r=r_n[xd2][yd2]; temp_q=q_n[xd2][yd2];
for(k=1; k<=length; k++)
{
temp_seq[k]=sq_seq[xd2][yd2][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
temp_seq[k]=0;
}
length=sq_length[xd7][yd7];
sq_length[xd2][yd2]=sq_length[xd7][yd7];
t_n[xd2][yd2]=t_n[xd7][yd7]; r_n[xd2][yd2]=r_n[xd7][yd7]; q_n[xd2][yd2]=q_n[xd7][yd7];
for(k=1; k<=length; k++)
{
sq_seq[xd2][yd2][k]=sq_seq[xd7][yd7][k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd2][yd2][k]=0;
}
length=temp_l;
sq_length[xd7][yd7]=temp_l;
t_n[xd7][yd7]=temp_t; r_n[xd7][yd7]=temp_r; q_n[xd7][yd7]=temp_q;
for(k=1; k<=length; k++)
{
sq_seq[xd7][yd7][k]=temp_seq[k];
}
for(k=(length+1); k<=MAXL_N; k++)
{
sq_seq[xd7][yd7][k]=0;
}
length=sq_length[xd4][yd4];
temp_l=length;
temp_t=t_n[xd4][yd4]; temp_r=r_n[xd4][yd4]; temp_q=q_n[xd4][yd4];
for(k=1; k<=length; k++)
{
temp_seq[k]=sq_seq[xd4][yd4][k];
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -