3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Green Red Orange Magenta Azure Cyan Skyblue
46 t_mat *init_mat(int n1, gmx_bool b1D)
59 m->mat = mk_matrix(n1, n1, b1D);
68 void enlarge_mat(t_mat *m, int deltan)
72 srenew(m->erow, m->nn+deltan);
73 srenew(m->m_ind, m->nn+deltan);
74 srenew(m->mat, m->nn+deltan);
76 /* Reallocate existing rows in the matrix, and set them to zero */
77 for (i = 0; (i < m->nn); i++)
79 srenew(m->mat[i], m->nn+deltan);
80 for (j = m->nn; (j < m->nn+deltan); j++)
85 /* Allocate new rows of the matrix, set energies to zero */
86 for (i = m->nn; (i < m->nn+deltan); i++)
89 snew(m->mat[i], m->nn+deltan);
94 void reset_index(t_mat *m)
98 for (i = 0; (i < m->n1); i++)
104 void set_mat_entry(t_mat *m, int i, int j, real val)
106 m->mat[i][j] = m->mat[j][i] = val;
107 m->maxrms = max(m->maxrms, val);
110 m->minrms = min(m->minrms, val);
113 m->nn = max(m->nn, max(j+1, i+1));
116 void done_mat(t_mat **m)
118 done_matrix((*m)->n1, &((*m)->mat));
125 real row_energy(int nn, int row, real *mat)
130 for (i = 0; (i < nn); i++)
132 re += abs(i-row)*mat[i];
137 real mat_energy(t_mat *m)
143 for (j = 0; (j < m->nn); j++)
146 re = row_energy(m->nn, jj, m->mat[j]);
150 m->emat = retot/m->nn;
154 void swap_rows(t_mat *m, int isw, int jsw)
161 m->mat[isw] = m->mat[jsw];
164 for (i = 0; (i < m->nn); i++)
166 ttt = m->mat[isw][i];
167 m->mat[isw][i] = m->mat[jsw][i];
168 m->mat[jsw][i] = ttt;
172 void swap_mat(t_mat *m)
177 tmp = init_mat(m->nn, FALSE);
178 for (i = 0; (i < m->nn); i++)
180 for (j = 0; (j < m->nn); j++)
182 tmp->mat[m->m_ind[i]][m->m_ind[j]] = m->mat[i][j];
185 /*tmp->mat[i][j] = m->mat[m->m_ind[i]][m->m_ind[j]]; */
186 for (i = 0; (i < m->nn); i++)
188 for (j = 0; (j < m->nn); j++)
190 m->mat[i][j] = tmp->mat[i][j];
196 void low_rmsd_dist(const char *fn, real maxrms, int nn, real **mat,
197 const output_env_t oenv)
205 for (i = 0; i < nn; i++)
207 for (j = i+1; j < nn; j++)
209 x = (int)(fac*mat[i][j]+0.5);
217 fp = xvgropen(fn, "RMS Distribution", "RMS (nm)", "a.u.", oenv);
218 for (i = 0; (i < 101); i++)
220 fprintf(fp, "%10g %10d\n", i/fac, histo[i]);
226 void rmsd_distribution(const char *fn, t_mat *rms, const output_env_t oenv)
228 low_rmsd_dist(fn, rms->maxrms, rms->nn, rms->mat, oenv);
231 t_clustid *new_clustid(int n1)
237 for (i = 0; (i < n1); i++)