text: """
  Test all geodesic functions here.
  
  This file compares results of geodesic functions with the
  Black Hole Perturbation Toolkit written in Mathematica.
  """
  import pytest
  import numpy as np
  from mpmath import mpf, mp, almosteq
  from functions import radial_roots
  
  digits = 100  # accuracy requested
  eps = 10 ** (-digits)  # number of digits which may be different
  mp.dps = digits  # set precision
  
  
  def test_radial_sc_polar_circular():
      aa = 0
      slr = 7
      ecc = 0
      x = 0
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r2 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r3 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_circular():
      aa = 0
      slr = 6
      ecc = 0
      x = 0.5
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r2 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r3 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_circular_equatorial():
      aa = 0
      slr = 6
      ecc = 0
      x = 1
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r2 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r3 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_kerr_polar_circular():
      aa = 0.9
      slr = 6
      ecc = 0
      x = 0
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r2 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r3 = mpf(
          "4.561339049348941915090005888955490497757108282190278250268"
          + "4625285077077054441036020702100329684574966701226"
      )
      r4 = mpf(
          "0.483045392174790095490549865184069452309707597442151923012"
          + "88456761367172626346831436319680699788709729280148"
      )
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_kerr_circular():
      aa = 0.9
      slr = 6
      ecc = 0
      x = 0.5
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r2 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r3 = mpf(
          "2.107201636766182477809036698351808785794321226845138552595"
          + "5383217548393926383743377377252142128804042090704"
      )
      r4 = mpf(
          "0.547821639683837105742296185097970225527479832053292059782"
          + "11748376227722039540247078504892434983412451743782"
      )
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_kerr_circular_equatorial():
      aa = 0.9
      slr = 6
      ecc = 0
      x = 1
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r2 = mpf(
          "6.000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000"
      )
      r3 = mpf(
          "1.439963351423020220550535199930285603136722337191912206300"
          + "85034338007837900996302977689884413080123991761707306462946"
          + "81767914649"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_polar_low_ecc():
      aa = 0
      slr = 6
      ecc = 0.1
      x = 0
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.66666666666666666666666666666666666666666666666666666666666"
          + "6666666666666666666666666666666666666666666666666666666666666"
          + "666666667"
      )
      r2 = mpf(
          "5.45454545454545454545454545454545454545454545454545454545454"
          + "5454545454545454545454545454545454545454545454545454545454545"
          + "454545455"
      )
      r3 = mpf(
          "6.00000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "0000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_inclined_low_ecc():
      aa = 0
      slr = 6
      ecc = 0.1
      x = 0.5
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.66666666666666666666666666666666666666666666666666666666666"
          + "6666666666666666666666666666666666666666666666666666666666666"
          + "666666667"
      )
      r2 = mpf(
          "5.45454545454545454545454545454545454545454545454545454545454"
          + "5454545454545454545454545454545454545454545454545454545454545"
          + "454545455"
      )
      r3 = mpf(
          "6.00000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "0000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_equatorial_low_ecc():
      aa = 0
      slr = 6
      ecc = 0.1
      x = 1
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6.66666666666666666666666666666666666666666666666666666666666"
          + "6666666666666666666666666666666666666666666666666666666666666"
          + "666666667"
      )
      r2 = mpf(
          "5.45454545454545454545454545454545454545454545454545454545454"
          + "5454545454545454545454545454545454545454545454545454545454545"
          + "454545455"
      )
      r3 = mpf(
          "6.00000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "0000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_equatorial_medium_ecc():
      aa = 0
      slr = 6
      ecc = 0.5
      x = 1
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "12.0000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "000000000"
      )
      r2 = mpf(
          "4.00000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          "000000000"
      )
      r3 = mpf(
          "6.00000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "0000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_polar_high_ecc():
      aa = 0
      slr = 6
      ecc = 0.9999
      x = 0
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "60000.0000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "00000"
      )
      r2 = mpf(
          "3.00015000750037501875093754687734386719335966798339916995849"
          + "7924896244812240612030601530076503825191259562978148907445372"
          + "268613431"
      )
      r3 = mpf(
          "6.00000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "0000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_high_ecc():
      aa = 0
      slr = 6
      ecc = 0.9999
      x = 0.5
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "60000.0000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "00000"
      )
      r2 = mpf(
          "3.00015000750037501875093754687734386719335966798339916995849"
          + "7924896244812240612030601530076503825191259562978148907445372"
          + "268613431"
      )
      r3 = mpf(
          "6.00000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "0000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_sc_equatorial_high_ecc():
      aa = 0
      slr = 6
      ecc = 0.9999
      x = 1
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "60000.0000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "00000"
      )
      r2 = mpf(
          "3.00015000750037501875093754687734386719335966798339916995849"
          + "7924896244812240612030601530076503825191259562978148907445372"
          + "268613431"
      )
      r3 = mpf(
          "6.00000000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "0000000"
      )
      r4 = mpf("0")
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_kerr_polar_high_ecc():
      aa = 0.9
      slr = 6
      ecc = 0.9999
      x = 0
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "60000.0000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "00000"
      )
      r2 = mpf(
          "3.00015000750037501875093754687734386719335966798339916995849"
          + "7924896244812240612030601530076503825191259562978148907445372"
          + "268613431"
      )
      r3 = mpf(
          "4.23160150405948518323888791308221551987498093592515288030115"
          + "945354089115861330121257783337210134827499637485386611824933"
      )
      r4 = mpf(
          "0.49286586348019811079088367661175558062941254382055084863641"
          + "0998598718916996297039941129876200032081561790988663968206728"
      )
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_kerr_high_ecc():
      aa = 0.9
      slr = 6
      ecc = 0.9999
      x = 0.5
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "60000.0000000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "00000"
      )
      r2 = mpf(
          "3.00015000750037501875093754687734386719335966798339916995849"
          + "7924896244812240612030601530076503825191259562978148907445372"
          + "268613431"
      )
      r3 = mpf(
          "2.09909631071719289113848288619813896273182343075993333860146"
          + "83275256205764755741564757727720651989485563951626108414"
      )
      r4 = mpf(
          "0.54443257307207486035051683847000616594344254063022881916146"
          + "92087614587124232498346668603991435731832877296511569358"
      )
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  
  
  def test_radial_kerr_even_higher_ecc():
      aa = 0.9
      slr = 6
      ecc = 0.999999
      x = 0.5
      r1_ch, r2_ch, r3_ch, r4_ch = radial_roots(aa, slr, ecc, x, digits)
      r1 = mpf(
          "6000000.00000000000000000000000000000000000000000000000000000"
          + "0000000000000000000000000000000000000000000000000000000000000"
          + "00000"
      )
      r2 = mpf(
          "3.00000150000075000037500018750009375004687502343751171875585"
          + "9377929688964844482422241211120605560302780151390075695037847"
          + "518923759"
      )
      r3 = mpf(
          "2.09909500849781985664788052724951917390754315839076840288258"
          + "30079641104679477801203514429959540281638544173538"
      )
      r4 = mpf(
          "0.54443186482151970276714708796651839555550722672128921192288"
          + "90618833999920935428791168827444874586729766166531"
      )
      assert almosteq(r1_ch, r1, eps)
      assert almosteq(r2_ch, r2, eps)
      assert almosteq(r3_ch, r3, eps)
      assert almosteq(r4_ch, r4, eps)
  