User: hinder Date: 2013/10/09 12:55 PM

Modified: /trunk/src/ utils.cc

Log: utils.cc: Treat all integration methods uniformly

All tests pass

File Changes:

Directory: /trunk/src/ ======================

File [modified]: utils.cc Delta lines: +10 -19 =================================================================== --- trunk/src/utils.cc 2013-10-09 17:55:39 UTC (rev 92) +++ trunk/src/utils.cc 2013-10-09 17:55:48 UTC (rev 93) @@ -345,25 +345,10 @@ iu = Multipole_Index(0,1,ntheta); CCTK_REAL dph = ph[iu] - ph[il];

- if (CCTK_Equals(integration_method, "midpoint")) + if (CCTK_Equals(integration_method, "Simpson") || + CCTK_Equals(integration_method, "DriscollHealy") || + CCTK_Equals(integration_method, "midpoint")) { - CCTK_REAL tempr = 0.0; - CCTK_REAL tempi = 0.0; - - for (int i=0; i < array_size; i++) { - // the below calculations take the integral of conj(array1)*array2 - tempr += ( array1r[i]*array2r[i] + array1i[i]*array2i[i] ) - *sin(th[i])*dth*dph; - tempi += ( array1r[i]*array2i[i] - array1i[i]*array2r[i] ) - *sin(th[i])*dth*dph; - - *outre = tempr; - *outim = tempi; - } - } - else if (CCTK_Equals(integration_method, "Simpson") || - CCTK_Equals(integration_method, "DriscollHealy")) - { static CCTK_REAL *fr = 0; static CCTK_REAL *fi = 0; static bool allocated_memory = false; @@ -376,6 +361,7 @@ allocated_memory = true; }

+ // the below calculations take the integral of conj(array1)*array2*sin(th) for (int i = 0; i < array_size; i++) { fr[i] = (array1r[i] * array2r[i] + @@ -384,8 +370,13 @@ array1i[i] * array2r[i] ) * sin(th[i]); }

- if (CCTK_Equals(integration_method, "Simpson")) + if (CCTK_Equals(integration_method, "midpoint")) { + *outre = Midpoint2DIntegral(fr, ntheta, nphi, dth, dph); + *outim = Midpoint2DIntegral(fi, ntheta, nphi, dth, dph); + } + else if (CCTK_Equals(integration_method, "Simpson")) + { if (nphi % 2 != 0 || ntheta % 2 != 0) { CCTK_WARN (CCTK_WARN_ABORT, "The Simpson integration method requires even ntheta and even nphi");