Merge pull request #55 from fangq/trinity

Making trinity branch new master

commons/Makefile_common.mk

@@ -20,14 +20,18 @@ MMCSRC :=$(MMCDIR)/src
 
 CXX        := g++
 AR         := $(CC)
+CUDACC     :=nvcc
 BIN        := bin
 BUILT      := built
 BINDIR     := $(BIN)
 OBJDIR 	   := $(BUILT)
 CCFLAGS    += -c -Wall -g -DMCX_EMBED_CL -fno-strict-aliasing#-pedantic -std=c99 -mfpmath=sse -ffast-math -mtune=core2
-INCLUDEDIR := $(MMCDIR)/src
+INCLUDEDIR := $(MMCDIR)/src -I$(MMCDIR)/src/zmat/easylzma -I$(MMCDIR)/src/ubj
 AROUTPUT   += -o
-MAKE       := make
+MAKE       ?= make
+
+ZMATLIB    :=libzmat.a
+USERARFLAGS?=$(ZMATLIB) -lz
 
 LIBOPENCLDIR ?= /usr/local/cuda/lib64
 LIBOPENCL  ?=-lOpenCL
@@ -36,6 +40,8 @@ EXTRALIB   += -lm -lstdc++ -L$(LIBOPENCLDIR)
 OPENMP     := -fopenmp
 OPENMPLIB  := -fopenmp
 FASTMATH   := #-ffast-math
+CUCCOPT    +=-Xcompiler $(OPENMP) -use_fast_math -Xptxas -O3,-v
+CUDA_STATIC=--cudart static -Xcompiler "-static-libgcc -static-libstdc++"
 
 ECHO	   := echo
 MKDIR      := mkdir
@@ -83,8 +89,27 @@ else ifeq ($(findstring Darwin,$(PLATFORM)), Darwin)
     LIBOPENCL=-framework OpenCL
     LIBOPENCLDIR=/System/Library/Frameworks/OpenCL.framework/Versions/A
     OPENMPLIB=-static-libgcc /usr/local/lib/libgomp.a
+    CUDA_STATIC=--cudart static
+endif
+
+ifeq ($(BACKEND),ocelot)
+  LINKOPT=-L/usr/local/lib `OcelotConfig -l` -ltinfo
+  CUCCOPT=-D__STRICT_ANSI__ -g #--maxrregcount 32
+else ifeq ($(BACKEND),cudastatic)
+  ifeq ($(findstring Darwin,$(PLATFORM)), Darwin)
+      CUDART=-lcudadevrt -lcudart_static -ldl -static-libgcc -static-libstdc++
+  else
+      CUDART=-lcudadevrt -lcudart_static -ldl -lrt -static-libgcc -static-libstdc++
+  endif
+  LINKOPT=-L/usr/local/cuda/lib -lm $(CUDART)
+  CUCCOPT+=-g -lineinfo -Xcompiler -Wall#-arch compute_20 #--maxrregcount 32
+else
+  LINKOPT=-L/usr/local/cuda/lib -lm $(CUDART)
+  CUCCOPT+=-g -lineinfo -Xcompiler -Wall#-arch compute_20 #--maxrregcount 32
 endif
 
+CUGENCODE?=-arch=sm_30
+
 INCLUDEDIR+=$(INCLUDEDIRS)
 EXTRALIB+=$(LIBOPENCL)
 
@@ -154,7 +179,10 @@ mex mexomp:     ARFLAGS+=mmclab.cpp -I$(INCLUDEDIR)
 oct:            BINARY=mmc.mex
 oct octomp:     ARFLAGS+=--mex mmclab.cpp -I$(INCLUDEDIR)
 oct octomp:     AR=CC=$(CC) CXX=$(CXX) LFLAGS='$(LFLAGS) $(OPENMPLIB) $(LIBOPENCL) $(MEXLINKOPT)' CPPFLAGS='$(CCFLAGS) $(USERCCFLAGS) -std=c++11' $(USEROCTOPT) $(MKOCT)
-oct octomp:     USERARFLAGS=-o $(BINDIR)/mmc
+oct octomp:     USERARFLAGS+=-o $(BINDIR)/mmc
+
+debug:     sse
+debug:     CUCCOPT+=-DMCX_DEBUG
 
 TARGETSUFFIX:=$(suffix $(BINARY))
 
@@ -172,20 +200,28 @@ ifeq ($(TARGETSUFFIX),.a)
 	OPENMPLIB  :=
 endif
 
-all release sse ssemath prof omp mex oct mexomp octomp pnacl web: $(SUBDIRS) $(BINDIR)/$(BINARY)
+cuda: sse
+
+all release sse ssemath prof omp mex oct mexomp octomp pnacl web debug cuda: $(SUBDIRS) $(BINDIR)/$(BINARY)
 
 $(SUBDIRS):
 	$(MAKE) -C $@ --no-print-directory
 
 makedirs:
 	@if test ! -d $(OBJDIR); then $(MKDIR) $(OBJDIR); fi
+	@if test ! -d $(OBJDIR)/ubj; then $(MKDIR) $(OBJDIR)/ubj; fi
 	@if test ! -d $(BINDIR); then $(MKDIR) $(BINDIR); fi
 
 makedocdir:
 	@if test ! -d $(DOCDIR); then $(MKDIR) $(DOCDIR); fi
 
 .SUFFIXES : $(OBJSUFFIX) .cpp
 
+
+##  Compile .cu files ##
+$(OBJDIR)/%$(OBJSUFFIX): %.cu
+	$(CUDACC) -c $(CUCCOPT) -o $@  $<
+
 ##  Compile .cpp files ##
 $(OBJDIR)/%$(OBJSUFFIX): %.cpp
 	@$(ECHO) Building $@
@@ -205,16 +241,21 @@ $(OBJDIR)/%$(OBJSUFFIX): %.c
 	xxd -i $(CLPROGRAM).cl | sed 's/\([0-9a-f]\)$$/\0, 0x00/' > $(CLPROGRAM).clh
 
 ##  Link  ##
-$(BINDIR)/$(BINARY): makedirs $(CLSOURCE) $(OBJS)
+$(BINDIR)/$(BINARY): makedirs $(CLSOURCE) $(ZMATLIB) $(OBJS)
 	@$(ECHO) Building $@
 	$(AR)  $(ARFLAGS) $(AROUTPUT) $@ $(OBJS) $(USERARFLAGS) $(EXTRALIB)
 
+
+$(ZMATLIB):
+	-$(MAKE) -C zmat lib AR=ar CPPOPT="$(DLLFLAG)" CCOPT="$(DLLFLAG)" USERLINKOPT=
+
 ##  Documentation  ##
 doc: makedocdir
 	$(DOXY) $(DOXYCFG)
 
 ## Clean
 clean:
+	-$(MAKE) -C zmat clean
 	rm -rf $(OBJS) $(OBJDIR) $(BINDIR) $(DOCDIR)
 ifdef SUBDIRS
 	for i in $(SUBDIRS); do $(MAKE) --no-print-directory -C $$i clean; done

deploy/buildmmc.sh

@@ -28,17 +28,17 @@ then
 fi
 DATE=`date +'%Y%m%d'`
 BUILDROOT=~/space/autobuild/$BUILD/mmc
-OS=`uname -s`
+OSID=`uname -s`
 MACHINE=`uname -m`
 
-if [ "$OS" == "Linux" ]
+if [ "$OSID" == "Linux" ]
 then
     OS=linux
     source ~/.bashrc
-elif [ "$OS" == "Darwin" ]; then
+elif [ "$OSID" == "Darwin" ]; then
     OS=osx
     source ~/.bash_profile
-elif [[ "$OS" == CYGWIN* ]] || [[ "$OS" == MINGW* ]] || [[ "$OS" == MSYS* ]]; then
+elif [[ "$OSID" == CYGWIN* ]] || [[ "$OSID" == MINGW* ]] || [[ "$OSID" == MSYS* ]]; then
     OS=win
 fi
 
@@ -65,6 +65,10 @@ cat <<EOF >> mmc/.git/config
 EOF
 
 cd mmc
+if [ ! -z "$2" ]
+then
+	git checkout $2
+fi
 rm -rf *
 git checkout *
 rm -rf .git
@@ -80,12 +84,27 @@ cd mmc/src
 
 rm -rf ../mmclab/AUTO_BUILD_*
 make clean
-make mex MEXLINKOPT="-static-libstdc++ -static-libgcc -fopenmp" &> ../mmclab/AUTO_BUILD_${DATE}.log
+if [ "$OS" == "win" ]
+then
+    if [[ "$OSID" == CYGWIN* ]]; then
+        PATH=/cygdrive/c/ProgramData/MATLAB/SupportPackages/R2017b/3P.instrset/mingw_w64.instrset/bin/:$PATH  make mex CC=gcc &> ../mmclab/AUTO_BUILD_${DATE}.log
+    else
+        PATH=/c/ProgramData/MATLAB/SupportPackages/R2017b/3P.instrset/mingw_w64.instrset/bin/:$PATH  make mex CC=gcc &> ../mmclab/AUTO_BUILD_${DATE}.log
+    fi
+else
+    make mex MEXLINKOPT="-static-libstdc++ -static-libgcc -fopenmp" &> ../mmclab/AUTO_BUILD_${DATE}.log
+fi
 
 make clean
 if [ "$OS" == "osx" ]
 then
 	make oct USEROCTOPT="CXXFLAGS='-pipe -Os -arch x86_64' DL_LD=g++ DL_LDFLAGS='-fopenmp -static-libgcc -static-libstdc++'" >>  ../mmclab/AUTO_BUILD_${DATE}.log 2>&1
+elif [ "$OS" == "win" ]; then
+    if [[ "$OSID" == CYGWIN* ]]; then
+        PATH=/cygdrive/c/Octave/Octave-4.2.1/bin:$PATH make oct CC=gcc LIBOPENCL='C:\Windows\System32\OpenCL.dll' MEXLINKOPT='-LC:\Octave\Octave-4.2.1\lib64 -LC:\Octave\Octave-4.2.1\lib\octave\4.2.1' &> ../mmclab/AUTO_BUILD_${DATE}.log
+    else
+        PATH=/c/Octave/Octave-4.2.1/bin:$PATH make oct CC=gcc LIBOPENCL='C:\Windows\System32\OpenCL.dll' MEXLINKOPT='-LC:\Octave\Octave-4.2.1\lib64 -LC:\Octave\Octave-4.2.1\lib\octave\4.2.1' &> ../mmclab/AUTO_BUILD_${DATE}.log
+    fi
 else
 	make oct  >>  ../mmclab/AUTO_BUILD_${DATE}.log 2>&1
 fi
@@ -102,7 +121,7 @@ then
 	rm -rf ../mmclab/AUTO_BUILD_${DATE}.log
 fi
 
-cp $BUILDROOT/dlls/*.dll ../mmclab
+#cp $BUILDROOT/dlls/*.dll ../mmclab
 cd ..
 zip -FSr $BUILDROOT/mmclab-${TAG}.zip mmclab
 cd src
@@ -113,8 +132,10 @@ make clean
 if [ "$OS" == "osx" ]
 then
 	make &> $BUILDROOT/mmc_buildlog_${DATE}.log
+elif [ "$OS" == "win" ]; then
+        make USERARFLAGS="-lwinmm -static-libgcc -static"
 else
-	make EXTRALIB="-static -lm" &> $BUILDROOT/mmc_buildlog_${DATE}.log
+	make AR=c++ EXTRALIB="-static-libstdc++ -static-libgcc -lOpenCL -lm" &> $BUILDROOT/mmc_buildlog_${DATE}.log
 fi
 
 if [ -f "bin/mmc" ]
@@ -130,6 +151,7 @@ fi
 cd ../
 rm -rf .git mmclab webmmc commons
 mv src/bin .
+cp $BUILDROOT/dlls/*.dll bin
 rm -rf src .git_filters .gitattributes
 mkdir -p src/bin
 cd src/bin
@@ -148,7 +170,12 @@ then
 	rm -rf mmc/AUTO_BUILD_${DATE}.log
 fi
 
-zip -FSr mmc-${TAG}.zip mmc
+if [ "$OS" == "win" ]
+then
+   zip -FSr mmc-${TAG}.zip mmc
+else
+   zip -FSry mmc-${TAG}.zip mmc
+fi
 
 #mv mmc-${TAG}.zip $BUILDROOT
 

mmclab/example/demo_example_replay.m

@@ -27,6 +27,7 @@
 cfg.tstep=5e-9;
 cfg.prop=[0 0 1 1;0.005 1.0 0.01 1.0];
 cfg.debuglevel='TP';
+cfg.method='elem';
 cfg.isreflect=0;
 cfg.detpos=[30. 20. 0. 1.
    30. 40. 0. 1.
@@ -44,16 +45,18 @@
 
 newcfg.replaydet=1;     % replay photons detected by det#1
 newcfg.seed=seeds.data(:,detp.data(1,:)==newcfg.replaydet);
+detp.ppath=detp.ppath(detp.data(1,:)==newcfg.replaydet,:);
+detp.w0=detp.w0(detp.data(1,:)==newcfg.replaydet,:);
 detp.data=detp.data(:,detp.data(1,:)==newcfg.replaydet);
 % calculate the detected photon weight using the partial path output and prop
-newcfg.replayweight=mmcdetweight(detp.data,newcfg.prop,detp.data(end,:));
-newcfg.replaytime=mmcdettime(detp.data,newcfg.prop);
+newcfg.replayweight=mmcdetweight(detp,newcfg.prop);
+newcfg.replaytime=mmcdettime(detp,newcfg.prop);
 newcfg.isnormalized=0;
 newcfg.outputtype='wl';    % replay and get wl
 
 % now replay the detected photons
 
-[cube2,detp2,~,~]=mmclab(newcfg);
+[cube2,detp2]=mmclab(newcfg);
 
 % the two detected photon arrays should be the same. however, because
 % the program uses multi-threading, the orders may be different

mmclab/example/demo_mmclab_slit.m

@@ -8,7 +8,7 @@
 
 % prepare simulation input
 
-cfg.nphoton=1e6;
+cfg.nphoton=1e7;
 [cfg.node face cfg.elem]=meshabox([0 0 0],[60 60 30],6);
 cfg.elemprop=ones(size(cfg.elem,1),1);
 cfg.srcpos=[20 30 -10];

mmclab/example/test_albedo.m

@@ -23,6 +23,7 @@
 cfg.issaveexit=1;  % save detected photon exit position and angles
 cfg.issaveseed=1;  % save detected photon seeds to replay later
 cfg.isnormalized=1;
+cfg.method='elem';
 
 % albedo method
 cfg.mcmethod=1; 

mmclab/mmclab.m

@@ -386,9 +386,9 @@
 
 if(useopencl==1)
     if(isfield(cfg,'gpuid') && ~ischar(cfg.gpuid) && cfg.gpuid<-1)
-	cfg.gpuid=1;
+	    cfg.gpuid=1;
     end
-    [varargout{1:nargout}]=mmc(cfg);
+    [varargout{1:mmcout}]=mmc(cfg);
 elseif(length(varargin)<2)
     [varargout{1:mmcout}]=mmc(cfg);
 elseif(strcmp(type,'omp'))
@@ -423,9 +423,9 @@
             end
             newdetp.w0=detp(end,:)';  % last column is the initial packet weight
             newdetp.prop=cfg(i).prop;
-	    if(isfield(cfg(i),'unitinmm'))
-		newdetp.unitinmm=cfg(i).unitinmm;
-	    end
+            if(isfield(cfg(i),'unitinmm'))
+                newdetp.unitinmm=cfg(i).unitinmm;
+            end
             newdetp.data=detp;      % enable this line for compatibility
             newdetpstruct(i)=newdetp;
         else

src/Makefile

@@ -1,7 +1,7 @@
 ROOTDIR = ..
 BINARY=mmc
 
-FILES=xorshift128p_rand simpmesh tettracing mcx_utils tictoc mmc cjson/cJSON mmc_host highordermesh mmc_cl_utils mmc_cl_host
+FILES=xorshift128p_rand simpmesh tettracing mcx_utils tictoc mmc cjson/cJSON mmc_host highordermesh mmc_cl_utils mmc_cl_host ubj/ubjw
 CLPROGRAM=mmc_core
 
 ifeq ($(findstring Darwin,$(PLATFORM)), Darwin)
@@ -15,5 +15,12 @@ USERCCFLAGS=-DUSE_OS_TIMER -DUSE_OPENCL -DMMC_XORSHIFT
 
 DUMMY:=$(shell mkdir -p built/cjson)
 
+ifeq ($(MAKECMDGOALS),cuda)
+    FILES+=mmc_cu_host
+    USERCCFLAGS+=-DUSE_CUDA
+    CUCCOPT= -DUSE_ATOMIC -DMCX_SAVE_DETECTORS -DMCX_DO_REFLECTION -DUSE_DMMC -DUSE_BLBADOUEL
+    EXTRALIB+=-lcudart
+endif
+
 include $(ROOTDIR)/commons/Makefile_common.mk
 

src/mcx_const.h

@@ -4,6 +4,8 @@
 #define DET_MASK           0xFFFF0000
 #define MED_MASK           0x0000FFFF
 
+#define MAX_PROP           4000
+
 #define MCX_DEBUG_REC_LEN  6  /**<  number of floating points per position saved when -D M is used for trajectory */
 
 #define MCX_SRC_PENCIL     0  /**<  default-Pencil beam src, no param */
@@ -23,4 +25,28 @@
 #define MCX_SRC_PENCILARRAY 14 /**<  a rectangular array of pencil beams */
 #define MCX_SRC_PATTERN3D  15  /**<  a 3D pattern source, starting from srcpos, srcparam1.{x,y,z} define the x/y/z dimensions */
 
+#define SAVE_DETID(a)         ((a)    & 0x1)   /**<  mask to save detector ID*/
+#define SAVE_NSCAT(a)         ((a)>>1 & 0x1)   /**<  output partial scattering counts */
+#define SAVE_PPATH(a)         ((a)>>2 & 0x1)   /**<  output partial path */
+#define SAVE_MOM(a)           ((a)>>3 & 0x1)   /**<  output momentum transfer */
+#define SAVE_PEXIT(a)         ((a)>>4 & 0x1)   /**<  save exit positions */
+#define SAVE_VEXIT(a)         ((a)>>5 & 0x1)   /**<  save exit vector/directions */
+#define SAVE_W0(a)            ((a)>>6 & 0x1)   /**<  save initial weight */
+
+#define SET_SAVE_DETID(a)     ((a) | 0x1   )   /**<  mask to save detector ID*/
+#define SET_SAVE_NSCAT(a)     ((a) | 0x1<<1)   /**<  output partial scattering counts */
+#define SET_SAVE_PPATH(a)     ((a) | 0x1<<2)   /**<  output partial path */
+#define SET_SAVE_MOM(a)       ((a) | 0x1<<3)   /**<  output momentum transfer */
+#define SET_SAVE_PEXIT(a)     ((a) | 0x1<<4)   /**<  save exit positions */
+#define SET_SAVE_VEXIT(a)     ((a) | 0x1<<5)   /**<  save exit vector/directions */
+#define SET_SAVE_W0(a)        ((a) | 0x1<<6)   /**<  save initial weight */
+
+#define UNSET_SAVE_DETID(a)     ((a) & ~(0x1)   )   /**<  mask to save detector ID*/
+#define UNSET_SAVE_NSCAT(a)     ((a) & ~(0x1<<1))   /**<  output partial scattering counts */
+#define UNSET_SAVE_PPATH(a)     ((a) & ~(0x1<<2))   /**<  output partial path */
+#define UNSET_SAVE_MOM(a)       ((a) & ~(0x1<<3))   /**<  output momentum transfer */
+#define UNSET_SAVE_PEXIT(a)     ((a) & ~(0x1<<4))   /**<  save exit positions */
+#define UNSET_SAVE_VEXIT(a)     ((a) & ~(0x1<<5))   /**<  save exit vector/directions */
+#define UNSET_SAVE_W0(a)        ((a) & ~(0x1<<6))   /**<  save initial weight */
+
 #endif