From 167b7353013f12c35d935a88a91c7c961ae5a2db Mon Sep 17 00:00:00 2001 From: eschnett Date: Sat, 16 Jun 2012 15:39:07 +0000 Subject: Implement Driscoll&Healy integration Implements a more accurate integration over the sphere using an algorithm by Driscoll & Healy. This algorithm uses Gaussian integration weights, leading (almost) to exponential convergence. git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/Multipole/trunk@83 4f5cb9a8-4dd8-4c2d-9bbd-173fa4467843 --- param.ccl | 5 +- src/integrate.cc | 49 +++++++++- src/integrate.hh | 3 + src/utils.cc | 31 +++++- test/test_driscollhealy.par | 89 ++++++++++++++++++ .../test_driscollhealy/mp_harmonic_im_r8.00.ph.asc | 104 +++++++++++++++++++++ .../test_driscollhealy/mp_harmonic_im_r8.00.th.asc | 55 +++++++++++ .../test_driscollhealy/mp_harmonic_l0_m0_r8.00.asc | 1 + .../mp_harmonic_l1_m-1_r8.00.asc | 1 + .../test_driscollhealy/mp_harmonic_l1_m0_r8.00.asc | 1 + .../test_driscollhealy/mp_harmonic_l1_m1_r8.00.asc | 1 + .../mp_harmonic_l2_m-1_r8.00.asc | 1 + .../mp_harmonic_l2_m-2_r8.00.asc | 1 + .../test_driscollhealy/mp_harmonic_l2_m0_r8.00.asc | 1 + .../test_driscollhealy/mp_harmonic_l2_m1_r8.00.asc | 1 + .../test_driscollhealy/mp_harmonic_l2_m2_r8.00.asc | 1 + .../test_driscollhealy/mp_harmonic_re_r8.00.ph.asc | 104 +++++++++++++++++++++ .../test_driscollhealy/mp_harmonic_re_r8.00.th.asc | 55 +++++++++++ test/test_driscollhealy/test_driscollhealy.par | 90 ++++++++++++++++++ 19 files changed, 588 insertions(+), 6 deletions(-) create mode 100644 test/test_driscollhealy.par create mode 100644 test/test_driscollhealy/mp_harmonic_im_r8.00.ph.asc create mode 100644 test/test_driscollhealy/mp_harmonic_im_r8.00.th.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l0_m0_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l1_m-1_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l1_m0_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l1_m1_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l2_m-1_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l2_m-2_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l2_m0_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l2_m1_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_l2_m2_r8.00.asc create mode 100644 test/test_driscollhealy/mp_harmonic_re_r8.00.ph.asc create mode 100644 test/test_driscollhealy/mp_harmonic_re_r8.00.th.asc create mode 100644 test/test_driscollhealy/test_driscollhealy.par diff --git a/param.ccl b/param.ccl index 9af3650..b8ed861 100644 --- a/param.ccl +++ b/param.ccl @@ -29,8 +29,9 @@ CCTK_STRING coord_system "What is the coord system?" KEYWORD integration_method "How to do surface integrals" STEERABLE=always { - "midpoint" :: "Midpoint rule (1st order)" - "Simpson" :: "Simpson's rule (4th order)" + "midpoint" :: "Midpoint rule (1st order)" + "Simpson" :: "Simpson's rule (4th order) [requires even ntheta and nphi]" + "DriscollHealy" :: "Driscoll & Healy (exponentially convergent) [requires even ntheta]" } "midpoint" CCTK_INT out_every "How often to output" \ diff --git a/src/integrate.cc b/src/integrate.cc index 0f1c392..1a9f44b 100644 --- a/src/integrate.cc +++ b/src/integrate.cc @@ -107,6 +107,52 @@ CCTK_REAL Simpson2DIntegral(CCTK_REAL const *f, int nx, int ny, CCTK_REAL hx, CC return ((double) 1 / (double) 9) * hx * hy * integrand_sum; } +// See: J.R. Driscoll and D.M. Healy Jr., Computing Fourier transforms +// and convolutions on the 2-sphere. Advances in Applied Mathematics, +// 15(2):202–250, 1994. +CCTK_REAL DriscollHealy2DIntegral(CCTK_REAL const *const f, + int const nx, int const ny, + CCTK_REAL const hx, CCTK_REAL const hy) +{ + assert(f); + assert(nx >= 0); + assert(ny >= 0); + assert(nx % 2 == 0); + + CCTK_REAL integrand_sum = 0.0; + + // Skip the poles (ix=0 and ix=nx), since the weight there is zero + // anyway +#pragma omp parallel for reduction(+: integrand_sum) + for (int ix = 1; ix < nx; ++ ix) + { + + // These weights lead to an almost spectral convergence + CCTK_REAL const theta = M_PI * ix / nx; + CCTK_REAL weight = 0.0; + for (int l = 0; l < nx/2; ++ l) + { + weight += sin((2*l+1)*theta)/(2*l+1); + } + weight *= 4.0 / M_PI; + + CCTK_REAL local_sum = 0.0; + // Skip the last point (iy=ny), since we assume periodicity and + // therefor it has the same value as the first point. We don't use + // weights in this direction, which leads to spectral convergence. + // (Yay periodicity!) + for (int iy = 0; iy < ny; ++ iy) + { + local_sum += f[idx(ix,iy)]; + } + + integrand_sum += weight * local_sum; + + } + + return hx * hy * integrand_sum; +} + // 1D integrals static CCTK_REAL Simpson1DIntegral(CCTK_REAL const *f, int n, CCTK_REAL h) @@ -114,7 +160,8 @@ static CCTK_REAL Simpson1DIntegral(CCTK_REAL const *f, int n, CCTK_REAL h) CCTK_REAL integrand_sum = 0; int i = 0; - assert(n > 0); assert(f); + assert(f); + assert(n > 0); assert(n % 2 == 0); int p = n / 2; diff --git a/src/integrate.hh b/src/integrate.hh index c45917e..5bcd43f 100644 --- a/src/integrate.hh +++ b/src/integrate.hh @@ -6,4 +6,7 @@ CCTK_REAL Simpson2DIntegral(CCTK_REAL const *f, int nx, int ny, CCTK_REAL hx, CCTK_REAL hy); +CCTK_REAL DriscollHealy2DIntegral(CCTK_REAL const *f, int nx, int ny, + CCTK_REAL hx, CCTK_REAL hy); + #endif diff --git a/src/utils.cc b/src/utils.cc index 74d9eef..cb25bdf 100644 --- a/src/utils.cc +++ b/src/utils.cc @@ -346,7 +346,8 @@ void Multipole_Integrate(int array_size, int nthetap, *outim = tempi; } } - else if (CCTK_Equals(integration_method, "simpson")) + else if (CCTK_Equals(integration_method, "Simpson") || + CCTK_Equals(integration_method, "DriscollHealy")) { static CCTK_REAL *fr = 0; static CCTK_REAL *fi = 0; @@ -368,7 +369,31 @@ void Multipole_Integrate(int array_size, int nthetap, array1i[i] * array2r[i] ) * sin(th[i]); } - *outre = Simpson2DIntegral(fr, ntheta, nphi, dth, dph); - *outim = Simpson2DIntegral(fi, ntheta, nphi, dth, dph); + 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"); + } + *outre = Simpson2DIntegral(fr, ntheta, nphi, dth, dph); + *outim = Simpson2DIntegral(fi, ntheta, nphi, dth, dph); + } + else if (CCTK_Equals(integration_method, "DriscollHealy")) + { + if (ntheta % 2 != 0) + { + CCTK_WARN (CCTK_WARN_ABORT, "The Driscoll&Healy integration method requires even ntheta"); + } + *outre = DriscollHealy2DIntegral(fr, ntheta, nphi, dth, dph); + *outim = DriscollHealy2DIntegral(fi, ntheta, nphi, dth, dph); + } + else + { + CCTK_WARN(CCTK_WARN_ABORT, "internal error"); + } + } + else + { + CCTK_WARN(CCTK_WARN_ABORT, "internal error"); } } diff --git a/test/test_driscollhealy.par b/test/test_driscollhealy.par new file mode 100644 index 0000000..c1fee6f --- /dev/null +++ b/test/test_driscollhealy.par @@ -0,0 +1,89 @@ + +ActiveThorns = "CoordBase SymBase Boundary CartGrid3d IOUtil Carpet CarpetLib CarpetInterp AEILocalInterp InitBase Multipole LoopControl" + +############################################################# +# Grid +############################################################# + +CartGrid3D::type = "coordbase" +CartGrid3D::domain = "full" +CartGrid3D::avoid_origin = "no" + +CoordBase::domainsize = minmax +CoordBase::xmin = -10 +CoordBase::ymin = -10 +CoordBase::zmin = -10 +CoordBase::xmax = 10 +CoordBase::ymax = 10 +CoordBase::zmax = 10 +CoordBase::dx = 0.2 +CoordBase::dy = 0.2 +CoordBase::dz = 0.2 +CoordBase::boundary_size_x_lower = 2 +CoordBase::boundary_size_y_lower = 2 +CoordBase::boundary_size_z_lower = 2 +CoordBase::boundary_shiftout_x_lower = 0 +CoordBase::boundary_shiftout_y_lower = 0 +CoordBase::boundary_shiftout_z_lower = 0 +CoordBase::boundary_size_x_upper = 2 +CoordBase::boundary_size_y_upper = 2 +CoordBase::boundary_size_z_upper = 2 +CoordBase::boundary_shiftout_x_upper = 0 +CoordBase::boundary_shiftout_y_upper = 0 +CoordBase::boundary_shiftout_z_upper = 0 + +############################################################# +# Carpet +############################################################# + +Carpet::ghost_size = 3 +Carpet::domain_from_coordbase = "yes" +Carpet::poison_new_timelevels = "yes" +Carpet::check_for_poison = "no" +CarpetLib::poison_value = 113 +Carpet::init_fill_timelevels = "yes" + +############################################################# +# CarpetLib +############################################################# + +CarpetLib::poison_new_memory = "yes" + +############################################################# +# Cactus +############################################################# + +Cactus::terminate = "iteration" +Cactus::cctk_itlast = 0 + +############################################################# +# Multipole +############################################################# + +Multipole::nradii = 1 +Multipole::radius[0] = 8.0 +Multipole::variables = "Multipole::harmonic_re{sw=-2 cmplx='Multipole::harmonic_im' name='harmonic'}" +Multipole::integration_method = "DriscollHealy" + +Multipole::enable_test = "yes" +Multipole::test_l = 2 +Multipole::test_m = 2 +Multipole::out_1d_every = 1 + +############################################################# +# Output +############################################################# + +IO::out_dir = $parfile +IO::out_fileinfo = "none" + +# Enabling 1D output for the test grid functions would be helpful for +# localising any failures but it makes the tests dependent on the +# number of processors, as CarpetIOASCII's output is dependent on +# this. + +# CarpetIOASCII::out1d_vars = "Multipole::harmonics" +# CarpetIOASCII::out1d_every = 1 +# CarpetIOASCII::out1d_x = yes +# CarpetIOASCII::out1d_y = yes +# CarpetIOASCII::out1d_z = yes diff --git a/test/test_driscollhealy/mp_harmonic_im_r8.00.ph.asc b/test/test_driscollhealy/mp_harmonic_im_r8.00.ph.asc new file mode 100644 index 0000000..b689c89 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_im_r8.00.ph.asc @@ -0,0 +1,104 @@ +"Time = 0 +0.000000 -6.907573430162997911e-17 +0.062832 0.0597709107978601703 +0.125664 0.1185991090836183148 +0.188496 0.1755571071177308706 +0.251327 0.22974631087048486 +0.314159 0.2803123461157627783 +0.376991 0.326457761226630172 +0.439823 0.3674546220110177508 +0.502655 0.4026566798676182168 +0.565487 0.4315084643619017113 +0.628319 0.4535551023743125687 +0.691150 0.4684489721639752968 +0.753982 0.4759551277630860722 +0.816814 0.4759551277630863497 +0.879646 0.4684489721639752968 +0.942478 0.4535551023743126797 +1.005310 0.4315084643619018778 +1.068142 0.4026566798676179393 +1.130973 0.3674546220110178063 +1.193805 0.3264577612266306161 +1.256637 0.2803123461157628893 +1.319469 0.2297463108704850265 +1.382301 0.1755571071177310649 +1.445133 0.1185991090836186479 +1.507964 0.05977091079786012173 +1.570796 -1.401964684891840811e-16 +1.633628 -0.05977091079786012867 +1.696460 -0.1185991090836181483 +1.759292 -0.1755571071177305376 +1.822124 -0.229746310870484749 +1.884956 -0.2803123461157626672 +1.947787 -0.326457761226630172 +2.010619 -0.3674546220110179173 +2.073451 -0.4026566798676175507 +2.136283 -0.4315084643619017668 +2.199115 -0.4535551023743126242 +2.261947 -0.4684489721639758519 +2.324779 -0.4759551277630864607 +2.387610 -0.4759551277630861277 +2.450442 -0.4684489721639749082 +2.513274 -0.4535551023743128463 +2.576106 -0.4315084643619019333 +2.638938 -0.4026566798676183834 +2.701770 -0.3674546220110177508 +2.764602 -0.3264577612266314488 +2.827433 -0.2803123461157625007 +2.890265 -0.2297463108704854706 +2.953097 -0.1755571071177309539 +3.015929 -0.1185991090836186895 +3.078761 -0.05977091079786037847 +3.141593 1.425743635182407722e-16 +3.204425 0.05977091079786019112 +3.267256 0.1185991090836186063 +3.330088 0.1755571071177305098 +3.392920 0.2297463108704850543 +3.455752 0.2803123461157625562 +3.518584 0.3264577612266305606 +3.581416 0.3674546220110170847 +3.644247 0.4026566798676178283 +3.707079 0.4315084643619019888 +3.769911 0.4535551023743124577 +3.832743 0.4684489721639751858 +3.895575 0.4759551277630860167 +3.958407 0.4759551277630860722 +4.021239 0.4684489721639746862 +4.084070 0.4535551023743122911 +4.146902 0.4315084643619018778 +4.209734 0.4026566798676178838 +4.272566 0.3674546220110170291 +4.335398 0.3264577612266304496 +4.398230 0.2803123461157624452 +4.461062 0.2297463108704845824 +4.523893 0.1755571071177310372 +4.586725 0.1185991090836189532 +4.649557 0.05977091079785964989 +4.712389 2.43532537710590483e-16 +4.775221 -0.05977091079785951111 +4.838053 -0.1185991090836182593 +4.900885 -0.1755571071177308429 +4.963716 -0.229746310870485082 +5.026548 -0.2803123461157622787 +5.089380 -0.3264577612266307272 +5.152212 -0.3674546220110169736 +5.215044 -0.4026566798676180503 +5.277876 -0.4315084643619020444 +5.340708 -0.4535551023743122911 +5.403539 -0.4684489721639749082 +5.466371 -0.4759551277630857391 +5.529203 -0.4759551277630865163 +5.592035 -0.4684489721639766291 +5.654867 -0.4535551023743119581 +5.717699 -0.4315084643619015448 +5.780530 -0.4026566798676182724 +5.843362 -0.3674546220110180839 +5.906194 -0.3264577612266302276 +5.969026 -0.2803123461157628338 +6.031858 -0.2297463108704853041 +6.094690 -0.1755571071177306763 +6.157522 -0.1185991090836190087 +6.220353 -0.059770910797860205 +6.283185 -4.802722211955578559e-16 + + diff --git a/test/test_driscollhealy/mp_harmonic_im_r8.00.th.asc b/test/test_driscollhealy/mp_harmonic_im_r8.00.th.asc new file mode 100644 index 0000000..41a6a88 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_im_r8.00.th.asc @@ -0,0 +1,55 @@ +"Time = 0 +0.000000 0.6307831305050407567 +0.061600 0.007024422862686013555 +0.123200 -8.818650445389855207e-16 +0.184800 2.386973079634007145e-16 +0.246399 3.090954546896031658e-16 +0.307999 1.88689857105083995e-16 +0.369599 1.787601751083378799e-16 +0.431199 6.42809675982812758e-17 +0.492799 -5.223354783714025765e-17 +0.554399 1.052874518799022369e-16 +0.615999 1.117372074496677921e-16 +0.677598 8.829044623742264276e-17 +0.739198 -6.907573430162997911e-17 +0.800798 -3.90179774949959511e-17 +0.862398 6.159252871967265588e-17 +0.923998 5.627127698640937574e-17 +0.985598 4.652774069628360505e-17 +1.047198 7.285297644863500048e-17 +1.108797 -6.542430094850337162e-18 +1.170397 6.465233305208766899e-17 +1.231997 1.779484214503371193e-17 +1.293597 -7.602576948937617888e-18 +1.355197 1.390238123941617177e-17 +1.416797 2.070976566529800156e-17 +1.478397 1.584136196322225341e-17 +1.539996 2.282432166526202462e-17 +1.601596 -3.480861916600505987e-17 +1.663196 4.246513194032316906e-18 +1.724796 8.325741691803758443e-18 +1.786396 8.248138457401257405e-18 +1.847996 1.165433598435465015e-17 +1.909596 6.66306970460519502e-18 +1.971195 5.865616208071797884e-18 +2.032795 4.518590889726565158e-18 +2.094395 5.639724593917154933e-19 +2.155995 8.462818399384888541e-19 +2.217595 7.220434734079161545e-18 +2.279195 -1.940108903329609695e-18 +2.340795 1.091953982035055254e-18 +2.402394 9.65339816474644947e-19 +2.463994 5.865854807296045742e-19 +2.525594 1.106952583858313206e-18 +2.587194 1.616933492592317859e-18 +2.648794 -8.596749700189567001e-19 +2.710394 3.992357306555509448e-19 +2.771994 1.582759388194540812e-19 +2.833593 2.217055110251093269e-20 +2.895193 2.827594302876126491e-20 +2.956793 2.496238261006691607e-20 +3.018393 7.633200804732358595e-21 +3.079993 1.074058560061705883e-21 +3.141593 -7.22781544709519842e-37 + + diff --git a/test/test_driscollhealy/mp_harmonic_l0_m0_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l0_m0_r8.00.asc new file mode 100644 index 0000000..23eaf79 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l0_m0_r8.00.asc @@ -0,0 +1 @@ +0.000000 0 0 diff --git a/test/test_driscollhealy/mp_harmonic_l1_m-1_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l1_m-1_r8.00.asc new file mode 100644 index 0000000..23eaf79 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l1_m-1_r8.00.asc @@ -0,0 +1 @@ +0.000000 0 0 diff --git a/test/test_driscollhealy/mp_harmonic_l1_m0_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l1_m0_r8.00.asc new file mode 100644 index 0000000..23eaf79 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l1_m0_r8.00.asc @@ -0,0 +1 @@ +0.000000 0 0 diff --git a/test/test_driscollhealy/mp_harmonic_l1_m1_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l1_m1_r8.00.asc new file mode 100644 index 0000000..23eaf79 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l1_m1_r8.00.asc @@ -0,0 +1 @@ +0.000000 0 0 diff --git a/test/test_driscollhealy/mp_harmonic_l2_m-1_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l2_m-1_r8.00.asc new file mode 100644 index 0000000..eab32f1 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l2_m-1_r8.00.asc @@ -0,0 +1 @@ +0.000000 -8.332830142340763872e-17 1.022973937129801201e-16 diff --git a/test/test_driscollhealy/mp_harmonic_l2_m-2_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l2_m-2_r8.00.asc new file mode 100644 index 0000000..2fa5550 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l2_m-2_r8.00.asc @@ -0,0 +1 @@ +0.000000 2.138002145330596485e-10 2.446208802313476189e-17 diff --git a/test/test_driscollhealy/mp_harmonic_l2_m0_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l2_m0_r8.00.asc new file mode 100644 index 0000000..08d5282 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l2_m0_r8.00.asc @@ -0,0 +1 @@ +0.000000 -1.931250071812291289e-16 2.833087430018189041e-07 diff --git a/test/test_driscollhealy/mp_harmonic_l2_m1_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l2_m1_r8.00.asc new file mode 100644 index 0000000..5b6f793 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l2_m1_r8.00.asc @@ -0,0 +1 @@ +0.000000 -8.269870904790199682e-17 -5.485377474579588873e-17 diff --git a/test/test_driscollhealy/mp_harmonic_l2_m2_r8.00.asc b/test/test_driscollhealy/mp_harmonic_l2_m2_r8.00.asc new file mode 100644 index 0000000..cfe52d1 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_l2_m2_r8.00.asc @@ -0,0 +1 @@ +0.000000 3.141510003733607892 -2.084879755699846063e-17 diff --git a/test/test_driscollhealy/mp_harmonic_re_r8.00.ph.asc b/test/test_driscollhealy/mp_harmonic_re_r8.00.ph.asc new file mode 100644 index 0000000..219f02c --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_re_r8.00.ph.asc @@ -0,0 +1,104 @@ +"Time = 0 +0.000000 0.4768960190949341205 +0.062832 0.4731354400496181056 +0.125664 0.4619133783331734744 +0.188496 0.443406739337559197 +0.251327 0.417907108617889822 +0.314159 0.3858170683438172932 +0.376991 0.3476422227979080182 +0.439823 0.3039851008795564358 +0.502655 0.255533712494027232 +0.565487 0.2030525812388497831 +0.628319 0.1473690189180452692 +0.691150 0.08936138578181727909 +0.753982 0.02994455985679102608 +0.816814 -0.02994455985679060975 +0.879646 -0.08936138578181719583 +0.942478 -0.1473690189180451304 +1.005310 -0.2030525812388500884 +1.068142 -0.255533712494027232 +1.130973 -0.3039851008795565468 +1.193805 -0.3476422227979084623 +1.256637 -0.3858170683438171267 +1.319469 -0.4179071086178900996 +1.382301 -0.443406739337558975 +1.445133 -0.461913378333173974 +1.507964 -0.4731354400496178836 +1.570796 -0.4768960190949332878 +1.633628 -0.4731354400496186607 +1.696460 -0.461913378333174196 +1.759292 -0.4434067393375594746 +1.822124 -0.4179071086178899885 +1.884956 -0.3858170683438170712 +1.947787 -0.3476422227979082957 +2.010619 -0.3039851008795566023 +2.073451 -0.255533712494027232 +2.136283 -0.2030525812388493112 +2.199115 -0.1473690189180452137 +2.261947 -0.08936138578181727909 +2.324779 -0.02994455985679057852 +2.387610 0.02994455985679041893 +2.450442 0.08936138578181734848 +2.513274 0.1473690189180449917 +2.576106 0.2030525812388492279 +2.638938 0.2555337124940273985 +2.701770 0.3039851008795561027 +2.764602 0.3476422227979088508 +2.827433 0.3858170683438172932 +2.890265 0.4179071086178900996 +2.953097 0.443406739337559308 +3.015929 0.4619133783331736964 +3.078761 0.4731354400496182722 +3.141593 0.4768960190949315114 +3.204425 0.4731354400496181611 +3.267256 0.4619133783331736409 +3.330088 0.4434067393375588639 +3.392920 0.4179071086178901551 +3.455752 0.3858170683438176263 +3.518584 0.3476422227979086288 +3.581416 0.3039851008795569354 +3.644247 0.2555337124940272875 +3.707079 0.203052581238849561 +3.769911 0.1473690189180455468 +3.832743 0.08936138578181757053 +3.895575 0.02994455985679117527 +3.958407 -0.02994455985679039117 +4.021239 -0.08936138578181766767 +4.084070 -0.1473690189180449084 +4.146902 -0.2030525812388496165 +4.209734 -0.2555337124940270654 +4.272566 -0.3039851008795565468 +4.335398 -0.3476422227979082402 +4.398230 -0.3858170683438171822 +4.461062 -0.4179071086178904326 +4.523893 -0.443406739337559086 +4.586725 -0.4619133783331737519 +4.649557 -0.4731354400496182722 +4.712389 -0.4768960190949321221 +4.775221 -0.4731354400496184387 +4.838053 -0.4619133783331735299 +4.900885 -0.4434067393375591415 +4.963716 -0.4179071086178902106 +5.026548 -0.3858170683438174597 +5.089380 -0.3476422227979087953 +5.152212 -0.3039851008795573795 +5.215044 -0.2555337124940266214 +5.277876 -0.2030525812388495333 +5.340708 -0.1473690189180454357 +5.403539 -0.08936138578181807013 +5.466371 -0.02994455985679158813 +5.529203 0.02994455985679015872 +5.592035 0.08936138578181786196 +5.654867 0.1473690189180450194 +5.717699 0.2030525812388488394 +5.780530 0.2555337124940263438 +5.843362 0.3039851008795563803 +5.906194 0.3476422227979081847 +5.969026 0.3858170683438172932 +6.031858 0.417907108617889711 +6.094690 0.4434067393375583643 +6.157522 0.4619133783331743626 +6.220353 0.4731354400496186052 +6.283185 0.476896019094934509 + + diff --git a/test/test_driscollhealy/mp_harmonic_re_r8.00.th.asc b/test/test_driscollhealy/mp_harmonic_re_r8.00.th.asc new file mode 100644 index 0000000..36c31b5 --- /dev/null +++ b/test/test_driscollhealy/mp_harmonic_re_r8.00.th.asc @@ -0,0 +1,55 @@ +"Time = 0 +0.000000 -4.885011147332195149e-16 +0.061600 0.6225628802689283647 +0.123200 0.626011174765530698 +0.184800 0.6200885652671891402 +0.246399 0.6118754000725465936 +0.307999 0.601448949309237535 +0.369599 0.588906848052695775 +0.431199 0.5743658428948120731 +0.492799 0.5579603039292125866 +0.554399 0.5398405168177836666 +0.615999 0.5201707733415493751 +0.677598 0.4991272890280229468 +0.739198 0.4768960190949341205 +0.800798 0.4536703256043045362 +0.862398 0.4296485914087362912 +0.923998 0.4050318059467630238 +0.985598 0.3800211327412000295 +1.047198 0.3548155100560987929 +1.108797 0.3296093187267329361 +1.170397 0.3045901516673371057 +1.231997 0.2799367090384122236 +1.293597 0.2558168423155763316 +1.355197 0.2323857922780937779 +1.416797 0.2097846387827586034 +1.478397 0.1881389819792761153 +1.539996 0.1675578340283152279 +1.601596 0.1481328320161413303 +1.663196 0.1299376592113383189 +1.724796 0.1130277740794584929 +1.786396 0.09744041746750978228 +1.847996 0.0831948475902995338 +1.909596 0.07029291427378291102 +1.971195 0.05871980508769748808 +2.032795 0.0484450861834871474 +2.094395 0.03942394482988391402 +2.155995 0.03159862529246817292 +2.217595 0.02490004252143376146 +2.279195 0.01924954067796907226 +2.340795 0.01456077249599862025 +2.402394 0.01074165773964111219 +2.463994 0.007696404127104090713 +2.525594 0.005327525361079997121 +2.587194 0.003537890131800171169 +2.648794 0.002232672092008527748 +2.710394 0.00132127104105327531 +2.771994 0.0007190932719948640991 +2.833593 0.000349210492446978024 +2.895193 0.0001438537891083224067 +2.956793 4.571824601967693921e-05 +3.018393 9.059240056068421681e-06 +3.079993 5.663901966311766602e-07 +3.141593 -1.077203552463963929e-22 + + diff --git a/test/test_driscollhealy/test_driscollhealy.par b/test/test_driscollhealy/test_driscollhealy.par new file mode 100644 index 0000000..4edac91 --- /dev/null +++ b/test/test_driscollhealy/test_driscollhealy.par @@ -0,0 +1,90 @@ + +ActiveThorns = "CoordBase SymBase Boundary CartGrid3d IOUtil Carpet CarpetLib CarpetInterp AEILocalInterp InitBase Multipole LoopControl" + +############################################################# +# Grid +############################################################# + +CartGrid3D::type = "coordbase" +CartGrid3D::domain = "full" +CartGrid3D::avoid_origin = "no" + +CoordBase::domainsize = minmax +CoordBase::xmin = -10 +CoordBase::ymin = -10 +CoordBase::zmin = -10 +CoordBase::xmax = 10 +CoordBase::ymax = 10 +CoordBase::zmax = 10 +CoordBase::dx = 0.2 +CoordBase::dy = 0.2 +CoordBase::dz = 0.2 +CoordBase::boundary_size_x_lower = 2 +CoordBase::boundary_size_y_lower = 2 +CoordBase::boundary_size_z_lower = 2 +CoordBase::boundary_shiftout_x_lower = 0 +CoordBase::boundary_shiftout_y_lower = 0 +CoordBase::boundary_shiftout_z_lower = 0 +CoordBase::boundary_size_x_upper = 2 +CoordBase::boundary_size_y_upper = 2 +CoordBase::boundary_size_z_upper = 2 +CoordBase::boundary_shiftout_x_upper = 0 +CoordBase::boundary_shiftout_y_upper = 0 +CoordBase::boundary_shiftout_z_upper = 0 + +############################################################# +# Carpet +############################################################# + +Carpet::ghost_size = 3 +Carpet::domain_from_coordbase = "yes" +Carpet::poison_new_timelevels = "yes" +Carpet::check_for_poison = "no" +CarpetLib::poison_value = 113 +Carpet::init_fill_timelevels = "yes" + +############################################################# +# CarpetLib +############################################################# + +CarpetLib::poison_new_memory = "yes" + +############################################################# +# Cactus +############################################################# + +Cactus::terminate = "iteration" +Cactus::cctk_itlast = 0 + +############################################################# +# Multipole +############################################################# + +Multipole::nradii = 1 +Multipole::ntheta = 51 +Multipole::radius[0] = 8.0 +Multipole::variables = "Multipole::harmonic_re{sw=-2 cmplx='Multipole::harmonic_im' name='harmonic'}" +Multipole::integration_method = "DriscollHealy" + +Multipole::enable_test = "yes" +Multipole::test_l = 2 +Multipole::test_m = 2 +Multipole::out_1d_every = 1 + +############################################################# +# Output +############################################################# + +IO::out_dir = $parfile +IO::out_fileinfo = "none" + +# Enabling 1D output for the test grid functions would be helpful for +# localising any failures but it makes the tests dependent on the +# number of processors, as CarpetIOASCII's output is dependent on +# this. + +# CarpetIOASCII::out1d_vars = "Multipole::harmonics" +# CarpetIOASCII::out1d_every = 1 +# CarpetIOASCII::out1d_x = yes +# CarpetIOASCII::out1d_y = yes +# CarpetIOASCII::out1d_z = yes -- cgit v1.2.3