/******************************************************************************* * * McCode, neutron/xray ray-tracing package * Copyright (C) 1997-2009, All rights reserved * Risoe National Laboratory, Roskilde, Denmark * Institut Laue Langevin, Grenoble, France * * Runtime: share/mccode-r.h * * %Identification * Written by: KN * Date: Aug 29, 1997 * Release: McStas 3.6.14 * Version: $Revision$ * * Runtime system header for McStas/McXtrace. * * In order to use this library as an external library, the following variables * and macros must be declared (see details in the code) * * struct mcinputtable_struct mcinputtable[]; * int numipar; * metadata_table_t metadata_table[]; * int num_metadata; * char instrument_name[], instrument_source[]; * int traceenabled, defaultmain; * extern MCNUM mccomp_storein[]; * extern MCNUM mcAbsorbProp[]; * extern MCNUM mcScattered; * #define MCCODE_STRING "the McStas/McXtrace version" * * Usage: Automatically embbeded in the c code. * * $Id$ * *******************************************************************************/ #ifndef MCCODE_R_H #define MCCODE_R_H "$Revision$" #include #include #include #include #include #include #include #ifndef _MSC_EXTENSIONS #include #endif #include #include #include #ifdef OPENACC #include #ifndef GCCOFFLOAD #include #else #include #endif #pragma acc routine int noprintf(); #pragma acc routine size_t str_len(const char *s); #else #include #endif /* In case of gcc / clang, ensure to use the built-in isnan/isinf functions */ #if defined(__GNUC__) || defined(__clang__) # ifdef isnan # undef isnan # endif # ifdef isinf # undef isinf # endif # define isnan(x) __builtin_isnan(x) # define isinf(x) __builtin_isinf(x) #endif #ifdef _MSC_EXTENSIONS #ifndef _TIMES_H #define _TIMES_H #if defined(WIN32) || defined(_WIN32) #include #include #include int gettimeofday(struct timeval* t,void* timezone); #define __need_clock_t #include /* Structure describing CPU time used by a process and its children. */ struct tms { clock_t tms_utime; /* User CPU time. */ clock_t tms_stime; /* System CPU time. */ clock_t tms_cutime; /* User CPU time of dead children. */ clock_t tms_cstime; /* System CPU time of dead children. */ }; /* Store the CPU time used by this process and all its dead children (and their dead children) in BUFFER. Return the elapsed real time, or (clock_t) -1 for errors. All times are in CLK_TCKths of a second. */ clock_t times (struct tms *__buffer); typedef long long suseconds_t ; int gettimeofday(struct timeval* t,void* timezone) { struct _timeb timebuffer; _ftime( &timebuffer ); t->tv_sec=timebuffer.time; t->tv_usec=1000*timebuffer.millitm; return 0; } clock_t times (struct tms *__buffer) { __buffer->tms_utime = clock(); __buffer->tms_stime = 0; __buffer->tms_cstime = 0; __buffer->tms_cutime = 0; return __buffer->tms_utime; } #endif #endif #endif /* If the runtime is embedded in the simulation program, some definitions can be made static. */ #ifdef MC_EMBEDDED_RUNTIME # define mcstatic #else # define mcstatic #endif #ifdef __dest_os # if (__dest_os == __mac_os) # define MAC # endif #endif #ifdef __FreeBSD__ # define NEED_STAT_H #endif #if defined(__APPLE__) && defined(__GNUC__) # define NEED_STAT_H #endif #if defined(WIN32) || defined(_WIN32) # define NEED_STAT_H # define NEED_TYPES_H #endif #ifdef NEED_STAT_H # include #endif #ifdef NEED_TYPES_H # include #endif #ifndef MC_PATHSEP_C #if defined(WIN32) || defined(_WIN32) # define MC_PATHSEP_C '\\' # define MC_PATHSEP_S "\\" # else /* !WIN32 */ # define MC_PATHSEP_C '/' # define MC_PATHSEP_S "/" # endif /* !WIN32 */ #endif /* MC_PATHSEP_C */ #if defined(WIN32) || defined(_WIN32) #if defined _MSC_VER #include #elif defined __GNUC__ #include #include #include #endif #define mkdir(a,b) mkdir(a) #define getpid() _getpid() #endif /* the version string is replaced when building distribution with mkdist */ #ifndef MCCODE_STRING # define MCCODE_STRING "McStas 3.6.14 - Apr. 12, 2026" #endif #ifndef MCCODE_DATE # define MCCODE_DATE "Apr. 12, 2026" #endif #ifndef MCCODE_VERSION # define MCCODE_VERSION "3.6.14" #endif #ifndef __MCCODE_VERSION__ #define __MCCODE_VERSION__ L #endif #ifndef MCCODE_NAME # define MCCODE_NAME "McStas" #endif #ifndef MCCODE_PARTICLE # define MCCODE_PARTICLE "neutron" #endif #ifndef MCCODE_PARTICLE_CODE # define MCCODE_PARTICLE_CODE 2112 #endif #ifndef MCCODE_LIBENV # define MCCODE_LIBENV "MCSTAS" #endif #ifndef FLAVOR_UPPER # define FLAVOR_UPPER MCCODE_NAME #endif #ifdef MC_PORTABLE # ifndef NOSIGNALS # define NOSIGNALS 1 # endif #endif #ifdef MAC # ifndef NOSIGNALS # define NOSIGNALS 1 # endif #endif #if (USE_MPI == 0) # undef USE_MPI #endif #ifdef USE_MPI /* default is to disable signals with MPI, as MPICH uses them to communicate */ # ifndef NOSIGNALS # define NOSIGNALS 1 # endif #endif #ifdef OPENACC /* default is to disable signals with PGI/OpenACC */ # ifndef NOSIGNALS # define NOSIGNALS 1 # endif #endif #ifndef OPENACC # ifndef USE_OFF /* default is to enable OFF when not using PGI/OpenACC */ # define USE_OFF # endif # ifndef CPUFUNNEL /* allow to enable FUNNEL-mode on CPU */ # ifdef FUNNEL /* by default disable FUNNEL-mode when not using PGI/OpenACC */ # undef FUNNEL # endif # endif #endif #if (NOSIGNALS == 0) # undef NOSIGNALS #endif /** Header information for metadata-r.c ----------------------------------------------------------------------------- */ struct metadata_table_struct { /* stores metadata strings from components */ char * source; // component name which provided the metadata char * name; // the name of the metadata char * type; // the MIME type of the metadata (free form, valid identifier) char * value; // the metadata string contents }; typedef struct metadata_table_struct metadata_table_t; char * metadata_table_key_component(char* key); char * metadata_table_key_literal(char * key); int metadata_table_defined(int, metadata_table_t *, char *); char * metadata_table_name(int, metadata_table_t *, char *); char * metadata_table_type(int, metadata_table_t *, char *); char * metadata_table_literal(int, metadata_table_t *, char *); void metadata_table_print_all_keys(int no, metadata_table_t * tab); int metadata_table_print_all_components(int no, metadata_table_t * tab); int metadata_table_print_component_keys(int no, metadata_table_t * tab, char * key); /* -------------------------------------------------------------------------- Header information for metadata-r.c --- */ /* Note: the enum instr_formal_types definition MUST be kept synchronized with the one in mccode.h and with the instr_formal_type_names array in cogen.c. */ enum instr_formal_types { instr_type_int, instr_type_string, instr_type_char, instr_type_vector, instr_type_double }; struct mcinputtable_struct { /* defines instrument parameters */ char *name; /* name of parameter */ void *par; /* pointer to instrument parameter (variable) */ enum instr_formal_types type; char *val; /* default value */ char *unit; /* expected unit for parameter; informational only */ }; #ifndef MCCODE_BASE_TYPES typedef double MCNUM; typedef struct {MCNUM x, y, z;} Coords; typedef MCNUM Rotation[3][3]; #endif /* the following variables are defined in the McStas generated C code but should be defined externally in case of independent library usage */ #ifndef DANSE extern struct mcinputtable_struct mcinputtable[]; /* list of instrument parameters */ extern int numipar; /* number of instrument parameters */ extern metadata_table_t metadata_table[]; /* list of component-defined string metadata */ extern int num_metadata; /* number of component-defined string metadata */ extern char instrument_name[], instrument_source[]; /* instrument name and filename */ extern char *instrument_exe; /* executable path = argv[0] or NULL */ extern char instrument_code[]; /* contains the initial 'instr' file */ #ifndef MC_ANCIENT_COMPATIBILITY extern int traceenabled, defaultmain; #endif #endif /* Useful macros ============================================================ */ /* SECTION: Dynamic Arrays */ typedef int* IArray1d; IArray1d create_iarr1d(int n); void destroy_iarr1d(IArray1d a); typedef int** IArray2d; IArray2d create_iarr2d(int nx, int ny); void destroy_iarr2d(IArray2d a); typedef int*** IArray3d; IArray3d create_iarr3d(int nx, int ny, int nz); void destroy_iarr3d(IArray3d a); typedef double* DArray1d; DArray1d create_darr1d(int n); void destroy_darr1d(DArray1d a); typedef double** DArray2d; DArray2d create_darr2d(int nx, int ny); void destroy_darr2d(DArray2d a); typedef double*** DArray3d; DArray3d create_darr3d(int nx, int ny, int nz); void destroy_darr3d(DArray3d a); /* MPI stuff */ #ifdef USE_MPI #include "mpi.h" #ifdef OMPI_MPI_H /* openmpi does not use signals: we may install our sighandler */ #ifndef OPENACC /* ... but only if we are not also running on GPU */ #undef NOSIGNALS #endif #endif /* * MPI_MASTER(i): * execution of i only on master node */ #define MPI_MASTER(statement) { \ if(mpi_node_rank == mpi_node_root)\ { statement; } \ } #ifndef MPI_REDUCE_BLOCKSIZE #define MPI_REDUCE_BLOCKSIZE 100000 #endif int mc_MPI_Sum(double* buf, long count); int mc_MPI_Send(void *sbuf, long count, MPI_Datatype dtype, int dest); int mc_MPI_Recv(void *rbuf, long count, MPI_Datatype dtype, int source); /* MPI_Finalize exits gracefully and should be preferred to MPI_Abort */ #define exit(code) do { \ MPI_Finalize(); \ exit(code); \ } while(0) #else /* !USE_MPI */ #define MPI_MASTER(instr) instr #endif /* USE_MPI */ #ifdef USE_MPI static int mpi_node_count; #endif #ifdef USE_THREADS /* user want threads */ #error Threading (USE_THREADS) support has been removed for very poor efficiency. Use MPI/SSH grid instead. #endif void mcset_ncount(unsigned long long count); /* wrapper to get mcncount */ #pragma acc routine unsigned long long int mcget_ncount(void); /* wrapper to set mcncount */ unsigned long long mcget_run_num(void); /* wrapper to get mcrun_num=0:mcncount-1 */ /* Following part is only embedded when not redundant with mccode.h ========= */ #ifndef MCCODE_H #ifndef NOSIGNALS #include char *mcsig_message; #define SIG_MESSAGE(msg) mcsig_message=(char *)(msg); #else #define SIG_MESSAGE(...) #endif /* !NOSIGNALS */ /* Useful macros and constants ============================================== */ #ifndef FLT_MAX #define FLT_MAX 3.40282347E+38F /* max decimal value of a "float" */ #endif #ifndef MIN #define MIN(a, b) (((a) < (b)) ? (a) : (b)) #endif #ifndef MAX #define MAX(a, b) (((a) > (b)) ? (a) : (b)) #endif #ifndef SQR #define SQR(x) ( (x) * (x) ) #endif #ifndef SIGN #define SIGN(x) (((x)>0.0)?(1):(-1)) #endif # ifndef M_E # define M_E 2.71828182845904523536 // e # endif # ifndef M_LOG2E # define M_LOG2E 1.44269504088896340736 // log2(e) # endif # ifndef M_LOG10E # define M_LOG10E 0.434294481903251827651 // log10(e) # endif # ifndef M_LN2 # define M_LN2 0.693147180559945309417 // ln(2) # endif # ifndef M_LN10 # define M_LN10 2.30258509299404568402 // ln(10) # endif # ifndef M_PI # define M_PI 3.14159265358979323846 // pi # endif # ifndef PI # define PI M_PI // pi - also used in some places # endif # ifndef M_PI_2 # define M_PI_2 1.57079632679489661923 // pi/2 # endif # ifndef M_PI_4 # define M_PI_4 0.785398163397448309616 // pi/4 # endif # ifndef M_1_PI # define M_1_PI 0.318309886183790671538 // 1/pi # endif # ifndef M_2_PI # define M_2_PI 0.636619772367581343076 // 2/pi # endif # ifndef M_2_SQRTPI # define M_2_SQRTPI 1.12837916709551257390 // 2/sqrt(pi) # endif # ifndef M_SQRT2 # define M_SQRT2 1.41421356237309504880 // sqrt(2) # endif # ifndef M_SQRT1_2 # define M_SQRT1_2 0.707106781186547524401 // 1/sqrt(2) # endif #define RAD2MIN ((180*60)/PI) #define MIN2RAD (PI/(180*60)) #define DEG2RAD (PI/180) #define RAD2DEG (180/PI) #define FWHM2RMS 0.424660900144 /* Convert between full-width-half-max and */ #define RMS2FWHM 2.35482004503 /* root-mean-square (standard deviation) */ #define HBAR 1.05457168e-34 /* [Js] h bar Planck constant CODATA 2002 */ #define MNEUTRON 1.67492728e-27 /* [kg] mass of neutron CODATA 2002 */ #define GRAVITY 9.81 /* [m/s^2] gravitational acceleration */ #define NA 6.02214179e23 /* [#atoms/g .mole] Avogadro's number*/ #define UNSET nan("0x6E6F74736574") int nans_match(double, double); int is_unset(double); int is_valid(double); int is_set(double); int all_unset(int n, ...); int all_set(int n, ...); int any_unset(int n, ...); int any_set(int n, ...); /* wrapper to get absolute and relative position of comp */ /* mccomp_posa and mccomp_posr are defined in McStas generated C code */ #define POS_A_COMP_INDEX(index) (instrument->_position_absolute[index]) #define POS_R_COMP_INDEX(index) (instrument->_position_relative[index]) /* setting parameters based COMP_GETPAR (returned as pointer) */ /* compname must be given as a string, type and par are symbols. */ #define COMP_GETPAR3(type, compname, par) \ &( ((_class_ ## type ##_parameters *) _getvar_parameters(compname))->par ) /* the body of this function depends on component instances, and is cogen'd */ void* _getvar_parameters(char* compname); int _getcomp_index(char* compname); /* Note: The two-stage approach to COMP_GETPAR is NOT redundant; without it, * after #define C sample, COMP_GETPAR(C,x) would refer to component C, not to * component sample. Such are the joys of ANSI C. * Anyway the usage of COMP_GETPAR requires that we use sometimes bare names... * NOTE: This can ONLY be used in instrument descriptions, not components. */ #define COMP_GETPAR2(comp, par) (_ ## comp ## _var._parameters.par) #define COMP_GETPAR(comp, par) COMP_GETPAR2(comp,par) #define INSTRUMENT_GETPAR(par) (_instrument_var._parameters.par) /* Current component name, index, position and orientation */ /* These macros work because, using class-based functions, "comp" is usually * the local variable of the active/current component. */ #define INDEX_CURRENT_COMP (_comp->_index) #define NAME_CURRENT_COMP (_comp->_name) #define TYPE_CURRENT_COMP (_comp->_type) #define POS_A_CURRENT_COMP (_comp->_position_absolute) #define POS_R_CURRENT_COMP (_comp->_position_relative) #define ROT_A_CURRENT_COMP (_comp->_rotation_absolute) #define ROT_R_CURRENT_COMP (_comp->_rotation_relative) #define NAME_INSTRUMENT (instrument->_name) /* MCDISPLAY/trace and debugging message sent to stdout */ #ifdef MC_TRACE_ENABLED #define DEBUG #endif #ifdef DEBUG #define DEBUG_INSTR() if(!mcdotrace); else { printf("INSTRUMENT:\n"); printf("Instrument '%s' (%s)\n", instrument_name, instrument_source); } #define DEBUG_COMPONENT(name,c,t) if(!mcdotrace); else {\ printf("COMPONENT: \"%s\"\n" \ "POS: %g, %g, %g, %g, %g, %g, %g, %g, %g, %g, %g, %g\n", \ name, c.x, c.y, c.z, t[0][0], t[0][1], t[0][2], \ t[1][0], t[1][1], t[1][2], t[2][0], t[2][1], t[2][2]); \ printf("Component %30s AT (%g,%g,%g)\n", name, c.x, c.y, c.z); } #define DEBUG_INSTR_END() if(!mcdotrace); else printf("INSTRUMENT END:\n"); #define DEBUG_ENTER() if(!mcdotrace); else printf("ENTER:\n"); #define DEBUG_COMP(c) if(!mcdotrace); else printf("COMP: \"%s\"\n", c); #define DEBUG_LEAVE() if(!mcdotrace); else printf("LEAVE:\n"); #define DEBUG_ABSORB() if(!mcdotrace); else printf("ABSORB:\n"); #else #define DEBUG_INSTR() #define DEBUG_COMPONENT(name,c,t) #define DEBUG_INSTR_END() #define DEBUG_ENTER() #define DEBUG_COMP(c) #define DEBUG_LEAVE() #define DEBUG_ABSORB() #endif // mcDEBUG_STATE and mcDEBUG_SCATTER are defined by mcstas-r.h and mcxtrace-r.h #ifdef TEST #define test_printf printf #else #define test_printf while(0) printf #endif /* send MCDISPLAY message to stdout to show gemoetry */ void mcdis_magnify(char *what); void mcdis_line(double x1, double y1, double z1, double x2, double y2, double z2); void mcdis_dashed_line(double x1, double y1, double z1, double x2, double y2, double z2, int n); void mcdis_multiline(int count, ...); void mcdis_rectangle(char* plane, double x, double y, double z, double width, double height); void mcdis_box(double x, double y, double z, double width, double height, double length, double thickness, double nx, double ny, double nz); void mcdis_circle(char *plane, double x, double y, double z, double r); void mcdis_Circle(double x, double y, double z, double r, double nx, double ny, double nz); void mcdis_cylinder( double x, double y, double z, double r, double height, double thickness, double nx, double ny, double nz); void mcdis_cone( double x, double y, double z, double r, double height, double nx, double ny, double nz); void mcdis_sphere(double x, double y, double z, double r); /* random number generation. ================================================ */ #if RNG_ALG == _RNG_ALG_MT // MT (currently not functional for GPU) # define MC_RAND_MAX ((uint32_t)0xffffffffUL) # define RANDSTATE_LEN 1 # define srandom(seed) mt_srandom_empty() # define random() mt_random() # define _random() mt_random() #elif RNG_ALG == _RNG_ALG_KISS // KISS # ifndef UINT64_MAX # define UINT64_MAX ((uint64_t)0xffffffffffffffffULL) # endif # define MC_RAND_MAX UINT64_MAX # define RANDSTATE_LEN 7 # define srandom(seed) kiss_srandom(_particle->randstate, seed) # define random() kiss_random(_particle->randstate) # define _random() kiss_random(state) #endif #pragma acc routine double _randnorm2(randstate_t* state); // Component writer interface #define randnorm() _randnorm2(_particle->randstate) // NOTE: can't use _randnorm on GPU #define rand01() _rand01(_particle->randstate) #define randpm1() _randpm1(_particle->randstate) #define rand0max(p1) _rand0max(p1, _particle->randstate) #define randminmax(p1, p2) _randminmax(p1, p2, _particle->randstate) #define randtriangle() _randtriangle(_particle->randstate) // Mersenne Twister rng uint32_t mt_random(void); void mt_srandom (uint32_t x); void mt_srandom_empty(); // KISS rng #pragma acc routine uint64_t *kiss_srandom(uint64_t state[7], uint64_t seed); #pragma acc routine uint64_t kiss_random(uint64_t state[7]); // Scrambler / hash function #pragma acc routine seq randstate_t _hash(randstate_t x); // internal RNG (transforms) interface #pragma acc routine double _rand01(randstate_t* state); #pragma acc routine double _randpm1(randstate_t* state); #pragma acc routine double _rand0max(double max, randstate_t* state); #pragma acc routine double _randminmax(double min, double max, randstate_t* state); #pragma acc routine double _randtriangle(randstate_t* state); #ifdef USE_OPENCL #include "opencl-lib.h" #include "opencl-lib.c" #endif #ifndef DANSE int init(void); int raytrace(_class_particle*); int save(FILE *); int finally(void); int display(void); #endif /* GPU related algorithms =================================================== */ /* * Divide-and-conquer strategy for parallel sort absorbed last. */ #ifdef FUNNEL long sort_absorb_last(_class_particle* particles, _class_particle* pbuffer, long len, long buffer_len, long flag_split, long* multiplier); #endif long sort_absorb_last_serial(_class_particle* particles, long len); /* simple vector algebra ==================================================== */ #define vec_prod(x, y, z, x1, y1, z1, x2, y2, z2) \ vec_prod_func(&x, &y, &z, x1, y1, z1, x2, y2, z2) #pragma acc routine seq mcstatic void vec_prod_func(double *x, double *y, double *z, double x1, double y1, double z1, double x2, double y2, double z2); #pragma acc routine seq mcstatic double scalar_prod( double x1, double y1, double z1, double x2, double y2, double z2); #pragma acc routine seq mcstatic void norm_func(double *x, double *y, double *z); #define NORM(x,y,z) norm_func(&x, &y, &z) #pragma acc routine seq void normal_vec(double *nx, double *ny, double *nz, double x, double y, double z); /** * Rotate the vector vx,vy,vz psi radians around the vector ax,ay,az * and put the result in x,y,z. */ #define rotate(x, y, z, vx, vy, vz, phi, ax, ay, az) \ do { \ double mcrt_tmpx = (ax), mcrt_tmpy = (ay), mcrt_tmpz = (az); \ double mcrt_vp, mcrt_vpx, mcrt_vpy, mcrt_vpz; \ double mcrt_vnx, mcrt_vny, mcrt_vnz, mcrt_vn1x, mcrt_vn1y, mcrt_vn1z; \ double mcrt_bx, mcrt_by, mcrt_bz; \ double mcrt_cos, mcrt_sin; \ NORM(mcrt_tmpx, mcrt_tmpy, mcrt_tmpz); \ mcrt_vp = scalar_prod((vx), (vy), (vz), mcrt_tmpx, mcrt_tmpy, mcrt_tmpz); \ mcrt_vpx = mcrt_vp*mcrt_tmpx; \ mcrt_vpy = mcrt_vp*mcrt_tmpy; \ mcrt_vpz = mcrt_vp*mcrt_tmpz; \ mcrt_vnx = (vx) - mcrt_vpx; \ mcrt_vny = (vy) - mcrt_vpy; \ mcrt_vnz = (vz) - mcrt_vpz; \ vec_prod(mcrt_bx, mcrt_by, mcrt_bz, \ mcrt_tmpx, mcrt_tmpy, mcrt_tmpz, mcrt_vnx, mcrt_vny, mcrt_vnz); \ mcrt_cos = cos((phi)); mcrt_sin = sin((phi)); \ mcrt_vn1x = mcrt_vnx*mcrt_cos + mcrt_bx*mcrt_sin; \ mcrt_vn1y = mcrt_vny*mcrt_cos + mcrt_by*mcrt_sin; \ mcrt_vn1z = mcrt_vnz*mcrt_cos + mcrt_bz*mcrt_sin; \ (x) = mcrt_vpx + mcrt_vn1x; \ (y) = mcrt_vpy + mcrt_vn1y; \ (z) = mcrt_vpz + mcrt_vn1z; \ } while(0) /** * Mirror (xyz) in the plane given by the point (rx,ry,rz) and normal (nx,ny,nz) * * TODO: This define is seemingly never used... */ #define mirror(x,y,z,rx,ry,rz,nx,ny,nz) \ do { \ double mcrt_tmpx= (nx), mcrt_tmpy = (ny), mcrt_tmpz = (nz); \ double mcrt_tmpt; \ NORM(mcrt_tmpx, mcrt_tmpy, mcrt_tmpz); \ mcrt_tmpt=scalar_prod((rx),(ry),(rz),mcrt_tmpx,mcrt_tmpy,mcrt_tmpz); \ (x) = rx -2 * mcrt_tmpt*mcrt_rmpx; \ (y) = ry -2 * mcrt_tmpt*mcrt_rmpy; \ (z) = rz -2 * mcrt_tmpt*mcrt_rmpz; \ } while (0) #pragma acc routine Coords coords_set(MCNUM x, MCNUM y, MCNUM z); #pragma acc routine Coords coords_get(Coords a, MCNUM *x, MCNUM *y, MCNUM *z); #pragma acc routine Coords coords_add(Coords a, Coords b); #pragma acc routine Coords coords_sub(Coords a, Coords b); #pragma acc routine Coords coords_neg(Coords a); #pragma acc routine Coords coords_scale(Coords b, double scale); #pragma acc routine double coords_sp(Coords a, Coords b); #pragma acc routine Coords coords_xp(Coords b, Coords c); #pragma acc routine double coords_len(Coords a); #pragma acc routine seq void coords_print(Coords a); #pragma acc routine seq mcstatic void coords_norm(Coords* c); #pragma acc routine seq void rot_set_rotation(Rotation t, double phx, double phy, double phz); #pragma acc routine seq int rot_test_identity(Rotation t); #pragma acc routine seq void rot_mul(Rotation t1, Rotation t2, Rotation t3); #pragma acc routine seq void rot_copy(Rotation dest, Rotation src); #pragma acc routine seq void rot_transpose(Rotation src, Rotation dst); #pragma acc routine seq Coords rot_apply(Rotation t, Coords a); #pragma acc routine seq void mccoordschange(Coords a, Rotation t, _class_particle *particle); #pragma acc routine seq void mccoordschange_polarisation(Rotation t, double *sx, double *sy, double *sz); double mcestimate_error(double N, double p1, double p2); void mcreadparams(void); /* this is now in mcstas-r.h and mcxtrace-r.h as the number of state parameters is no longer equal */ _class_particle mcgenstate(void); // trajectory/shape intersection routines #pragma acc routine seq int inside_rectangle(double, double, double, double); #pragma acc routine seq int box_intersect(double *dt_in, double *dt_out, double x, double y, double z, double vx, double vy, double vz, double dx, double dy, double dz); #pragma acc routine seq int cylinder_intersect(double *t0, double *t1, double x, double y, double z, double vx, double vy, double vz, double r, double h); #pragma acc routine seq int sphere_intersect(double *t0, double *t1, double x, double y, double z, double vx, double vy, double vz, double r); // second order equation roots #pragma acc routine seq int solve_2nd_order(double *t1, double *t2, double A, double B, double C); // random vector generation to shape // defines silently introducing _particle as the last argument #define randvec_target_circle(xo, yo, zo, solid_angle, xi, yi, zi, radius) \ _randvec_target_circle(xo, yo, zo, solid_angle, xi, yi, zi, radius, _particle) #define randvec_target_rect_angular(xo, yo, zo, solid_angle, xi, yi, zi, height, width, A) \ _randvec_target_rect_angular(xo, yo, zo, solid_angle, xi, yi, zi, height, width, A, _particle) #define randvec_target_rect_real(xo, yo, zo, solid_angle, xi, yi, zi, height, width, A, lx, ly, lz, order) \ _randvec_target_rect_real(xo, yo, zo, solid_angle, xi, yi, zi, height, width, A, lx, ly, lz, order, _particle) // defines forwarding to "inner" functions #define randvec_target_sphere randvec_target_circle #define randvec_target_rect(p0,p1,p2,p3,p4,p5,p6,p7,p8,p9) \ randvec_target_rect_real(p0,p1,p2,p3,p4,p5,p6,p7,p8,p9,0,0,0,1) // headers for randvec #pragma acc routine seq void _randvec_target_circle(double *xo, double *yo, double *zo, double *solid_angle, double xi, double yi, double zi, double radius, _class_particle* _particle); #pragma acc routine seq void _randvec_target_rect_angular(double *xo, double *yo, double *zo, double *solid_angle, double xi, double yi, double zi, double height, double width, Rotation A, _class_particle* _particle); #pragma acc routine seq void _randvec_target_rect_real(double *xo, double *yo, double *zo, double *solid_angle, double xi, double yi, double zi, double height, double width, Rotation A, double lx, double ly, double lz, int order, _class_particle* _particle); // this is the main() int mccode_main(int argc, char *argv[]); #endif /* !MCCODE_H */ #ifndef MCCODE_R_IO_H #define MCCODE_R_IO_H "$Revision$" #if (USE_NEXUS == 0) #undef USE_NEXUS #endif #ifndef CHAR_BUF_LENGTH #define CHAR_BUF_LENGTH 1024 #endif /* I/O section part ========================================================= */ /* ========================================================================== */ /* MCCODE_R_IO_C */ /* ========================================================================== */ /* main DETECTOR structure which stores most information to write to data files */ struct mcdetector_struct { char filename[CHAR_BUF_LENGTH]; /* file name of monitor */ double Position[3]; /* position of detector component*/ char position[CHAR_BUF_LENGTH]; /* position of detector component (string)*/ Rotation Rotation; /* position of detector component*/ char options[CHAR_BUF_LENGTH]; /* Monitor_nD style list-mode'options' (string)*/ char component[CHAR_BUF_LENGTH]; /* component instance name */ char nexuscomp[CHAR_BUF_LENGTH]; /* component naming in NeXus/HDF case */ char instrument[CHAR_BUF_LENGTH]; /* instrument name */ char type[CHAR_BUF_LENGTH]; /* data type, e.g. 0d, 1d, 2d, 3d */ char user[CHAR_BUF_LENGTH]; /* user name, e.g. HOME */ char date[CHAR_BUF_LENGTH]; /* date of simulation end/write time */ char title[CHAR_BUF_LENGTH]; /* title of detector */ char xlabel[CHAR_BUF_LENGTH]; /* X axis label */ char ylabel[CHAR_BUF_LENGTH]; /* Y axis label */ char zlabel[CHAR_BUF_LENGTH]; /* Z axis label */ char xvar[CHAR_BUF_LENGTH]; /* X variable name */ char yvar[CHAR_BUF_LENGTH]; /* Y variable name */ char zvar[CHAR_BUF_LENGTH]; /* Z variable name */ char ncount[CHAR_BUF_LENGTH]; /* number of events initially generated */ char limits[CHAR_BUF_LENGTH]; /* X Y Z limits, e.g. [xmin xmax ymin ymax zmin zmax] */ char variables[CHAR_BUF_LENGTH]; /* variables written into data block */ char statistics[CHAR_BUF_LENGTH]; /* center, mean and half width along axis */ char signal[CHAR_BUF_LENGTH]; /* min max and mean of signal (data block) */ char values[CHAR_BUF_LENGTH]; /* integrated values e.g. [I I_err N] */ double xmin,xmax; /* min max of axes */ double ymin,ymax; double zmin,zmax; double intensity; /* integrated values for data block */ double error; double events; double min; /* statistics for data block */ double max; double mean; double centerX; /* statistics for axes */ double halfwidthX; double centerY; double halfwidthY; int rank; /* dimensionaly of monitor, e.g. 0 1 2 3 */ char istransposed; /* flag to transpose matrix for some formats */ long m,n,p; /* dimensions of data block and along axes */ long date_l; /* same as date, but in sec since 1970 */ double *p0, *p1, *p2; /* pointers to saved data, NULL when freed */ char format[CHAR_BUF_LENGTH]; /* format for file generation */ }; typedef struct mcdetector_struct MCDETECTOR; static char *dirname = NULL; /* name of output directory */ static char *siminfo_name = "mccode"; /* default output sim file name */ char *mcformat = NULL; /* NULL (default) or a specific format */ /* file I/O definitions and function prototypes */ #ifndef MC_EMBEDDED_RUNTIME /* the mcstatic variables (from mccode-r.c) */ extern FILE * siminfo_file; /* handle to the output siminfo file */ extern int mcgravitation; /* flag to enable gravitation */ extern int mcdotrace; /* flag to print MCDISPLAY messages */ #else mcstatic FILE *siminfo_file = NULL; #endif /* I/O function prototypes ================================================== */ // from msysgit: https://code.google.com/p/msysgit/source/browse/compat/strcasestr.c char *strcasestr(const char *haystack, const char *needle); /* output functions */ MCDETECTOR mcdetector_out_0D(char *t, double p0, double p1, double p2, char *c, Coords pos, Rotation rot, int index); MCDETECTOR mcdetector_out_1D(char *t, char *xl, char *yl, char *xvar, double x1, double x2, long n, double *p0, double *p1, double *p2, char *f, char *c, Coords pos, Rotation rot, int index); MCDETECTOR mcdetector_out_2D(char *t, char *xl, char *yl, double x1, double x2, double y1, double y2, long m, long n, double *p0, double *p1, double *p2, char *f, char *c, Coords pos, Rotation rot, int index); MCDETECTOR mcdetector_out_list(char *t, char *xl, char *yl, long m, long n, double *p1, char *f, char *c, Coords posa, Rotation rot,char* options, int index); /* wrappers to output functions, that automatically set NAME and POSITION */ #define DETECTOR_OUT(p0,p1,p2) mcdetector_out_0D(NAME_CURRENT_COMP,p0,p1,p2,NAME_CURRENT_COMP,POS_A_CURRENT_COMP,ROT_A_CURRENT_COMP,INDEX_CURRENT_COMP) #define DETECTOR_OUT_0D(t,p0,p1,p2) mcdetector_out_0D(t,p0,p1,p2,NAME_CURRENT_COMP,POS_A_CURRENT_COMP,ROT_A_CURRENT_COMP,INDEX_CURRENT_COMP) #define DETECTOR_OUT_1D(t,xl,yl,xvar,x1,x2,n,p0,p1,p2,f) \ mcdetector_out_1D(t,xl,yl,xvar,x1,x2,n,p0,p1,p2,f,NAME_CURRENT_COMP,POS_A_CURRENT_COMP,ROT_A_CURRENT_COMP,INDEX_CURRENT_COMP) #define DETECTOR_OUT_2D(t,xl,yl,x1,x2,y1,y2,m,n,p0,p1,p2,f) \ mcdetector_out_2D(t,xl,yl,x1,x2,y1,y2,m,n,p0,p1,p2,f,NAME_CURRENT_COMP,POS_A_CURRENT_COMP,ROT_A_CURRENT_COMP,INDEX_CURRENT_COMP) #ifdef USE_NEXUS #include "napi.h" NXhandle nxhandle; #endif #endif /* ndef MCCODE_R_IO_H */ #endif /* MCCODE_R_H */