support photon numbers over 2^31, new limit 2^63-1 after fixing fangq#41

src/mcx_core.cu

@@ -1528,9 +1528,10 @@ void mcx_run_simulation(Config *cfg,GPUInfo *gpu){
      float4 c0=cfg->srcdir;
      float3 maxidx=float3(cfg->dim.x,cfg->dim.y,cfg->dim.z);
      float *energy;
-     int timegate=0, totalgates, gpuid, gpuphoton=0,threadid=0;
+     int timegate=0, totalgates, gpuid, threadid=0;
+     size_t gpuphoton=0, photoncount=0;
 
-     unsigned int photoncount=0,printnum;
+     unsigned int printnum;
      unsigned int tic,tic0,tic1,toc=0,debuglen=MCX_DEBUG_REC_LEN;
      size_t fieldlen;
      uint3 cp0=cfg->crop0,cp1=cfg->crop1;
@@ -1858,7 +1859,7 @@ void mcx_run_simulation(Config *cfg,GPUInfo *gpu){
 }
 #pragma omp barrier
 
-     MCX_FPRINTF(cfg->flog,"\nGPU=%d (%s) threadph=%d extra=%d np=%d nthread=%d maxgate=%d repetition=%d\n",gpuid+1,gpu[gpuid].name,param.threadphoton,param.oddphotons,
+     MCX_FPRINTF(cfg->flog,"\nGPU=%d (%s) threadph=%d extra=%d np=%ld nthread=%d maxgate=%d repetition=%d\n",gpuid+1,gpu[gpuid].name,param.threadphoton,param.oddphotons,
            gpuphoton,gpu[gpuid].autothread,gpu[gpuid].maxgate,ABS(cfg->respin));
      MCX_FPRINTF(cfg->flog,"initializing streams ...\t");
      fflush(cfg->flog);
@@ -2115,7 +2116,7 @@ is more than what your have specified (%d), please use the -H option to specify
                    int detid;
 		   for(detid=1;detid<=(int)cfg->detnum;detid++){
 	               scale=0.f; // the cfg->normalizer and cfg.his.normalizer are inaccurate in this case, but this is ok
-		       for(i=0;i<cfg->nphoton;i++)
+		       for(size_t i=0;i<cfg->nphoton;i++)
 		           if(cfg->replay.detid[i]==detid)
 	                       scale+=cfg->replay.weight[i];
 	               if(scale>0.f)
@@ -2126,7 +2127,7 @@ is more than what your have specified (%d), please use the -H option to specify
 		   isnormalized=1;
 	       }else{
 	           scale=0.f;
-	           for(i=0;i<cfg->nphoton;i++)
+	           for(size_t i=0;i<cfg->nphoton;i++)
 	               scale+=cfg->replay.weight[i];
 	           if(scale>0.f)
                        scale=cfg->unitinmm/scale;

src/mcx_utils.c

@@ -697,16 +697,17 @@ void mcx_prepdomain(char *filename, Config *cfg){
 
 void mcx_loadconfig(FILE *in, Config *cfg){
      uint i,gates,itmp;
+     size_t count;
      float dtmp;
      char filename[MAX_PATH_LENGTH]={'\0'}, comment[MAX_PATH_LENGTH],strtypestr[MAX_SESSION_LENGTH]={'\0'},*comm;
 
      if(in==stdin)
      	fprintf(stdout,"Please specify the total number of photons: [1000000]\n\t");
-     MCX_ASSERT(fscanf(in,"%d", &(i) )==1);
-     if(cfg->nphoton==0) cfg->nphoton=i;
+     MCX_ASSERT(fscanf(in,"%ld", &(count) )==1);
+     if(cfg->nphoton==0) cfg->nphoton=count;
      comm=fgets(comment,MAX_PATH_LENGTH,in);
      if(in==stdin)
-     	fprintf(stdout,"%d\nPlease specify the random number generator seed: [1234567]\n\t",cfg->nphoton);
+     	fprintf(stdout,"%ld\nPlease specify the random number generator seed: [1234567]\n\t",cfg->nphoton);
      if(cfg->seed==0){
         MCX_ASSERT(fscanf(in,"%d", &(cfg->seed) )==1);
      }else{
@@ -1220,7 +1221,7 @@ int mcx_loadjson(cJSON *root, Config *cfg){
 void mcx_saveconfig(FILE *out, Config *cfg){
      uint i;
 
-     fprintf(out,"%d\n", (cfg->nphoton) );
+     fprintf(out,"%ld\n", (cfg->nphoton) );
      fprintf(out,"%d\n", (cfg->seed) );
      fprintf(out,"%f %f %f\n", (cfg->srcpos.x),(cfg->srcpos.y),(cfg->srcpos.z) );
      fprintf(out,"%f %f %f\n", (cfg->srcdir.x),(cfg->srcdir.y),(cfg->srcdir.z) );
@@ -1670,9 +1671,7 @@ void mcx_parsecmd(int argc, char* argv[], Config *cfg){
 				break;
 		     case 'n':
 		     	        i=mcx_readarg(argc,argv,i,&(np),"float");
-				if((uint)np>0x7FFFFFFF)
-				    mcx_error(-2,"the maximum photon number per session is 2^31 = 2.1e9, please use -r if more photons are needed.",__FILE__,__LINE__);
-				cfg->nphoton=(int)np;
+				cfg->nphoton=(size_t)np;
 		     	        break;
 		     case 't':
 		     	        i=mcx_readarg(argc,argv,i,&(cfg->nthread),"int");
@@ -2003,7 +2002,7 @@ where possible parameters include (the first value in [*|*] is the default)\n\
 == MC options ==\n\
 \n\
  -n [0|int]    (--photon)      total photon number (exponential form accepted)\n\
-                               max value:2.14e9, use -r to run larger workload\n\
+                               max accepted value:9.2234e+18 on 64bit systems\n\
  -r [1|+/-int] (--repeat)      if positive, repeat by r times,total= #photon*r\n\
                                if negative, divide #photon into r subsets\n\
  -b [1|0]      (--reflect)     1 to reflect photons at ext. boundary;0 to exit\n\

src/mcx_utils.h

@@ -123,7 +123,7 @@ typedef struct MCXGPUInfo {
  */
 
 typedef struct MCXConfig{
-	int nphoton;                  /**<total simulated photon number*/
+	size_t nphoton;               /**<total simulated photon number*/
         unsigned int nblocksize;      /**<thread block size*/
 	unsigned int nthread;         /**<num of total threads, multiple of 128*/
 	int seed;                     /**<random number generator seed*/

src/mcxlab.cpp

@@ -625,7 +625,7 @@ void mcx_set_field(const mxArray *root,const mxArray *item,int idx, Config *cfg)
             memcpy(cfg->replay.seed,mxGetData(item),arraydim[0]*arraydim[1]);
             cfg->seed=SEED_FROM_FILE;
             cfg->nphoton=arraydim[1];
-            printf("mcx.nphoton=%d;\n",cfg->nphoton);
+            printf("mcx.nphoton=%ld;\n",cfg->nphoton);
         }
     }else if(strcmp(name,"gpuid")==0){
         int len=mxGetNumberOfElements(item);