taylor_expand function for Taylor1 and TaylorN #121

blas-ko · 2017-08-31T05:56:19Z

This makes reference to #68. A taylor_expand function is added for Taylor1 and TaylorN expansions around a Number or Vector, respectively.

Obs: The way "common" functions are expanded for arrays such as taylor_expand(cos.[0.0,1.0]) for julia 0.6 output a warning about a method deprecation for arrays. They suggest to use now the dot operator instead, but this generates some problems with user defined functions.

…tions.

`taylor_expand(f,x0::Int64)` and `taylor_expand(f,order)` clashed.

coveralls · 2017-08-31T06:07:59Z

Coverage increased (+0.5%) to 94.692% when pulling 65b5596 on blas-ko:taylor_expand into 2381ad4 on JuliaDiff:master.

lbenet · 2017-08-31T18:21:39Z

I like the idea of this PR, though I am still not sure about the name of the function (and perhaps other things), see below.

My proposal here is to twist a little the essence, and address #80. That is, given x::Taylor1 (for concreteness I'll focus on Taylor1s), the idea is to return x (probably update it in place), such that the new x is the polynomial evaluated at a different location. Then, update{T}(x::Taylor1{T}, t0::T) would be equivalent to evaluate(x, Taylor1( [t0,one(t0)], x.order )). The tricky part (and I am not sure if this is possible or not to accomplish) is to only update x; maybe the following works:

function update!(x::Taylor1, t0::TaylorSeries.NumberNotSeries)
    x.coeffs .= evaluate(x, Taylor1( [t0,one(t0)], x.order ))[:]
    return x
end

If this works, I guess this could also be extended to HomogeneousPolynomial and TaylorN, but I am not sure.

Surely, the original idea of this PR is more general, in the sense that it returns the evaluated expansion of a given function f around a point. This could still be considered.

The distinction of returning a Taylor1 or a TaylorN is by dispatch on x0. That's ok. Yet, I don't think that taylor_evaluate should change the number of variables nor the order for TaylorN; I think this should be done by the user explicitly, and taylor_evaluate should rely on get_order() and get_numvars().

lbenet · 2017-08-31T23:03:22Z

I think the idea above works

julia> using TaylorSeries

julia> function update!(x::Taylor1, t0::TaylorSeries.NumberNotSeries)
           x.coeffs .= evaluate(x, Taylor1( [t0,one(t0)], x.order ))[:]
           return nothing
       end
update! (generic function with 1 method)

julia> x = Taylor1([1.0,2.0,3.0])
 1.0 + 2.0 t + 3.0 t² + 𝒪(t³)

julia> update!(x, 1.0)

julia> x
 6.0 + 8.0 t + 3.0 t² + 𝒪(t³)

julia> update!(x, -1.0)

julia> x
 1.0 + 2.0 t + 3.0 t² + 𝒪(t³)

lbenet · 2017-08-31T23:04:11Z

Incidentally, the fact that tests are not passing may be related to changing the order and the number of variables; see above.

blas-ko · 2017-08-31T23:18:39Z

given x::Taylor1 (for concreteness I'll focus on Taylor1s), the idea is to return x (probably update it in place), such that the new x is the polynomial evaluated at a different location.

I think this is a great feature (that I think should be called something like shift! or displace!), but it's, in essence, different from what taylor_expand does. I can incorporate the update! feature with the name that we decide in this PR, but it will be defined via evaluate as you suggest and not via taylor_expand. What do you think @lbenet ?

Methods for evaluating HomogeneousPolynomial and TaylorN in HomogenousPolynomial and TaylorN are not yet defined. This should be constructed in order to implement the update! function accordingly.

Yet, I don't think that taylor_evaluate should change the number of variables nor the order for TaylorN; I think this should be done by the user explicitly, and taylor_evaluate should rely on get_order() and get_numvars().

The reason for this was that it took the dimension of the function f to expand it accordingly, but maybe it's a good practice for the user for him to change the set_variables environment a priori of the expansion.

Incidentally, the fact that tests are not passing may be related to changing the order and the number of variables; see above.

Output for travis (both Linux nightly & Mac nightly) is:

while loading /home/travis/.julia/v0.7/TaylorSeries/test/manyvariables.jl, in expression starting on line 7
Tests for HomogeneousPolynomial and TaylorN: Error During Test
  Test threw an exception of type MethodError
  Expression: taylor_expand(exp, [0, 0.0]) == exp.(set_variables("x", numvars=2))
  MethodError: no method matching exp(::Array{TaylorSeries.TaylorN{Float64},1})
  Closest candidates are:
    exp(!Matched::Float16) at math.jl:993
    exp(!Matched::Complex{Float16}) at math.jl:994
    exp(!Matched::BigFloat) at mpfr.jl:527
    ...

I think it has to do to broadcasting as exp.(x) still has some compatibility issues, but I'm not sure.

lbenet · 2017-08-31T23:46:45Z

I agree that there is a difference between update! (or displace!, but not shift!, since it has a different meaning in Base) and taylor_evaluate; the former overwrittes the given Taylor1, while the latter creates a new one, but centered at some given value. In any case, implementing both is maybe the best.

For the tests, I think it's better to check the actual series constructed directly, either using taylor_expand(f, x0, order) or as f(Taylor1([zero(x0),one(x0)],order)); they should be identic.

Regarding the TaylorN case, we should perhaps leave it aside, to begin, and then latter address it.

lbenet · 2017-09-01T09:20:36Z

I was checking this PR, and there is something odd in the many-variable case. For simplicity I will consider two variables. I hope the following comments are helpful.

My expectation is that taylor_expand should correspond the proper translation around the new point. That is, the expansion of f(x,y) around [x0,y0] should be f(x0+xT, y0+yT) where xT and yT are TaylorNs.

Yet, this is not the case currently.

julia> using TaylorSeries

julia> xT, yT = set_variables("x y");

julia> f(x,y) = (x+y)^2
f (generic function with 1 method)

julia> f(xT,yT)
 1.0 x² + 2.0 x y + 1.0 y² + 𝒪(‖x‖⁷)

julia> taylor_expand(f, [0.0,0.0])
ERROR: MethodError: no method matching f(::Array{TaylorSeries.TaylorN{Float64},1})
Closest candidates are:
  f(::Any, ::Any) at REPL[3]:1
Stacktrace:
 [1] #taylor_expand#33(::Int64, ::Function, ::#f, ::Array{Float64,1}) at /Users/benet/.julia/v0.6/TaylorSeries/src/other_functions.jl:146
 [2] taylor_expand(::Function, ::Array{Float64,1}) at /Users/benet/.julia/v0.6/TaylorSeries/src/other_functions.jl:143

The error is because in taylor_expand, X is actually a vector of the TaylorN independent variables, instead of being some sort of tuple. I guess that splatting (...) the argument of the function could be helpful here.

Currently, the function f is broadcasted, which means that f(X) returns a vector of the form [f(xT), f(yT)], which is not useful for the example above. This explains that

julia> taylor_expand(exp,[0,0.])
...
2-element Array{TaylorSeries.TaylorN{Float64},1}:
  1.0 + 1.0 x₁ + 0.5 x₁² + 0.16666666666666666 x₁³ + 0.041666666666666664 x₁⁴ + 0.008333333333333333 x₁⁵ + 0.0013888888888888887 x₁⁶ + 𝒪(‖x‖⁷)
  1.0 + 1.0 x₂ + 0.5 x₂² + 0.16666666666666666 x₂³ + 0.041666666666666664 x₂⁴ + 0.008333333333333333 x₂⁵ + 0.0013888888888888887 x₂⁶ + 𝒪(‖x‖⁷)

is a vector, and each entry depends only on one of the TaylorN variables. (I supressed some irrelevant warnings.)

I think we should consider the function to be defined as in the first example above. Broadcasting (f.(x,y)) may be considered to extend this to the case where f(x,y) returns a vector.

…es.jl into taylor_expand

…bles is changed.

blas-ko · 2017-09-05T18:10:38Z

I've covered the case you're telling about @lbenet ! It works as

julia> using TaylorSeries

julia> xT, yT = set_variables("x y");

julia> f(x,y) = (x+y)^2
f (generic function with 1 method)

julia> f(xT,yT)
 1.0 x² + 2.0 x y + 1.0 y² + 𝒪(‖x‖⁷)

julia> taylor_expand(f, 0.0, 0.0)
1.0 x² + 2.0 x y + 1.0 y² + 𝒪(‖x‖⁷)

julia> @which taylor_expand(f, 0.0, 0.0)
taylor_expand(f::Function, x0...) at .../src/other_functions.jl:152

Additionally, taylor_expand for TaylorN throws a warning if the set variables are note equal to the required for the expansion.

EDIT: I think I'll address #80 in a new PR if you agree with me in this.

EDIT 2: Rethinking this issue, I feel that taylor_expand! would make sense for shifting an existing Taylor1 or TaylorN around a new point x0. In this fashion, taylor_expand would take a Function and x0 as arguments to make the expansion, while taylor_expand! would take an existing Taylor{1,N} and update it around its second argument x0. What do you think, @lbenet?

By the way, taylor_expand! for TaylorN would need an evaluate method for TaylorN in a Vector{TaylorN}, as @PerezHz points out in #119.

coveralls · 2017-09-05T20:40:52Z

Coverage decreased (-0.4%) to 93.814% when pulling 3c914ae on blas-ko:taylor_expand into 2381ad4 on JuliaDiff:master.

blas-ko · 2017-09-08T00:08:48Z

I edited my last comment, @lbenet. I'm not sure if Github notifies this...
Cheers!

@lbenet

Followed suggestion from @lbenet.

@blas-ko

JuliaDiff#118) * Add function-like behavior for Taylor1 * Relocate new code * Add function-like behavior for TaylorN * Fix TaylorN functor methods * Add tests for Taylor1 * Add more Taylor1 tests * Add TaylorN tests; more Taylor1 tests; add missing evaluate methods * Add function-like behavior for HomogeneousPolynomial and corresponding tests * Add missing tests for HomogeneousPolynomial * Add another test for Taylor1s * Fix test * Add extra evaluate method for Taylor1 (suggested by @blas-ko) * Add an evaluate test for mixtures (more to come) * Add tests for mixtures; add/fix evaluate methods * A small fix * Add missing evaluate methods for mixtures and tests * Update docs * Add evaluate method and tests * Fix new method * Update docstrings * Changes suggested by @lbenet 's review

coveralls · 2017-10-02T02:27:09Z

Coverage decreased (-1.0%) to 94.943% when pulling 5d978eb on blas-ko:taylor_expand into 3b68fb3 on JuliaDiff:master.

lbenet · 2017-10-06T16:50:03Z

@blas-ko Can you rebase this branch to current master? This will also restart travis check and then I'll review it.

…tions.

`taylor_expand(f,x0::Int64)` and `taylor_expand(f,order)` clashed.

…es.jl into taylor_expand

blas-ko · 2017-10-06T17:32:39Z

Rebased! (Hope I did it correctly since I don't have experience rebasing yet).

They use set_variables internally.

lbenet

In general, LGTM. Yet, i have some objection with respect to the use of set_variables, for which I suggest an alternative. Also I suggest that taylor_expand! is renamed as update! (see #80) and is also exported.

Cc: @dpsanders

lbenet · 2017-10-06T19:31:59Z

src/evaluate.jl

+        x0::Array{TaylorN{T},1})
+    @assert length(x) == length(x0)
+    @inbounds for i in eachindex(x)
+        x0[i] = evaluate( x[i], δx )


Perhaps using here (as well in similar occurrences in other methods) evaluate!(x[i], δx, x0[i]) instead, yields better performance.

Agreed for evaluate methods but not for evaluate! methods as evaluate!, as the docs say, just work for arrays.

Good point!

lbenet · 2017-10-06T19:32:08Z

src/evaluate.jl

+    @assert length(vals) == get_numvars()
+
+    num_vars = get_numvars()
+    ct = coeff_table


This creates a copy of coeff_table. Perhaps using a view of it, or directly coeff_table saves this.

I didn't notice this before, sorry, but I think this should also be changed in other occurrences.

lbenet · 2017-10-06T19:32:13Z

src/other_functions.jl

+
+#taylor_expand function for Taylor1
+doc"""
+    taylor_expand(f,x0;order)


Can you add white spaces? I mean taylor_expand(f, x0; order).

As commented below, I am not in favor of changing the number of variables inside the function.

lbenet · 2017-10-06T19:32:16Z

src/other_functions.jl

+computed, changing the number of `TaylorN` variables
+according to the dimension of `f`.
+"""
+function taylor_expand(f::Function; order::Int64=15)


Note that this will return a Taylor1{Float64} expansion. My point focuses in the Float64 part.

Do you think taylor_expand(T::Type,f::Function;order::Int64) is a good idea? In this case, should the same be done when passing an explicit x0 to taylor_expand?

lbenet · 2017-10-06T19:32:18Z

src/other_functions.jl

+end
+
+#taylor_expand function for TaylorN
+function taylor_expand{T<:Number}(f::Function, x0::Vector{T}; order::Int64=get_order()) #a Taylor expansion around x0


I agree that, for Taylor1 expansions it makes sense to allow the user to set the order of the expansion; I am not in favor of this for TaylorN expansions, at least as it is now. My objection is related to the use of set_variables, since it changes some internals, in particular, the maximum order of the TaylorN expansions. This costs time; also, the idea of the current implementation is to set the parameters related to TaylorN internals only once at the very beginning.

I think this method should assert that the the length of x0 corresponds to get_numvars(), as you do it, and only when this is fulfilled evaluate f(X .+ x0) (notice the spaces). In this case, using set_variables will not change anything internally.

An alternative that avoids using set_variables is for example:

function taylor_expand{T<:Number}(f::Function, x0::Vector{T}; order::Int64=get_order()) ll = length(x0) # check that x0 has the correct length and `order` does not exceed the maximum order of the expansions @assert ll == get_numvars() && order <= get_order() X = Array{TaylorN{T}}(length(x0) for i in eachindex(X) X[i] = TaylorN(T, i, order=order) end return f( X .+ x0 ) end

This alternative permits to have expansions of given order, smaller or equal to get_order(), e.g., order=2 when get_order() is 6.

Really in favor of this. By the way, really good "type preserving strategy" doing the X[i] = TaylorN(T, i, order).

Regarding the problem with travis in Julia 0.5 (and considering your comment about this), I think a possible solution is to have X as the shifted vector. I mean, define X[i] = x0[i] + TaylorN(T, i, order=order) inside the for-loop, and then only return the evaluation of f( X ). This could be used in each method of taylor_expand.

I think this could actually solve the problem with travis, and yield a tiny performance improvement, since only one loop would be needed.

lbenet · 2017-10-06T19:32:23Z

src/other_functions.jl

+function taylor_expand(f::Function, x0...; order::Int64=get_order()) #a Taylor expansion around x0
+    ll = length(x0)
+    ll == get_numvars() ? X = get_variables() : begin
+        X = set_variables("x",order=order,numvars=ll)


Same comment as above.

lbenet · 2017-10-06T19:32:25Z

src/other_functions.jl

+
+Takes `a <: Union{Taylo1,TaylorN}` and expands it around the coordinate `x0`.
+"""
+function taylor_expand!(a::Taylor1)


What about naming this function update! (as proposed in #80), and export it? Since this method is some sort of fallback, i think it can be written so a::Union{Taylor1, TaylorN}, which also covers another method below.

lbenet · 2017-10-06T19:33:21Z

src/evaluate.jl

@@ -224,7 +233,6 @@ function evaluate{T<:Number,S<:NumberNotSeries}(a::TaylorN{T},
    suma = Taylor1(zeros(R, ord))

    for homPol in 1:length(a)
-        sun = zero(R)


Good catch!

blas-ko · 2017-10-06T19:50:39Z

Great! I've read all the changes and I'll make them as you propose. I wanted you to see the warning for set_variables in the TaylorN case, but I agree it's better to use an @assert. I'm gonna push this changes in a moment, so you can re-review.

lbenet · 2017-10-06T21:14:21Z

Incidentally, I just noticed this line, which seems to me unnecessary, since there is a method real(T::Type) (in Base at complex.jl:97) which does exactly this.

PerezHz · 2017-10-06T22:55:42Z

Incidentally, I just noticed this line, which seems to me unnecessary, since there is a method real(T::Type) (in Base at complex.jl:97) which does exactly this.

I added it in order to keep support for 0.5; once we drop 0.5 support, we can erase this line

* It doesn't use set_variables() anymore. * typeof(x0) is preserved if possible.

* Added 1 more test...

coveralls · 2017-10-07T21:40:08Z

Coverage decreased (-0.7%) to 95.288% when pulling 9d71757 on blas-ko:taylor_expand into 2111846 on JuliaDiff:master.

blas-ko · 2017-10-07T21:50:45Z

Seems that for Julia 0.5 broadcasting between Arrays and Tuples is not compatible...

Tests for HomogeneousPolynomial and TaylorN: Error During Test
  Test threw an exception of type MethodError
  Expression: taylor_expand(f11,1.0,2.0) == taylor_expand(f22,[1,2.0])
  MethodError: no method matching .+(::Tuple{Float64,Float64}, ::Array{TaylorSeries.TaylorN{Float64},1})

This is, in particular, for line 216 of src/other_functions.jl for function taylor_expand(f::Function, x0...; order::Int64=get_order()) method.
What can be done here?

coveralls · 2017-10-07T21:57:45Z

Coverage decreased (-0.6%) to 95.366% when pulling 33a23a1 on blas-ko:taylor_expand into 2111846 on JuliaDiff:master.

lbenet · 2017-10-08T01:16:57Z

src/other_functions.jl

+end
+
+function taylor_expand(f::Function, x0...; order::Int64=get_order()) #a Taylor expansion around x0
+    T = eltype(x0[1])


Since x0 is a tuple, each entry could have a different type. I think that T should be defined as T = eltype( promote(x0...)[1] ).

coveralls · 2017-10-08T01:53:42Z

Coverage decreased (-0.4%) to 95.516% when pulling f8046ba on blas-ko:taylor_expand into 2111846 on JuliaDiff:master.

blas-ko · 2017-10-08T01:57:24Z

Tests are passing!

lbenet · 2017-10-08T01:58:51Z

Great! Thanks for all! I'm merging right away!

PerezHz · 2017-10-08T02:18:52Z

Really cool stuff! Thanks for fixing #119 @blas-ko!

blas-ko added 4 commits August 30, 2017 11:29

Added function and macro taylor_expand for expanding arbitraty func…

37fd143

…tions.

Deleted macros and set order as keyword argument so methods doesnt clash

61d4f70

`taylor_expand(f,x0::Int64)` and `taylor_expand(f,order)` clashed.

Added some tests for taylor_expand.

2333386

un-exported taylor_expand macros.

65b5596

blas-ko added 3 commits September 5, 2017 11:53

Merge branch 'taylor_expand' of https://github.com/blas-ko/TaylorSeri…

13617c5

…es.jl into taylor_expand

Updated taylor_expand tests.

af194ab

Added taylor_expand method for TaylorN and warning if number of varia…

3c914ae

…bles is changed.

blas-ko mentioned this pull request Sep 7, 2017

Add method for evaluating a TaylorN on a Vector{TaylorN} #119

Closed

Fix travis issue for taylor_expand (hopefully).

49d2e4e

blas-ko and others added 7 commits October 1, 2017 18:36

Added taylor_expand! and tests for Taylor1.

29082dd

Followed suggestion from @lbenet.

Added method for evaluating a TaylorN with an array of TaylorNs.

097034d

Added taylor_expand! method for TaylorN

3cbe55a

Added some test for taylor_expand!

0abd32c

Documentation for taylor_expand and taylor_expand!

fe7d8b7

Merge branch 'master' into taylor_expand

5d978eb

blas-ko added 2 commits October 6, 2017 11:58

Added function and macro taylor_expand for expanding arbitraty func…

8e46714

…tions.

Deleted macros and set order as keyword argument so methods doesnt clash

db185d4

`taylor_expand(f,x0::Int64)` and `taylor_expand(f,order)` clashed.

blas-ko added 4 commits October 6, 2017 12:05

Added some test for taylor_expand!

7291204

Documentation for taylor_expand and taylor_expand!

fdf6cd5

Merge branch 'taylor_expand' of https://github.com/blas-ko/TaylorSeri…

f0b9518

…es.jl into taylor_expand

Corrected silly mistake from rebasing.

8bd2bff

lbenet mentioned this pull request Oct 6, 2017

Update stepsize function using TaylorSeries' norm PerezHz/TaylorIntegration.jl#40

Merged

blas-ko force-pushed the taylor_expand branch 2 times, most recently from 696c375 to 9cf86ab Compare October 6, 2017 18:25

Rearanged taylor_expand tests to another place.

f1bf4c6

They use set_variables internally.

lbenet reviewed Oct 6, 2017

View reviewed changes

blas-ko added 3 commits October 7, 2017 15:40

called coeff_table directly so it doesn't make a copy.

a90808c

changed docs for taylor_expand

3acf771

Changed taylor_expand for TaylorN.

9d71757

* It doesn't use set_variables() anymore. * typeof(x0) is preserved if possible.

blas-ko force-pushed the taylor_expand branch from 9cf86ab to 9d71757 Compare October 7, 2017 21:19

Changed taylor_expand! for update!

33a23a1

* Added 1 more test...

lbenet reviewed Oct 8, 2017

View reviewed changes

Little performance and compatibility fix for taylor_expand

f8046ba

lbenet merged commit 909689b into JuliaDiff:master Oct 8, 2017

This was referenced Oct 8, 2017

update! functionality #80

Closed

Implement an expand function? #68

Closed

taylor_expand function for Taylor1 and TaylorN #121

taylor_expand function for Taylor1 and TaylorN #121

Conversation

blas-ko commented Aug 31, 2017

coveralls commented Aug 31, 2017 • edited Loading

lbenet commented Aug 31, 2017

lbenet commented Aug 31, 2017

lbenet commented Aug 31, 2017

blas-ko commented Aug 31, 2017

lbenet commented Aug 31, 2017

lbenet commented Sep 1, 2017

blas-ko commented Sep 5, 2017 • edited Loading

coveralls commented Sep 5, 2017 • edited Loading

blas-ko commented Sep 8, 2017

coveralls commented Oct 2, 2017 • edited Loading

lbenet commented Oct 6, 2017

blas-ko commented Oct 6, 2017

lbenet left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

blas-ko commented Oct 6, 2017

lbenet commented Oct 6, 2017

PerezHz commented Oct 6, 2017

coveralls commented Oct 7, 2017 • edited Loading

blas-ko commented Oct 7, 2017

coveralls commented Oct 7, 2017 • edited Loading

Choose a reason for hiding this comment

coveralls commented Oct 8, 2017 • edited Loading

blas-ko commented Oct 8, 2017

lbenet commented Oct 8, 2017

PerezHz commented Oct 8, 2017

coveralls commented Aug 31, 2017 •

edited

Loading

blas-ko commented Sep 5, 2017 •

edited

Loading

coveralls commented Sep 5, 2017 •

edited

Loading

coveralls commented Oct 2, 2017 •

edited

Loading

coveralls commented Oct 7, 2017 •

edited

Loading

coveralls commented Oct 7, 2017 •

edited

Loading

coveralls commented Oct 8, 2017 •

edited

Loading