Compare C and C++ speed

compare.cc

@@ -216,32 +216,35 @@ int main(int argc, char *argv[]) {
 
     long cpptime_nsec, best_cpp;
     long ctime_nsec, best_c;
-
+    long cpptime_sec, ctime_sec;
+    long repeats = 209;
     //HeapProfilerStart("heapprof");    
-    for (int i = 0; i<150000; ++i) {
-      //for (int i = 0; i<100; ++i) {
-      // clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time1);
-      // nMie(L, x, m, nt, Theta, &Qext, &Qsca, &Qabs, &Qbk, &Qpr, &g, &Albedo, S1, S2);
-      // clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time2);
-      // //ctime_nsec = std::min(ctime_nsec,diff(time1,time2).tv_nsec);
-      // ctime_nsec = diff(time1,time2).tv_nsec;
-      // if (ctime_nsec < best_c) best_c = ctime_nsec;
-      
-      // clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time1);
-      nmie::nMie_wrapper(L, x, m, nt, Theta, &Qextw, &Qscaw, &Qabsw, &Qbkw, &Qprw, &gw, &Albedow, S1w, S2w);
-    // clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time2);
-      // //cpptime_nsec = std::min(cpptime_nsec, diff(time1,time2).tv_nsec);
-      // cpptime_nsec = diff(time1,time2).tv_nsec;
-      // if (cpptime_nsec < best_cpp) best_cpp = cpptime_nsec;
-      //printf("-- C++ time consumed %ld sec : %ld nsec\n",diff(time1,time2).tv_sec,cpptime_nsec);
-      
-      
-      //printf("-- C/C++ time ratio: %Lg\n", static_cast<long double>(ctime_nsec)/static_cast<long double>(cpptime_nsec));
-      //printf("--best C/C++ time ratio: %Lg\n", static_cast<long double>(best_c)/static_cast<long double>(best_cpp));
-    }  
-    //HeapProfilerStop();
-
-    //printf("--best C/C++ time ratio: %Lg\n", static_cast<long double>(best_c)/static_cast<long double>(best_cpp));
+    do {
+      clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time1);
+      for (int i = 0; i<repeats; ++i) {
+	nmie::nMie_wrapper(L, x, m, nt, Theta, &Qextw, &Qscaw,
+			   &Qabsw, &Qbkw, &Qprw, &gw, &Albedow, S1w, S2w);
+      }
+      clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time2);
+      cpptime_nsec = diff(time1,time2).tv_nsec;
+      cpptime_sec = diff(time1,time2).tv_sec;
+      clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time1);
+      for (int i = 0; i<repeats; ++i) {      
+	nMie(L, x, m, nt, Theta, &Qext, &Qsca, &Qabs, &Qbk, &Qpr, &g, &Albedo, S1, S2);
+      }  
+      clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time2);
+      ctime_nsec = diff(time1,time2).tv_nsec;
+      ctime_sec = diff(time1,time2).tv_sec;
+      if ( ctime_sec == 0 && cpptime_sec == 0) {
+	printf("-- total repeats: %ld\n", repeats);
+	printf("-- C++ time consumed %lg sec\n", (cpptime_nsec/1e9));
+	printf("-- C time consumed %lg sec\n", (ctime_nsec/1e9));
+	printf("-- C/C++ time ratio: %Lg\n", static_cast<long double>(ctime_nsec)/static_cast<long double>(cpptime_nsec));
+      } else {
+        printf("==Test is too long!\n");
+      }
+      repeats *= 10;
+    } while (cpptime_nsec < 1e8 && ctime_nsec < 1e8);
 
     // clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &time1);
     nMie(L, x, m, nt, Theta, &Qext, &Qsca, &Qabs, &Qbk, &Qpr, &g, &Albedo, S1, S2);

nmie-wrapper.cc

@@ -838,46 +838,28 @@ c    MM       + 1  and - 1, alternately
   //**********************************************************************************//
   void MultiLayerMie::calcPiTau(std::vector< std::vector<double> >& Pi,
 				std::vector< std::vector<double> >& Tau) {
-    
-  int n;
-  //****************************************************//
-  // Equations (26a) - (26c)                            //
-  //****************************************************//
-  // Initialize Pi and Tau
-  Pi[0] = 1.0;
-  Tau[0] = cos(Theta);
-  // Calculate the actual values
-  if (nmax > 1) {
-    Pi[1] = 3*Tau[0]*Pi[0];
-    Tau[1] = 2*Tau[0]*Pi[1] - 3*Pi[0];
-    for (n = 2; n < nmax; n++) {
-      Pi[n] = ((n + n + 1)*Tau[0]*Pi[n - 1] - (n + 1)*Pi[n - 2])/n;
-      Tau[n] = (n + 1)*Tau[0]*Pi[n] - (n + 2)*Pi[n - 1];
+    //****************************************************//
+    // Equations (26a) - (26c)                            //
+    //****************************************************//
+    std::vector<double> costheta(theta_.size(), 0.0);
+    for (int t = 0; t < theta_.size(); t++) {	
+      costheta[t] = cos(theta_[t]);
+    }
+    for (int n = 0; n < nmax_; n++) {
+      for (int t = 0; t < theta_.size(); t++) {	
+  	if (n == 0) {
+  	  // Initialize Pi and Tau
+  	  Pi[n][t] = 1.0;
+  	  Tau[n][t] = (n + 1)*costheta[t]; 
+  	} else {
+  	  // Calculate the actual values
+  	  Pi[n][t] = ((n == 1) ? ((n + n + 1)*costheta[t]*Pi[n - 1][t]/n)
+  		   : (((n + n + 1)*costheta[t]*Pi[n - 1][t]
+  		       - (n + 1)*Pi[n - 2][t])/n));
+  	  Tau[n][t] = (n + 1)*costheta[t]*Pi[n][t] - (n + 2)*Pi[n - 1][t];
+  	}
+      }
     }
-  }
-
-  // //****************************************************//
-  //   // Equations (26a) - (26c)                            //
-  //   //****************************************************//
-  //   std::vector<double> costheta(theta_.size(), 0.0);
-  //   for (int t = 0; t < theta_.size(); t++) {	
-  //     costheta[t] = cos(theta_[t]);
-  //   }
-  //   for (int n = 0; n < nmax_; n++) {
-  //     for (int t = 0; t < theta_.size(); t++) {	
-  // 	if (n == 0) {
-  // 	  // Initialize Pi and Tau
-  // 	  Pi[n][t] = 1.0;
-  // 	  Tau[n][t] = (n + 1)*costheta[t]; 
-  // 	} else {
-  // 	  // Calculate the actual values
-  // 	  Pi[n][t] = ((n == 1) ? ((n + n + 1)*costheta[t]*Pi[n - 1][t]/n)
-  // 		   : (((n + n + 1)*costheta[t]*Pi[n - 1][t]
-  // 		       - (n + 1)*Pi[n - 2][t])/n));
-  // 	  Tau[n][t] = (n + 1)*costheta[t]*Pi[n][t] - (n + 2)*Pi[n - 1][t];
-  // 	}
-  //     }
-  //   }
   }  // end of void MultiLayerMie::calcPiTau(...)
   //**********************************************************************************//
   // This function calculates the scattering coefficients required to calculate       //