@@ -8,6 +8,7 @@ export Pade, padeval
8
8
9
9
import Base: length, endof, getindex, setindex!, copy, zero, one, convert
10
10
import Base: show, print, * , / , // , - , + , == , divrem, rem, eltype
11
+ import Base: promote_rule
11
12
12
13
eps {T} (:: Type{T} ) = zero (T)
13
14
eps {F<:FloatingPoint} (x:: Type{F} ) = Base. eps (F)
24
25
Poly {T<:Number} (a:: Vector{T} , var:: Union(String,Symbol,Char) = :x ) = Poly {T} (a, var)
25
26
26
27
convert {T} (:: Type{Poly{T}} , p:: Poly ) = Poly (convert (Vector{T}, p. a), p. var)
28
+ convert {T, S<:Number} (:: Type{Poly{T}} , x:: S ) = Poly (promote_type (T, S)[x])
29
+ convert {T, S<:Number,n} (:: Type{Poly{T}} , x:: Array{S,n} ) = map (el-> convert (Poly{promote_type (T,S)},el),x)
27
30
promote_rule {T, S} (:: Type{Poly{T}} , :: Type{Poly{S}} ) = Poly{promote_type (T, S)}
28
31
eltype {T} (:: Poly{T} ) = T
29
32
30
33
length (p:: Poly ) = length (p. a)
31
34
endof (p:: Poly ) = length (p) - 1
32
35
33
36
getindex {T} (p:: Poly{T} , i) = (i+ 1 > length (p. a) ? zero (T) : p. a[i+ 1 ])
34
- function setindex! (p:: Poly , v, i)
37
+ function setindex! (p:: Poly , v, i)
35
38
n = length (p. a)
36
39
if n < i+ 1
37
40
resize! (p. a,i+ 1 )
43
46
44
47
copy (p:: Poly ) = Poly (copy (p. a), p. var)
45
48
46
- zero {T} (p:: Poly{T} ) = zero (Poly{T} )
49
+ zero {T} (p:: Poly{T} ) = Poly ([ zero (T)], p . var )
47
50
zero {T} (:: Type{Poly{T}} ) = Poly (T[])
48
- one {T} (p:: Poly{T} ) = one (Poly{T} )
51
+ one {T} (p:: Poly{T} ) = Poly ([ one (T)], p . var )
49
52
one {T} (:: Type{Poly{T}} ) = Poly ([one (T)])
50
53
51
54
function show (io:: IO , p:: Poly )
@@ -91,10 +94,10 @@ function printterm{T<:Complex}(io::IO,p::Poly{T},j,first)
91
94
end
92
95
93
96
# We show a negative sign either for any complex number with negative
94
- # real part (and then negate the immaginary part) of for complex
97
+ # real part (and then negate the immaginary part) of for complex
95
98
# numbers that are pure imaginary with negative imaginary part
96
99
97
- neg = ((abs_repj > 2 * eps (T)) && real (pj) < 0 ) ||
100
+ neg = ((abs_repj > 2 * eps (T)) && real (pj) < 0 ) ||
98
101
((abs_impj > 2 * eps (T)) && imag (pj) < 0 )
99
102
100
103
if first
@@ -130,16 +133,16 @@ function print{T}(io::IO, p::Poly{T})
130
133
printed_anything || print (io,zero (T))
131
134
end
132
135
133
- * {T<: Number ,S}(c:: T , p:: Poly{S} ) = Poly (c * p. a)
134
- * {T<: Number ,S}(p:: Poly{S} , c:: T ) = Poly (p. a * c)
135
- / (p:: Poly , c:: Number ) = Poly (p. a / c)
136
- - (p:: Poly ) = Poly (- p. a)
136
+ * {T<: Number ,S}(c:: T , p:: Poly{S} ) = Poly (c * p. a, p . var )
137
+ * {T<: Number ,S}(p:: Poly{S} , c:: T ) = Poly (p. a * c, p . var )
138
+ / (p:: Poly , c:: Number ) = Poly (p. a / c, p . var )
139
+ - (p:: Poly ) = Poly (- p. a, p . var )
137
140
138
141
- (p:: Poly , c:: Number ) = + (p, - c)
139
142
+ (c:: Number , p:: Poly ) = + (p, c)
140
143
function + (p:: Poly , c:: Number )
141
144
if length (p) < 1
142
- return Poly ([c,])
145
+ return Poly ([c,], p . var )
143
146
else
144
147
p2 = copy (p)
145
148
p2. a[1 ] += c;
@@ -148,7 +151,7 @@ function +(p::Poly, c::Number)
148
151
end
149
152
function - {T}(c:: Number , p:: Poly{T} )
150
153
if length (p) < 1
151
- return Poly (T[c,])
154
+ return Poly (T[c,], p . var )
152
155
else
153
156
p2 = - p;
154
157
p2. a[1 ] += c;
@@ -160,13 +163,13 @@ function +{T,S}(p1::Poly{T}, p2::Poly{S})
160
163
if p1. var != p2. var
161
164
error (" Polynomials must have same variable" )
162
165
end
163
- Poly ([p1[i] + p2[i] for i = 0 : max (length (p1),length (p2))])
166
+ Poly ([p1[i] + p2[i] for i = 0 : max (length (p1),length (p2))], p1 . var )
164
167
end
165
168
function - {T,S}(p1:: Poly{T} , p2:: Poly{S} )
166
169
if p1. var != p2. var
167
170
error (" Polynomials must have same variable" )
168
171
end
169
- Poly ([p1[i] - p2[i] for i = 0 : max (length (p1),length (p2))])
172
+ Poly ([p1[i] - p2[i] for i = 0 : max (length (p1),length (p2))], p1 . var )
170
173
end
171
174
172
175
@@ -177,7 +180,7 @@ function *{T,S}(p1::Poly{T}, p2::Poly{S})
177
180
R = promote_type (T,S)
178
181
n = degree (p1)
179
182
m = degree (p2)
180
- a = Poly (zeros (R,m+ n+ 1 ))
183
+ a = Poly (zeros (R,m+ n+ 1 ), p1 . var )
181
184
for i = 0 : n
182
185
for j = 0 : m
183
186
a[i+ j] += p1[i] * p2[j]
@@ -233,22 +236,22 @@ function ==(p1::Poly, p2::Poly)
233
236
else
234
237
for i = 1 : max (length (p1),length (p2))
235
238
if p1[i] != p2[i]
236
- return false
239
+ return false
237
240
end
238
241
end
239
242
return true
240
243
end
241
244
end
242
245
243
- function polyval {T} (p:: Poly{T} , x:: Number )
244
- R = promote_type (T, typeof (x) )
246
+ function polyval {T,S } (p:: Poly{T} , x:: S )
247
+ R = promote_type (T,S )
245
248
lenp = length (p)
246
249
if lenp == 0
247
250
return zero (R)
248
251
else
249
252
y = convert (R, p[end ])
250
253
for i = (endof (p)- 1 ): - 1 : 0
251
- y = p[i] + x. * y
254
+ y = p[i] + x* y
252
255
end
253
256
return y
254
257
end
0 commit comments