module Math.AngerJ
( AngerResult(..), angerJ )
where
import Data.Complex ( imagPart, realPart, Complex(..) )
import Numerical.Integration ( integration, IntegralResult(..) )
import Foreign.C ( CDouble )
data AngerResult = AngerResult {
AngerResult -> Complex Double
_result :: Complex Double
, AngerResult -> (Double, Double)
_errors :: (Double, Double)
, AngerResult -> (Int, Int)
_codes :: (Int, Int)
} deriving Int -> AngerResult -> ShowS
[AngerResult] -> ShowS
AngerResult -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [AngerResult] -> ShowS
$cshowList :: [AngerResult] -> ShowS
show :: AngerResult -> String
$cshow :: AngerResult -> String
showsPrec :: Int -> AngerResult -> ShowS
$cshowsPrec :: Int -> AngerResult -> ShowS
Show
cpxdbl2cpxcdbl :: Complex Double -> Complex CDouble
cpxdbl2cpxcdbl :: Complex Double -> Complex CDouble
cpxdbl2cpxcdbl Complex Double
z = forall a b. (Real a, Fractional b) => a -> b
realToFrac (forall a. Complex a -> a
realPart Complex Double
z) forall a. a -> a -> Complex a
:+ forall a b. (Real a, Fractional b) => a -> b
realToFrac (forall a. Complex a -> a
imagPart Complex Double
z)
angerJ :: Complex Double
-> Complex Double
-> Double
-> Int
-> IO AngerResult
angerJ :: Complex Double -> Complex Double -> Double -> Int -> IO AngerResult
angerJ Complex Double
nu Complex Double
z Double
err Int
subdiv = do
let z' :: Complex CDouble
z' = Complex Double -> Complex CDouble
cpxdbl2cpxcdbl Complex Double
z
x' :: CDouble
x' = forall a. Complex a -> a
realPart Complex CDouble
z'
y' :: CDouble
y' = forall a. Complex a -> a
imagPart Complex CDouble
z'
nu' :: Complex CDouble
nu' = Complex Double -> Complex CDouble
cpxdbl2cpxcdbl Complex Double
nu
a' :: CDouble
a' = forall a. Complex a -> a
realPart Complex CDouble
nu'
b' :: CDouble
b' = forall a. Complex a -> a
imagPart Complex CDouble
nu'
reintegrand :: CDouble -> CDouble
reintegrand CDouble
t = forall a. Floating a => a -> a
cosh(CDouble
b' forall a. Num a => a -> a -> a
* CDouble
t forall a. Num a => a -> a -> a
- CDouble
y' forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
sin CDouble
t) forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
cos(CDouble
a' forall a. Num a => a -> a -> a
* CDouble
t forall a. Num a => a -> a -> a
- CDouble
x' forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
sin CDouble
t)
imintegrand :: CDouble -> CDouble
imintegrand CDouble
t = -forall a. Floating a => a -> a
sinh(CDouble
b' forall a. Num a => a -> a -> a
* CDouble
t forall a. Num a => a -> a -> a
- CDouble
y' forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
sin CDouble
t) forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
sin(CDouble
a' forall a. Num a => a -> a -> a
* CDouble
t forall a. Num a => a -> a -> a
- CDouble
x' forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
sin CDouble
t)
IntegralResult
re <- (CDouble -> CDouble)
-> Double -> Double -> Double -> Double -> Int -> IO IntegralResult
integration CDouble -> CDouble
reintegrand Double
0 forall a. Floating a => a
pi Double
0.0 Double
err Int
subdiv
IntegralResult
im <- (CDouble -> CDouble)
-> Double -> Double -> Double -> Double -> Int -> IO IntegralResult
integration CDouble -> CDouble
imintegrand Double
0 forall a. Floating a => a
pi Double
0.0 Double
err Int
subdiv
forall (m :: * -> *) a. Monad m => a -> m a
return AngerResult {
_result :: Complex Double
_result = (IntegralResult -> Double
_value IntegralResult
re forall a. a -> a -> Complex a
:+ IntegralResult -> Double
_value IntegralResult
im) forall a. Fractional a => a -> a -> a
/ forall a. Floating a => a
pi
, _errors :: (Double, Double)
_errors = (IntegralResult -> Double
_error IntegralResult
re, IntegralResult -> Double
_error IntegralResult
im)
, _codes :: (Int, Int)
_codes = (IntegralResult -> Int
_code IntegralResult
re, IntegralResult -> Int
_code IntegralResult
im)
}