diff --git a/polynomial.cc b/polynomial.cc
index ab3ae59..2fbf055 100644
--- a/polynomial.cc
+++ b/polynomial.cc
@@ -264,6 +264,8 @@ template class Polynomial<int>;
 template class Polynomial<unsigned int>;
 template class Polynomial<double>;
 template class Polynomial<std::complex<double> >;
+template class Polynomial<float>;
+template class Polynomial<std::complex<float> >;
 
 #define INSTANTIZE(T) \
 template NRMat<T> Sylvester(const Polynomial<T> &p, const Polynomial<T> &q); \
@@ -277,12 +279,15 @@ INSTANTIZE(int)
 INSTANTIZE(unsigned int)
 INSTANTIZE(double)
 INSTANTIZE(std::complex<double>)
+INSTANTIZE(float)
+INSTANTIZE(std::complex<float>)
 
 #define INSTANTIZE2(T) \
 template Polynomial<T> hermite_polynomial(int); \
 
 INSTANTIZE2(int)
 INSTANTIZE2(double)
+INSTANTIZE2(float)
 
 
 
diff --git a/polynomial.h b/polynomial.h
index dd33006..f67f36d 100644
--- a/polynomial.h
+++ b/polynomial.h
@@ -86,7 +86,7 @@ public:
 		while(n>0 && MYABS((*this)[n])<=thr) --n;
 		resize(n,true);
 		};
-	void normalize() {if((*this)[degree()]==(T)0) laerror("zero coefficient at highest degree - simplify first"); *this /= (*this)[degree()];};
+	void normalize() {const T lead=((const Polynomial &)*this)[degree()]; if(lead==(T)0) laerror("zero coefficient at highest degree - simplify first"); *this /= lead;};
 	Polynomial shifted(const int shift) const
 		{
 		if(shift==0) return *this;
@@ -154,7 +154,7 @@ public:
 	Polynomial operator/(const Polynomial &rhs) const {Polynomial q,r; polydiv(rhs,q,r); return q;};
 	Polynomial operator%(const Polynomial &rhs) const {Polynomial q,r; polydiv(rhs,q,r); return r;};
 	NRMat<T> companion() const; //matrix which has this characteristic polynomial
-	NRVec<typename LA_traits<T>::complextype> roots() const; //implemented for complex<double> and double only
+	NRVec<typename LA_traits<T>::complextype> roots() const; //implemented for double, float and their complex types only
 	NRVec<T> realroots(const typename LA_traits<T>::normtype thr) const;
 	T newton(const T x0, const typename LA_traits<T>::normtype thr=LA_precision<typename LA_traits<T>::normtype>::newton_threshold(), const int maxit=1000) const; //solve root from the guess
 	Polynomial pow(const int i) const; //integer power