reworked on RVec

diff --git a/src/domains.cpp b/src/domains.cpp
index 6ecb70b9d61d8ae151306f7ae112fbd51ff71249..f8efe64f944f508dfa1b1474104de1f071f51a21 100644 (file)
--- a/src/domains.cpp
+++ b/src/domains.cpp
@@ -63,7 +63,7 @@ struct node {
     bool yep;
 };
 
-void make_graph(int mgwi_natoms, rvec *mgwi_x, std::vector< std::vector< node > > &mgwi_graph)
+void make_graph(int mgwi_natoms, std::vector < RVec > mgwi_x, std::vector< std::vector< node > > &mgwi_graph)
 {
     mgwi_graph.resize(mgwi_natoms);
     for (int i = 0; i < mgwi_natoms; i++) {
@@ -71,39 +71,19 @@ void make_graph(int mgwi_natoms, rvec *mgwi_x, std::vector< std::vector< node >
     }
     for (int i = 0; i < mgwi_natoms; i++) {
         for (int j = 0; j < mgwi_natoms; j++) {
-            rvec_sub(mgwi_x[i], mgwi_x[j], mgwi_graph[i][j].r);
+            rvec_sub(mgwi_x[i].as_vec(), mgwi_x[j].as_vec(), mgwi_graph[i][j].r.as_vec());
             mgwi_graph[i][j].n = 0;
         }
     }
 }
 
-void update_graph(std::vector< std::vector< node > > &ugwi_graph, rvec *ugwi_x, long double ugwi_epsi) {
-    rvec ugwi_temp;
+void update_graph(std::vector< std::vector< node > > &ugwi_graph, std::vector < RVec > ugwi_x, long double ugwi_epsi) {
+    RVec ugwi_temp;
     int ugwi_for = ugwi_graph.size();
     for (int i = 0; i < ugwi_for; i++) {
         for (int j = i; j < ugwi_for; j++) {
-            rvec_sub(ugwi_x[i], ugwi_x[j], ugwi_temp);
-            rvec_dec(ugwi_temp, ugwi_graph[i][j].r.as_vec());
-            if (norm(ugwi_temp) <= ugwi_epsi) {
-                if (i == j) {
-                    ugwi_graph[i][j].n++;
-                }
-                else {
-                    ugwi_graph[i][j].n++;
-                    ugwi_graph[j][i].n++;
-                }
-            }
-        }
-    }
-}
-
-void update_graph_v2(std::vector< std::vector< node > > &ugwi_graph, std::vector < RVec > ugwi_x, long double ugwi_epsi) {
-    rvec ugwi_temp;
-    int ugwi_for = ugwi_graph.size();
-    for (int i = 0; i < ugwi_for; i++) {
-        for (int j = i; j < ugwi_for; j++) {
-            rvec_sub(ugwi_x[i].as_vec(), ugwi_x[j].as_vec(), ugwi_temp);
-            rvec_dec(ugwi_temp, ugwi_graph[i][j].r.as_vec());
+            rvec_sub(ugwi_x[i].as_vec(), ugwi_x[j].as_vec(), ugwi_temp.as_vec());
+            rvec_dec(ugwi_temp.as_vec(), ugwi_graph[i][j].r.as_vec());
             if (norm(ugwi_temp) <= ugwi_epsi) {
                 if (i == j) {
                     ugwi_graph[i][j].n++;
@@ -293,7 +273,6 @@ class Domains : public TrajectoryAnalysisModule
         int                                                         frames              = 0;
         int                                                         domain_min_size     = 5; // selectable
 
-        real                                                        **w_rls;
         std::vector< std::vector< std::pair< int, int > > >         w_rls2;
         int                                                         bone;
         double                                                      delta               = 0.90; // selectable
@@ -371,34 +350,15 @@ Domains::initAfterFirstFrame(const TrajectoryAnalysisSettings       &settings,
         trajectory[0][i] = fr.x[index[i]];
     }
 
-
-
-
     bone = index.size() - domain_min_size + 1;
     graph.resize(bone);
-    snew(w_rls, bone);
     w_rls2.resize(bone);
     for (int i = 0; i < bone; i++) {
-        snew(w_rls[i], index.size());
-        for (int j = 0; j < index.size(); j++) {
-            if (j >= i && j <= i + domain_min_size - 1) {
-                w_rls[i][j] = 1;
-            } else {
-                w_rls[i][j] = 0;
-            }
-        }
         w_rls2[i].resize(0);
         for (int j = 0; j < domain_min_size; j++) {
             w_rls2[i].push_back(std::make_pair(j + i, j + i));
         }
-        rvec *etalon;
-        snew(etalon, index.size());
-        for (int j = 0; j < index.size(); j++) {
-            copy_rvec(trajectory[0][j].as_vec(), etalon[j]);
-        }
-        reset_x(index.size(), NULL, index.size(), NULL, etalon, w_rls[i]);
-        make_graph(index.size(), etalon, graph[i]);
-        sfree(etalon);
+        make_graph(index.size(), trajectory[0], graph[i]);
     }
     trajectory[1].resize(index.size());
 }
@@ -410,34 +370,17 @@ Domains::analyzeFrame(int frnr, const t_trxframe &fr, t_pbc *pbc,
     for (int i = 0; i < index.size(); i++) {
         trajectory[1][i] = fr.x[index[i]];
     }
-    std::vector < RVec > temp;
     frames++;
 
-    #pragma omp parallel
-    {
-        #pragma omp for schedule(dynamic)
+    //#pragma omp parallel
+    //{
+        //#pragma omp for schedule(dynamic)
         for (int j = 0; j < bone; j++) {
-            /*rvec *etalon, *traj;
-            snew(etalon, index.size());
-            for (int k = 0; k < index.size(); k++) {
-                copy_rvec(trajectory[0][k].as_vec(), etalon[k]);
-            }
-            snew(traj, index.size());
-            for (int k = 0; k < index.size(); k++) {
-                copy_rvec(trajectory[1][k].as_vec(), traj[k]);
-            }
-            reset_x(index.size(), NULL, index.size(), NULL, etalon, w_rls[j]);
-            reset_x(index.size(), NULL, index.size(), NULL, traj, w_rls[j]);
-            do_fit(index.size(), w_rls[j], etalon, traj);
-            update_graph(graph[j], traj, epsi);
-            sfree(etalon);
-            sfree(traj);*/
-
             std::vector < RVec > temp = trajectory[1];
             MyFitNew(trajectory[0], temp, w_rls2[j]);
-            update_graph_v2(graph[j], temp, epsi);
+            update_graph(graph[j], temp, epsi);
         }
-    }
+    //}
     std::cout << "frame: " << frames << " analyzed\n";
 }
 
@@ -467,10 +410,6 @@ Domains::finishAnalysis(int nframes)
         domsizes.resize(0);
         find_domain_sizes(graph, domsizes);
     }
-    for (int i = 0; i < bone; i++) {
-        sfree(w_rls[i]);
-    }
-    sfree(w_rls);
 }
 
 void
@@ -489,33 +428,3 @@ main(int argc, char *argv[])
 {
     return gmx::TrajectoryAnalysisCommandLineRunner::runAsMain<Domains>(argc, argv);
 }
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-