Skip to content

Commit f1ee586

Browse files
authored
Arnoldifit (#632)
* update reference; stub for fourier extension
1 parent 1f96881 commit f1ee586

1 file changed

Lines changed: 30 additions & 3 deletions

File tree

src/polynomials/standard-basis/standard-basis.jl

Lines changed: 30 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -665,9 +665,7 @@ polynomial evaluations afterwards for higher-degree polynomials.
665665
666666
The two main functions are translations from example code in:
667667
668-
*VANDERMONDE WITH ARNOLDI*;
669-
PABLO D. BRUBECK, YUJI NAKATSUKASA, AND LLOYD N. TREFETHEN;
670-
[arXiv:1911.09988](https://people.maths.ox.ac.uk/trefethen/vander_revised.pdf)
668+
Pablo D. Brubeck, Yuji Nakatsukasa, and Lloyd N. Trefethen, “Vandermonde with Arnoldi”, SIAM Review 63, doi.org/10.1137/19M130100X (2021).
671669
672670
For more details, see also:
673671
@@ -744,6 +742,35 @@ function polyvalA(d, H::AbstractMatrix{S}, s::T) where {T, S}
744742
sum(W[i]*d[i] for i in eachindex(d))
745743
end
746744

745+
#=
746+
# Finds a + bim following formula 5.3 of paper, which solves
747+
# Re(Vₙ) [a₀ ⋯ aₙ] - Im(Vₙ[2:n]) * [b₁ ⋯ bₙ] = [y₁ ⋯ yₘ]
748+
function polyfit_fourier_extension(x, y, n::Int=length(x) - 1; var=Var(:x))
749+
m = length(x)
750+
T = eltype(y)
751+
Q = ones(T, m, n+1)
752+
H = zeros(T, n+1, n)
753+
754+
q = zeros(T, m)
755+
756+
@inbounds for k = 1:n
757+
q .= x .* Q[:,k]
758+
for j in 1:k
759+
λ = dot(Q[:,j], q)/m
760+
H[j,k] = λ
761+
q .-= λ * Q[:,j]
762+
end
763+
H[k+1,k] = norm(q)/sqrt(m)
764+
Q[:,k+1] .= q/H[k+1,k]
765+
end
766+
767+
d = [real(Q) imag(view(Q, :, 2:n+1))] \ y
768+
d = d[1:(n+1)] - vcat(0, view(d, (n+2):(2n+1)))
769+
ArnoldiFit{eltype(d),Symbol(var)}(d, H)
770+
end
771+
=#
772+
773+
747774
# Polynomial Interface
748775
"""
749776
ArnoldiFit

0 commit comments

Comments
 (0)