jyapayne/Nim
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Fixed some minor stuff in module poly

Browse source 
Removed the stupid initPolyFromDegree which only served ro re-allocate
results. Also fixed some minor stuff with nil return values in roots.
This commit is contained in:
Robert Persson 2013-07-07 01:50:10 +02:00
commit 5a18100816
1 changed files with 9 additions and 15 deletions
Download patch file Download diff file Expand all files Collapse all files

24
lib/pure/poly.nim
View file

@ -16,13 +16,6 @@ type
cofs:seq[float] cofs:seq[float]
proc initPolyFromDegree(n:int):TPoly=
## internal usage only
## caller must initialize coefficients of poly
## and possibly `clean` away zero exponents
var numcof=n+1 #num. coefficients is one more than degree
result.cofs=newSeq[float](numcof)
proc degree*(p:TPoly):int= proc degree*(p:TPoly):int=
## Returns the degree of the polynomial, ## Returns the degree of the polynomial,
## that is the number of coefficients-1 ## that is the number of coefficients-1
@ -130,7 +123,7 @@ proc diff*(p:TPoly,x:float):float=
proc integral*(p:TPoly):TPoly= proc integral*(p:TPoly):TPoly=
## Returns a new polynomial which is the indefinite ## Returns a new polynomial which is the indefinite
## integral of `p`. The constant term is set to 0.0 ## integral of `p`. The constant term is set to 0.0
result=initPolyFromDegree(p.degree+1) newSeq(result.cofs,p.cofs.len+1)
result.cofs[0]=0.0 #constant arbitrary term, use 0.0 result.cofs[0]=0.0 #constant arbitrary term, use 0.0
for i in 1..high(result.cofs): for i in 1..high(result.cofs):
result.cofs[i]=p.cofs[i-1]/float(i) result.cofs[i]=p.cofs[i-1]/float(i)
@ -227,7 +220,7 @@ proc `*` *(p1:TPoly,p2:TPoly):TPoly=
proc `*` *(p:TPoly,f:float):TPoly= proc `*` *(p:TPoly,f:float):TPoly=
## Multiplies the polynomial `p` with a real number ## Multiplies the polynomial `p` with a real number
result=initPolyFromDegree(p.degree) newSeq(result.cofs,p.cofs.len)
for i in 0..high(p.cofs): for i in 0..high(p.cofs):
result[i]=p.cofs[i]*f result[i]=p.cofs[i]*f
result.clean result.clean
@ -254,7 +247,7 @@ proc `-` *(p1:TPoly,p2:TPoly):TPoly=
proc `/`*(p:TPoly,f:float):TPoly= proc `/`*(p:TPoly,f:float):TPoly=
## Divides polynomial `p` with a real number `f` ## Divides polynomial `p` with a real number `f`
result=initPolyFromDegree(p.degree) newSeq(result.cofs,p.cofs.len)
for i in 0..high(p.cofs): for i in 0..high(p.cofs):
result[i]=p.cofs[i]/f result[i]=p.cofs[i]/f
result.clean result.clean
@ -351,13 +344,12 @@ proc roots*(p:TPoly,tol=1.0e-9,zerotol=1.0e-6,mergetol=1.0e-12,maxiter=1000):seq
## `maxiter` can be used to limit the number of iterations for each root. ## `maxiter` can be used to limit the number of iterations for each root.
## Returns a (possibly empty) sorted sequence with the solutions. ## Returns a (possibly empty) sorted sequence with the solutions.
var deg=p.degree var deg=p.degree
result= @[] if deg<=0: #constant only => no roots
if deg<=0: return @[]
return nil elif p.degree==1: #linear
elif p.degree==1:
var linrt= -p.cofs[0]/p.cofs[1] var linrt= -p.cofs[0]/p.cofs[1]
if linrt==inf or linrt==neginf: if linrt==inf or linrt==neginf:
return nil #constant only => no roots return @[] #constant only => no roots
return @[linrt] return @[linrt]
elif p.degree==2: elif p.degree==2:
return solveQuadric(p.cofs[2],p.cofs[1],p.cofs[0],zerotol) return solveQuadric(p.cofs[2],p.cofs[1],p.cofs[0],zerotol)
@ -366,8 +358,10 @@ proc roots*(p:TPoly,tol=1.0e-9,zerotol=1.0e-6,mergetol=1.0e-12,maxiter=1000):seq
# derivative and do a numerical search for root between each min/max # derivative and do a numerical search for root between each min/max
var range=p.getRangeForRoots() var range=p.getRangeForRoots()
var minmax=p.derivative.roots(tol,zerotol,mergetol) var minmax=p.derivative.roots(tol,zerotol,mergetol)
result= @[]
if minmax!=nil: #ie. we have minimas/maximas in this function if minmax!=nil: #ie. we have minimas/maximas in this function
for x in minmax.items: for x in minmax.items:
addRoot(p,result,range.xmin,x,tol,zerotol,mergetol,maxiter) addRoot(p,result,range.xmin,x,tol,zerotol,mergetol,maxiter)
range.xmin=x range.xmin=x
addRoot(p,result,range.xmin,range.xmax,tol,zerotol,mergetol,maxiter) addRoot(p,result,range.xmin,range.xmax,tol,zerotol,mergetol,maxiter)