speedup, internalize libop

4 files changed, 70 insertions, 26 deletions

diff --git a/CMakeLists.txt b/CMakeLists.txt
index c579851..4ae8e51 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -6,7 +6,7 @@ INCLUDE_DIRECTORIES ( "/usr/include/eigen3" )
 
 # Use Random Library for Random Number Generation
 
-set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11")
+set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11 -fopenmp")
 
 FIND_PACKAGE( Boost COMPONENTS program_options REQUIRED )
 INCLUDE_DIRECTORIES( ${Boost_INCLUDE_DIR})
diff --git a/model/ranksvm.cpp b/model/ranksvm.cpp
index a468729..acc2f4e 100644
--- a/model/ranksvm.cpp
+++ b/model/ranksvm.cpp
@@ -15,7 +15,7 @@ int maxiter = 1; // max iteration count
 double prec=1e-10; // precision
 // conjugate gradient
 double cg_prec=1e-3; // precision
-int cg_maxiter = 500;
+int cg_maxiter = 1000;
 int ls_maxiter = 20;
 // line search
 double line_prec=1e-10; // precision
diff --git a/model/ranksvmtn.cpp b/model/ranksvmtn.cpp
index 7f292f5..ba2d5a3 100644
--- a/model/ranksvmtn.cpp
+++ b/model/ranksvmtn.cpp
@@ -9,23 +9,26 @@ using namespace Eigen;
 void cal_Dw(RidList &D,const VectorXd &w, VectorXd &Dw)
 {
     int n = D.getSize();
-    #pragma omp parallel for
+    // static chunk size of 1 to interleave the iterations
+    #pragma omp parallel for schedule(static,1)
     for (int i=0;i<n;++i)
-        Dw(i) = D.getVec(i).dot(w);
+        Dw(i) = D.getVecDot(i,w);
 }
 
 void cal_Dtw(RidList &D,const VectorXd &w, VectorXd &Dw)
 {
     int n = D.getSize();
-    #pragma omp parallel
+    int fsize = D.getfSize();
+    #pragma omp parallel shared(D,Dw)
     {
         VectorXd Dw_private = VectorXd::Zero(D.getfSize());
         #pragma omp for nowait
         for (int i=0;i<n;++i)
-            Dw_private = Dw_private + D.getVec(i)*w(i);
+            D.addVecw(i,w(i),Dw_private);
         #pragma omp critical
         {
-            Dw=Dw+Dw_private;
+            for (int i=0;i<fsize;++i)
+                Dw(i) = Dw(i) + Dw_private(i);
         }
     }
 }
@@ -55,8 +58,10 @@ int cal_Hs(RidList &D,const vector<int> &rank,const VectorXd &corr,const VectorX
             else
                 g+=Ds[rank[i+j]];
     }
-    VectorXd tmp = alpha.cwiseProduct(Ds)-gamma;
-    VectorXd res = VectorXd::Zero(D.getSize());
+    VectorXd tmp(n);
+    for (int i=0;i<n;++i)
+        tmp(i) = alpha(i)*Ds(i)-gamma(i);
+    VectorXd res = VectorXd::Zero(D.getfSize());
     cal_Dtw(D,tmp,res);
     Hs = s + C*res;
     return 0;
diff --git a/tools/dataProvider.h b/tools/dataProvider.h
index 348d15c..4dc9c41 100644
--- a/tools/dataProvider.h
+++ b/tools/dataProvider.h
@@ -4,6 +4,7 @@
 #include<Eigen/Dense>
 #include "../tools/easylogging++.h"
 #include<vector>
+#include<math.h>
 
 // TODO decide how to construct training data
 // One possible way for training data:
@@ -93,8 +94,8 @@ public:
     }
     inline std::string getQid(int x)
     {
-        int a,b,n=getqSize();
-        a=x/n;
+        int a,b,q=getqSize();
+        a=x/q;
         return getU(a)->qid;
     }
     inline int getqSize()
@@ -112,25 +113,63 @@ public:
         return getuSize()*getqSize();
     }
     inline Eigen::VectorXd getVec(int x){
-        int a,b,n=getqSize();
-        a=x/n;
-        b=x%n;
-        Eigen::VectorXd vec;
+        int a,b,q=getqSize();
+        a=x/q;
+        b=x%q;
+        Eigen::VectorXd *id,*oth;
+        id = &(uniq[a]->feature);
         if (single)
-            return (uniq[a]->feature-other[b]->feature).cwiseAbs();
-        if (b<a)
-            vec=uniq[a]->feature-uniq[b]->feature;
+            oth = &(other[b]->feature);
+        else if (b<a)
+            oth = &(uniq[b]->feature);
+        else if (b<uniq.size()-1)
+            oth = &(uniq[b+1]->feature);
         else
-        if (b<uniq.size()-1)
-            vec=uniq[a]->feature-uniq[b+1]->feature;
-        else
-            vec=uniq[a]->feature-other[b-uniq.size()+1]->feature;
-        return vec.cwiseAbs();
+            oth = &(other[b-uniq.size()+1]->feature);
+        return (*id-*oth).cwiseAbs();
     };
+    inline double getVecDot(int x,const Eigen::VectorXd &w)
+    {
+        int a,b,q=getqSize();
+        a=x/q;
+        b=x%q;
+        double res = 0;
+        Eigen::VectorXd *id,*oth;
+        id = &(uniq[a]->feature);
+        if (single)
+            oth = &(other[b]->feature);
+        else if (b<a)
+            oth = &(uniq[b]->feature);
+        else if (b<uniq.size()-1)
+            oth = &(uniq[b+1]->feature);
+        else
+            oth = &(other[b-uniq.size()+1]->feature);
+        for (int i=0;i<n;++i)
+            res += fabs((*id)[i] - (*oth)[i])*w[i];
+        return res;
+    }
+    inline void addVecw(int x,double w,Eigen::VectorXd &X)
+    {
+        int a,b,q=getqSize();
+        a=x/q;
+        b=x%q;
+        Eigen::VectorXd *id,*oth;
+        id = &(uniq[a]->feature);
+        if (single)
+            oth = &(other[b]->feature);
+        else if (b<a)
+            oth = &(uniq[b]->feature);
+        else if (b<uniq.size()-1)
+            oth = &(uniq[b+1]->feature);
+        else
+            oth = &(other[b-uniq.size()+1]->feature);
+        for (int i=0;i<n;++i)
+            X[i] += fabs((*id)[i] - (*oth)[i])*w;
+    }
     inline double getL(int x){
-        int a,b,n=getqSize();
-        a=x/n;
-        b=x%n;
+        int a,b,q=getqSize();
+        a=x/q;
+        b=x%q;
         if (single)
         {
             if (std::fabs(other[b]->rank - a) < 1e-5)