四参数随机生长法生成多孔介质代码

//四参数随机生长法，生成三维多孔介质//voidQuartetStructureGenerationSet(){ doubled1_6=0.02; doubled7_18=d1_6/2.; doubled19_26=d1_6/8.; doublen=0.8; doublecdd=0.001; doublenumtotal_need=(1-n)*NX*NY*NZ; intnumsoild=0; intTnumsoild; intsoild[M+1][3]; srand((unsigned)time(NULL)); //第一步，遍历所有网格，生成固相内核 for(i=0;i<=NX;i++) { for(j=0;j<=NY;j++) { for(k=0;k<=NZ;k++) { arrgrid[i][j][k]=0; soild[numsoild][0]=0; soild[numsoild][1]=0; soild[numsoild][2]=0; if(((rand()%1000)/1000.0)<cdd) { arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; numsoild=numsoild+1; } } } } Tnumsoild=numsoild; //第二步，从固相核向周围26个方向生长 while(Tnumsoild<numtotal_need) { for(intindex_soild=0;index_soild<Tnumsoild;index_soild++) { intindex_i=soild[index_soild][0]; intindex_j=soild[index_soild][1]; intindex_k=soild[index_soild][2]; //1向右方向生长 if(index_j<NY-1) { i=index_i; j=index_j+1; k=index_k; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d1_6) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //2向后方向生长 if(index_i>0) { i=index_i-1; j=index_j; k=index_k; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d1_6) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //3向左方向生长 if(index_j>0) { i=index_i; j=index_j-1; k=index_k; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d1_6) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //4向前方向生长 if(index_i<NX-1) { i=index_i+1; j=index_j; k=index_k; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d1_6) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //5向上方向生长 if(index_k<NZ-1) { i=index_i; j=index_j; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d1_6) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //6向下方向生长 if(index_k>0) { i=index_i; j=index_j; k=index_k-1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d1_6) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //7向水平方向右前生长 if(index_i<NX-1&&index_j<NY-1) { i=index_i+1; j=index_j+1; k=index_k; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //8向水平方向左前生长 if(index_i<NX-1&&index_j>0) { i=index_i+1; j=index_j-1; k=index_k; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //9向水平方向右后生长 if(index_i>0&&index_j<NY-1) { i=index_i-1; j=index_j+1; k=index_k; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //10向水平方向左后生长 if(index_i>0&&index_j>0) { i=index_i-1; j=index_j-1; k=index_k; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //11向右上方向生长 if(index_j<NY-1&&index_k<NZ-1) { i=index_i; j=index_j+1; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //12向右下方向生长 if(index_j<NY-1&&index_k>0) { i=index_i; j=index_j+1; k=index_k-1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //13向左上方向生长 if(index_j>0&&index_k<NZ-1) { i=index_i; j=index_j-1; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //14向左下方向生长 if(index_j>0&&index_k>0) { i=index_i; j=index_j-1; k=index_k-1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //15向前上方向生长 if(index_i<NX-1&&index_k<NZ-1) { i=index_i+1; j=index_j; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //16向前下方向生长 if(index_i<NX-1&&index_k>0) { i=index_i+1; j=index_j; k=index_k-1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //17向后上方向生长 if(index_i>0&&index_k<NZ-1) { i=index_i-1; j=index_j; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //18向后下方向生长 if(index_i>0&&index_k>0) { i=index_i-1; j=index_j; k=index_k-1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d7_18) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //19向右前上对角线方向生长 if(index_i<NX-1&&index_j<NY-1&&index_k<NZ-1) { i=index_i+1; j=index_j+1; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d19_26) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //20向右后上对角线方向生长 if(index_i>0&&index_j<NY-1&&index_k<NZ-1) { i=index_i-1; j=index_j+1; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d19_26) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //21向左后上对角线方向生长 if(index_i>0&&index_j>0&&index_k<NZ-1) { i=index_i-1; j=index_j-1; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d19_26) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //22向左前上对角线方向生长 if(index_i<NX-1&&index_j>0&&index_k<NZ-1) { i=index_i+1; j=index_j-1; k=index_k+1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d19_26) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //23向右前下对角线方向生长 if(index_i<NX-1&&index_j<NY-1&&index_k>0) { i=index_i+1; j=index_j+1; k=index_k-1; if(arrgrid[i][j][k]==0&&((rand()%1000)/1000.0)<d19_26) { numsoild=numsoild+1; arrgrid[i][j][k]=1; soild[numsoild][0]=i; soild[numsoild][1]=j; soild[numsoild][2]=k; } } //24向右后下对角线方向生长 if(index_i>0&&index_j<NY-1&&index_k>0) { i=index_i-1; j=index_j+1; k