diff options
-rw-r--r-- | interface.ccl | 8 | ||||
-rw-r--r-- | schedule.ccl | 6 | ||||
-rw-r--r-- | src/tests.cc | 46 |
3 files changed, 60 insertions, 0 deletions
diff --git a/interface.ccl b/interface.ccl index 2e893ab..26d7ef7 100644 --- a/interface.ccl +++ b/interface.ccl @@ -20,3 +20,11 @@ CCTK_REAL test_integration_convergence_orders type=SCALAR test_trapezoidal_convergence_order, test_simpson_convergence_order } "Test integration convergence orders" + +CCTK_REAL test_integration_symmetries type=SCALAR +{ + test_midpoint_pi_symmetry, + test_trapezoidal_pi_symmetry, + test_simpson_pi_symmetry, + test_driscollhealy_pi_symmetry +} "Test integration symmetries" diff --git a/schedule.ccl b/schedule.ccl index 328e1de..3b46cae 100644 --- a/schedule.ccl +++ b/schedule.ccl @@ -4,6 +4,7 @@ if (enable_test) { STORAGE: harmonics[1] STORAGE: test_integration_convergence_orders + STORAGE: test_integration_symmetries } schedule Multipole_Calc at CCTK_ANALYSIS after (calc_np,PsiKadelia,Accelerator_CopyBack) @@ -34,4 +35,9 @@ if (enable_test) { LANG: C } "Test convergence of integration" + + schedule Multipole_TestIntegrationSymmetry at CCTK_PARAMCHECK + { + LANG: C + } "Test symmetry of integration" } diff --git a/src/tests.cc b/src/tests.cc index 42e79b3..de711eb 100644 --- a/src/tests.cc +++ b/src/tests.cc @@ -40,6 +40,38 @@ static CCTK_REAL test_integral(int n, CCTK_REAL (*integration_fn) (const CCTK_RE return result; } + +static CCTK_REAL test_pi_symmetric_sphere_integral(CCTK_REAL (*integration_fn) (const CCTK_REAL *, int, int, CCTK_REAL, CCTK_REAL)) +{ + const int n = 100; + const int nth = n; + const int nph = n; + const int array_size=(nth+1)*(nph+1); + const CCTK_REAL PI = acos(-1.0); + + CCTK_REAL *f = new CCTK_REAL[array_size]; + + const CCTK_REAL dth = PI/nth; + const CCTK_REAL dph = 2*PI/nph; + + for (int ith = 0; ith <= nth; ith++) + { + for (int iph = 0; iph <= nph; iph++) + { + const int i = Multipole_Index(ith, iph, nth); + + const CCTK_REAL th = ith*dth; + const CCTK_REAL ph = iph*dph; + + f[i] = -(cos(ph)*sqrt(5/PI)*pow(cos(th/2.),3)*sin(th/2.)); + } + } + + const CCTK_REAL result = integration_fn(f, nth, nph, dth, dph); + delete [] f; + return result; +} + CCTK_REAL integration_convergence_order(CCTK_REAL (*integration_fn) (const CCTK_REAL *, int, int, CCTK_REAL, CCTK_REAL)) { const int n1 = 100; @@ -61,6 +93,20 @@ void Multipole_TestIntegrationConvergence(CCTK_ARGUMENTS) *test_midpoint_convergence_order = integration_convergence_order(&Midpoint2DIntegral); } +void Multipole_TestIntegrationSymmetry(CCTK_ARGUMENTS) +{ + DECLARE_CCTK_ARGUMENTS; + + *test_simpson_pi_symmetry = test_pi_symmetric_sphere_integral(&Simpson2DIntegral); + *test_midpoint_pi_symmetry = test_pi_symmetric_sphere_integral(&Midpoint2DIntegral); + *test_trapezoidal_pi_symmetry = test_pi_symmetric_sphere_integral(&Trapezoidal2DIntegral); + *test_driscollhealy_pi_symmetry = test_pi_symmetric_sphere_integral(&DriscollHealy2DIntegral); + printf("Pi symmetry Simpson integral: %.19g\n", *test_simpson_pi_symmetry); + printf("Pi symmetry midpoint integral: %.19g\n", *test_midpoint_pi_symmetry); + printf("Pi symmetry trapezoidal integral: %.19g\n", *test_trapezoidal_pi_symmetry); + printf("Pi symmetry Driscoll and Healy integral: %.19g\n", *test_driscollhealy_pi_symmetry); +} + // void Multipole_TestIntegrate(CCTK_ARGUMENTS) // { // const int n = 100; |