/
FresnelIntegrals.jl
95 lines (87 loc) · 2.54 KB
/
FresnelIntegrals.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
module FresnelIntegrals
using SpecialFunctions
using IrrationalConstants: sqrtπ
export fresnelc, fresnels, fresnel
"""
fresnelc(z::Number)
Calculate the normalized Fresnel cosine integral
```math
C(z) = \\int_{0}^{z} \\cos{\\left(\\frac{\\pi t^2}{2}\\right)} \\, \\mathrm{d}t
```
for the number ``z``.
"""
function fresnelc(z::Number)
x = (z * sqrtπ) / 2
re_x, im_x = reim(x)
a = (re_x + im_x) + (im_x - re_x) * im
b = (re_x - im_x) + (im_x + re_x) * im
re_erf_a, im_erf_a = reim(erf(a))
re_erf_b, im_erf_b = reim(erf(b))
re_y = (re_erf_a - im_erf_a + re_erf_b + im_erf_b) / 4
im_y = (im_erf_a + re_erf_a - re_erf_b + im_erf_b) / 4
y = re_y + im_y * im
return y
end
function fresnelc(z::Real)
x = (z * sqrtπ) / 2
a = x + x * im
re_erf_a, im_erf_a = reim(erf(a))
y = (re_erf_a + im_erf_a) / 2
return y
end
"""
fresnels(z::Number)
Calculate the normalized Fresnel sine integral
```math
S(z) = \\int_{0}^{z} \\sin{\\left(\\frac{\\pi t^2}{2}\\right)} \\, \\mathrm{d}t
```
for the number ``z``.
"""
function fresnels(z::Number)
x = (z * sqrtπ) / 2
re_x, im_x = reim(x)
a = (re_x + im_x) + (im_x - re_x) * im
b = (re_x - im_x) + (im_x + re_x) * im
re_erf_a, im_erf_a = reim(erf(a))
re_erf_b, im_erf_b = reim(erf(b))
re_y = (re_erf_a + im_erf_a + re_erf_b - im_erf_b) / 4
im_y = (im_erf_a - re_erf_a + re_erf_b + im_erf_b) / 4
y = re_y + im_y * im
return y
end
function fresnels(z::Real)
x = (z * sqrtπ) / 2
a = x + x * im
re_erf_a, im_erf_a = reim(erf(a))
y = (re_erf_a - im_erf_a) / 2
return y
end
"""
fresnel(z::Number)
Calculate the normalized cosine and sine fresnel integrals.
See also [`fresnels`](@ref), [`fresnelc`](@ref).
"""
function fresnel(z::Number)
x = (z * sqrtπ) / 2
re_x, im_x = reim(x)
a = (re_x + im_x) + (im_x - re_x) * im
b = (re_x - im_x) + (im_x + re_x) * im
re_erf_a, im_erf_a = reim(erf(a))
re_erf_b, im_erf_b = reim(erf(b))
re_y_sin = (re_erf_a + im_erf_a + re_erf_b - im_erf_b) / 4
im_y_sin = (im_erf_a - re_erf_a + re_erf_b + im_erf_b) / 4
re_y_cos = (re_erf_a - im_erf_a + re_erf_b + im_erf_b) / 4
im_y_cos = (im_erf_a + re_erf_a - re_erf_b + im_erf_b) / 4
y_sin = re_y_sin + im_y_sin * im
y_cos = re_y_cos + im_y_cos * im
return (y_cos, y_sin)
end
function fresnel(z::Real)
x = (z * sqrtπ) / 2
a = x + x * im
re_erf_a, im_erf_a = reim(erf(a))
y_sin = (re_erf_a - im_erf_a) / 2
y_cos = (re_erf_a + im_erf_a) / 2
return (y_cos, y_sin)
end
end # module