diff options
author | schnetter <schnetter@90eea020-d82d-4da5-bf6e-4ee79ff7632f> | 2003-11-05 20:02:40 +0000 |
---|---|---|
committer | schnetter <schnetter@90eea020-d82d-4da5-bf6e-4ee79ff7632f> | 2003-11-05 20:02:40 +0000 |
commit | b86fe77ff66c9c880fa2dfa012a8dc564a26c838 (patch) | |
tree | 865d024bef260802b4b9993ad8500834748eaf7d /src | |
parent | 5a10b056bae14fb944fdbdf5a025594d1c089fad (diff) |
Undo my accidental changes.
git-svn-id: http://svn.cactuscode.org/arrangements/CactusWave/WaveToyC/trunk@94 90eea020-d82d-4da5-bf6e-4ee79ff7632f
Diffstat (limited to 'src')
-rw-r--r-- | src/WaveToy.c | 208 |
1 files changed, 27 insertions, 181 deletions
diff --git a/src/WaveToy.c b/src/WaveToy.c index 5ecbdd8..7c98b1e 100644 --- a/src/WaveToy.c +++ b/src/WaveToy.c @@ -8,15 +8,8 @@ @version $Header$ @@*/ -#include <assert.h> - -#if 0 -#include <sys/time.h> -#include <sys/resource.h> -#include <unistd.h> -#endif - -#include "cctk.h" + +#include "cctk.h" #include "cctk_Arguments.h" #include "cctk_Parameters.h" @@ -27,44 +20,7 @@ static const char *rcsid = "$Header$"; CCTK_FILEVERSION(CactusWave_WaveToyC_WaveToy_c); void WaveToyC_Boundaries(CCTK_ARGUMENTS); -void WaveToyC_Evolution(CCTK_ARGUMENTS); - -#define DEBUG 0 - - - -static inline int min (int const a, int const b) -{ - return a<b ? a : b; -} - -static inline int max (int const a, int const b) -{ - return a>b ? a : b; -} - -#if 0 -static inline long long timeval2longlong (const struct timeval * const tv) -{ - return tv->tv_sec * 1000000LL + tv->tv_usec; -} - -static long long currenttime (void) -{ - struct rusage ru; - getrusage (RUSAGE_SELF, &ru); - return timeval2longlong (ru.ru_utime) + timeval2longlong (ru.ru_stime); -} -#endif - -static long long currenttime (void) -{ - unsigned long eax, edx; - asm volatile ("rdtsc" : "=a" (eax), "=d" (edx)); - return (unsigned long long)edx << 32 | (unsigned long long)eax; -} - - +void WaveToyC_Evolution(CCTK_ARGUMENTS); /*@@ @routine WaveToyC_Evolution @@ -81,14 +37,14 @@ static long long currenttime (void) void WaveToyC_Evolution(CCTK_ARGUMENTS) { + DECLARE_CCTK_ARGUMENTS; - DECLARE_CCTK_PARAMETERS; int i,j,k; int index; int istart, jstart, kstart, iend, jend, kend; CCTK_REAL dx,dy,dz,dt,dx2,dy2,dz2,dt2; - CCTK_REAL dt2dx2i,dt2dy2i,dt2dz2i; + CCTK_REAL dx2i,dy2i,dz2i; CCTK_REAL factor; @@ -103,9 +59,9 @@ void WaveToyC_Evolution(CCTK_ARGUMENTS) dz2 = dz*dz; dt2 = dt*dt; - dt2dx2i = dt2/dx2; - dt2dy2i = dt2/dy2; - dt2dz2i = dt2/dz2; + dx2i = 1.0/dx2; + dy2i = 1.0/dy2; + dz2i = 1.0/dz2; istart = 1; jstart = 1; @@ -116,139 +72,29 @@ void WaveToyC_Evolution(CCTK_ARGUMENTS) kend = cctk_lsh[2]-1; /* Do the evolution */ - factor = 2 * (1 - (dt2dx2i + dt2dy2i + dt2dz2i)); + factor = 2*(1 - (dt2)*(dx2i + dy2i + dz2i)); + for (k=kstart; k<kend; k++) { - /* user-defined values */ - const char * const loopthorn = CCTK_THORNSTRING; - const char * const loopfile = __FILE__; - int const loopline = __LINE__; - const char * const loopname = "Evolution"; - - int const dim = 3; - - int const imin = istart; - int const jmin = jstart; - int const kmin = kstart; - int const imax = iend; - int const jmax = jend; - int const kmax = kend; - int const ishape = cctk_lsh[0]; - int const jshape = cctk_lsh[1]; - int const kshape = cctk_lsh[2]; - - /* tuning parameters */ - if (DEBUG) printf ("level1_stride = [%d,%d,%d]\n", level1_stride[0], level1_stride[1], level1_stride[2]); - if (DEBUG) printf ("level2_stride = [%d,%d,%d]\n", level2_stride[0], level2_stride[1], level2_stride[2]); - - int const i1step = level1_stride[0]<0 ? imax-imin : level1_stride[0]; - int const j1step = level1_stride[1]<0 ? jmax-jmin : level1_stride[1]; - int const k1step = level1_stride[2]<0 ? kmax-kmin : level1_stride[2]; - - int const i2step = level2_stride[0]<0 ? imax-imin : level2_stride[0]; - int const j2step = level2_stride[1]<0 ? jmax-jmin : level2_stride[1]; - int const k2step = level2_stride[2]<0 ? kmax-kmin : level2_stride[2]; - - /* internal stuff */ - int rep; - for (rep = 0; rep < repetitions; ++rep) { - - long long l0time = 0, l1time = 0, l2time = 0; - long long l0count = 0, l1count = 0, l2count = 0; - - /* level 2 loop */ - l2time = currenttime() - l2time; - ++ l2count; - - int const i2min = imin; - int const j2min = jmin; - int const k2min = kmin; - int const i2max = imax; - int const j2max = jmax; - int const k2max = kmax; - - int i2, j2, k2; - - for (k2 = k2min; k2 < k2max; k2 += k2step) { - for (j2 = j2min; j2 < j2max; j2 += j2step) { - for (i2 = i2min; i2 < i2max; i2 += i2step) { - - /* level 1 loop */ - l1time = currenttime() - l1time; - ++ l1count; - - int const i1min = i2; - int const j1min = j2; - int const k1min = k2; - int const i1max = min (imax, i2 + i2step); - int const j1max = min (jmax, j2 + j2step); - int const k1max = min (kmax, k2 + k2step); - - int i1, j1, k1; - - for (k1 = k1min; k1 < k1max; k1 += k1step) { - for (j1 = j1min; j1 < j1max; j1 += j1step) { - for (i1 = i1min; i1 < i1max; i1 += i1step) { - - /* level 0 loop */ - l0time = currenttime() - l0time; - ++ l0count; - - int const i0min = i1; - int const j0min = j1; - int const k0min = k1; - int const i0max = min (imax, i1 + i1step); - int const j0max = min (jmax, j1 + j1step); - int const k0max = min (kmax, k1 + k1step); - - int i0, j0, k0; - - for (k0 = k0min; k0 < k0max; ++k0) { - for (j0 = j0min; j0 < j0max; ++j0) { - for (i0 = i0min; i0 < i0max; ++i0) { - - int const i = i0; - int const j = j0; - int const k = k0; - - int const di = 1; - int const dj = ishape; - int const dk = ishape * jshape; - - int const index = i *di + j * dj + k * dk; - + for (j=jstart; j<jend; j++) { - phi[index] = factor * phi_p[index] - phi_p_p[index] - + (phi_p[index+di] + phi_p[index-di]) * dt2dx2i - + (phi_p[index+dj] + phi_p[index-dj]) * dt2dy2i - + (phi_p[index+dk] + phi_p[index-dk]) * dt2dz2i; - } - - } } } /* end level 0 loop */ - l0time = currenttime() - l0time; - - } } } /* end level 1 loop */ - l1time = currenttime() - l1time; - - } } } /* end level 2 loop */ - l2time = currenttime() - l2time; - - if (DEBUG) - CCTK_VInfo (CCTK_THORNSTRING, - "Timing information for loop %s::%s\n" - " (file \"%s\", line %d):\n" - " Level 0: shape %d %d %d, count %lld, time %lld\n" - " Level 1: shape %d %d %d, count %lld, time %lld\n" - " Level 2: shape %d %d %d, count %lld, time %lld", - loopthorn, loopname, loopfile, loopline, - i1step , j1step , k1step , l0count, l0time, - i2step , j2step , k2step , l1count, l1time, - imax-imin, jmax-jmin, kmax-kmin, l2count, l2time); - - } /* end repetition loop */ - + for (i=istart; i<iend; i++) + { + index = CCTK_GFINDEX3D(cctkGH,i,j,k); + + phi[index] = factor* + phi_p[index] - phi_p_p[index] + + (dt2) * + ( ( phi_p[CCTK_GFINDEX3D(cctkGH,i+1,j ,k )] + +phi_p[CCTK_GFINDEX3D(cctkGH,i-1,j ,k )] )*dx2i + +( phi_p[CCTK_GFINDEX3D(cctkGH,i ,j+1,k )] + +phi_p[CCTK_GFINDEX3D(cctkGH,i ,j-1,k )] )*dy2i + +( phi_p[CCTK_GFINDEX3D(cctkGH,i ,j ,k+1)] + +phi_p[CCTK_GFINDEX3D(cctkGH,i ,j, k-1)] )*dz2i); + } + } } - + return; } |