Updated g_cluster

[alexxy/gromacs.git] / src / gromacs / gmxana / cmat.c

diff --git a/src/gromacs/gmxana/cmat.c b/src/gromacs/gmxana/cmat.c
index 5fba0be20760857b83dbdbf4ec8e05691b55e236..788d4e791d2c8b53d4b235c6d33e5a51f4f32f1a 100644 (file)
--- a/src/gromacs/gmxana/cmat.c
+++ b/src/gromacs/gmxana/cmat.c
@@ -39,6 +39,7 @@
 #include "cmat.h"
 #include "smalloc.h"
 #include "macros.h"
+#include "vec.h"
 #include "xvgr.h"
 #include "matio.h"
 #include "futil.h"
@@ -51,11 +52,9 @@ t_mat *init_mat(int n1, gmx_bool b1D)
     m->n1     = n1;
     m->nn     = 0;
     m->b1D    = b1D;
-    m->emat   = 0;
     m->maxrms = 0;
     m->minrms = 1e20;
     m->sumrms = 0;
-    m->nn     = 0;
     m->mat    = mk_matrix(n1, n1, b1D);
 
     snew(m->erow, n1);
@@ -65,6 +64,29 @@ t_mat *init_mat(int n1, gmx_bool b1D)
     return m;
 }
 
+void copy_t_mat(t_mat *dst, t_mat *src)
+{
+    int i, j;
+
+    if (dst->nn != src->nn)
+    {
+        fprintf(stderr, "t_mat structures not identical in size dst %d src %d\n",dst->nn,src->nn);
+        return;
+    }
+    dst->maxrms = src->maxrms;
+    dst->minrms = src->minrms;
+    dst->sumrms = src->sumrms;
+    for(i = 0; (i < src->nn); i++)
+    {
+        for(j = 0; (j < src->nn); j++)
+        {
+            dst->mat[i][j] = src->mat[i][j];
+        }
+        dst->erow[i] = src->erow[i];
+        dst->m_ind[i] = src->m_ind[i];
+    }
+}
+
 void enlarge_mat(t_mat *m, int deltan)
 {
     int i, j;
@@ -122,50 +144,39 @@ void done_mat(t_mat **m)
     *m = NULL;
 }
 
-real row_energy(int nn, int row, real *mat)
-{
-    real re = 0;
-    int  i;
-
-    for (i = 0; (i < nn); i++)
-    {
-        re += abs(i-row)*mat[i];
-    }
-    return re/nn;
-}
-
 real mat_energy(t_mat *m)
 {
-    real re, retot;
-    int  j, jj;
+    int  j;
+    real emat = 0;
 
-    retot = 0;
-    for (j = 0; (j < m->nn); j++)
+    for (j = 0; (j < m->nn-1); j++)
     {
-        jj         = m->m_ind[j];
-        re         = row_energy(m->nn, jj, m->mat[j]);
-        m->erow[j] = re;
-        retot     += re;
+        emat += sqr(m->mat[j][j+1]);
     }
-    m->emat = retot/m->nn;
-    return m->emat;
+    return emat;
 }
 
-void swap_rows(t_mat *m, int isw, int jsw)
+void swap_rows(t_mat *m, int iswap, int jswap)
 {
     real *tmp, ttt;
-    int   i;
+    int   i, itmp;
+
+    /* Swap indices */
+    itmp            = m->m_ind[iswap];
+    m->m_ind[iswap] = m->m_ind[jswap];
+    m->m_ind[jswap] = itmp;
+
+    /* Swap rows (since the matrix is an array of pointers) */
+    tmp           = m->mat[iswap];
+    m->mat[iswap] = m->mat[jswap];
+    m->mat[jswap] = tmp;
 
-    /* Swap rows */
-    tmp         = m->mat[isw];
-    m->mat[isw] = m->mat[jsw];
-    m->mat[jsw] = tmp;
     /* Swap columns */
     for (i = 0; (i < m->nn); i++)
     {
-        ttt            = m->mat[isw][i];
-        m->mat[isw][i] = m->mat[jsw][i];
-        m->mat[jsw][i] = ttt;
+        ttt              = m->mat[i][iswap];
+        m->mat[i][iswap] = m->mat[i][jswap];
+        m->mat[i][jswap] = ttt;
     }
 }