Merge branch 'master' of github.com:ovidiopr/scattnlay

nmie.cc

@@ -422,21 +422,22 @@ void calcD1D3(std::complex<double> z, int nmax,
 		      std::vector<std::complex<double> >& D3) {
 
   int n;
-  std::vector<std::complex<double> > PsiZeta;
-  PsiZeta.resize(nmax + 1);
+  std::complex<double> nz, PsiZeta;
 
   // Downward recurrence for D1 - equations (16a) and (16b)
   D1[nmax] = std::complex<double>(0.0, 0.0);
   for (n = nmax; n > 0; n--) {
-    D1[n - 1] = double(n)/z - 1.0/(D1[n] + double(n)/z);
+    nz = double(n)/z;
+    D1[n - 1] = nz - 1.0/(D1[n] + nz);
   }
 
   // Upward recurrence for PsiZeta and D3 - equations (18a) - (18d)
-  PsiZeta[0] = 0.5*(1.0 - std::complex<double>(cos(2.0*z.real()), sin(2.0*z.real()))*exp(-2.0*z.imag()));
+  PsiZeta = 0.5*(1.0 - std::complex<double>(cos(2.0*z.real()), sin(2.0*z.real()))*exp(-2.0*z.imag()));
   D3[0] = std::complex<double>(0.0, 1.0);
   for (n = 1; n <= nmax; n++) {
-    PsiZeta[n] = PsiZeta[n - 1]*(double(n)/z - D1[n - 1])*(double(n)/z- D3[n - 1]);
-    D3[n] = D1[n] + std::complex<double>(0.0, 1.0)/PsiZeta[n];
+    nz = double(n)/z;
+    PsiZeta = PsiZeta*(nz - D1[n - 1])*(nz - D3[n - 1]);
+    D3[n] = D1[n] + std::complex<double>(0.0, 1.0)/PsiZeta;
   }
 }