Sort all includes in src/gromacs
[alexxy/gromacs.git] / src / gromacs / linearalgebra / sparsematrix.c
index a6d78873a5d36ac6e10e7ab2e74b53ebdf4de44a..5b2405c3de99bcefac6c790027f0cd23cfd49e39 100644 (file)
@@ -1,68 +1,69 @@
 /*
- * 
- *                This source code is part of
- * 
- *                 G   R   O   M   A   C   S
- * 
- *          GROningen MAchine for Chemical Simulations
- * 
- *                        VERSION 3.2.0
- * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
+ * This file is part of the GROMACS molecular simulation package.
+ *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2004, The GROMACS development team,
- * check out http://www.gromacs.org for more information.
- * This program is free software; you can redistribute it and/or
- * modify it under the terms of the GNU General Public License
- * as published by the Free Software Foundation; either version 2
+ * Copyright (c) 2001-2004, The GROMACS development team.
+ * Copyright (c) 2012,2014, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
  * of the License, or (at your option) any later version.
- * 
- * If you want to redistribute modifications, please consider that
- * scientific software is very special. Version control is crucial -
- * bugs must be traceable. We will be happy to consider code for
- * inclusion in the official distribution, but derived work must not
- * be called official GROMACS. Details are found in the README & COPYING
- * files - if they are missing, get the official version at www.gromacs.org.
- * 
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
  * To help us fund GROMACS development, we humbly ask that you cite
- * the papers on the package - you can find them in the top README file.
- * 
- * For more info, check our website at http://www.gromacs.org
- * 
- * And Hey:
- * Green Red Orange Magenta Azure Cyan Skyblue
+ * the research papers on the package. Check out http://www.gromacs.org.
  */
+#include "gmxpre.h"
+
 #include "sparsematrix.h"
 
 #include <assert.h>
-#include <stdlib.h>
 #include <stdio.h>
+#include <stdlib.h>
 
-#include "gromacs/legacyheaders/smalloc.h"
-
-
-
+#include "gromacs/utility/smalloc.h"
 
 gmx_sparsematrix_t *
 gmx_sparsematrix_init(int                    nrow)
 {
-    int i;
+    int                 i;
     gmx_sparsematrix_t *A;
-    
-    snew(A,1);
-    
-    A->nrow=nrow;
-    snew(A->ndata,nrow);
-    snew(A->nalloc,nrow);
-    snew(A->data,nrow);
-        
-    for(i=0;i<nrow;i++) 
+
+    snew(A, 1);
+
+    A->nrow = nrow;
+    snew(A->ndata, nrow);
+    snew(A->nalloc, nrow);
+    snew(A->data, nrow);
+
+    for (i = 0; i < nrow; i++)
     {
         A->ndata[i]    = 0;
         A->nalloc[i]   = 0;
         A->data[i]     = NULL;
     }
-    return A;    
+    return A;
 }
 
 
@@ -71,19 +72,21 @@ void
 gmx_sparsematrix_destroy(gmx_sparsematrix_t *   A)
 {
     int i;
-    
+
     /* Release each row */
-    for(i=0;i<A->nrow;i++) 
+    for (i = 0; i < A->nrow; i++)
     {
-        if(A->data[i]!=NULL)
+        if (A->data[i] != NULL)
+        {
             sfree(A->data[i]);
+        }
     }
     /* Release the rowdata arrays */
     sfree(A->ndata);
     sfree(A->nalloc);
     sfree(A->data);
     /* Release matrix structure itself */
-    sfree(A);    
+    sfree(A);
 }
 
 
@@ -92,57 +95,63 @@ void
 gmx_sparsematrix_print(FILE *                 stream,
                        gmx_sparsematrix_t *   A)
 {
-    int i,j,k;
-    
-    for(i=0;i<A->nrow;i++)
+    int i, j, k;
+
+    for (i = 0; i < A->nrow; i++)
     {
-        if(A->ndata[i]==0) 
+        if (A->ndata[i] == 0)
         {
-            for(j=0;j<A->nrow;j++)
-                fprintf(stream," %6.3f",0.0);
+            for (j = 0; j < A->nrow; j++)
+            {
+                fprintf(stream, " %6.3f", 0.0);
+            }
         }
         else
         {
-            k=0;
-            j=0;
-            for(j=0;j<A->ndata[i];j++) 
+            k = 0;
+            j = 0;
+            for (j = 0; j < A->ndata[i]; j++)
             {
-                while(k++<A->data[i][j].col) 
-                    fprintf(stream," %6.3f",0.0);
-                fprintf(stream," %6.3f",A->data[i][j].value);
+                while (k++ < A->data[i][j].col)
+                {
+                    fprintf(stream, " %6.3f", 0.0);
+                }
+                fprintf(stream, " %6.3f", A->data[i][j].value);
+            }
+            while (k++ < A->nrow)
+            {
+                fprintf(stream, " %6.3f", 0.0);
             }
-            while(k++<A->nrow)
-                fprintf(stream," %6.3f",0.0);
         }
-        fprintf(stream,"\n");
+        fprintf(stream, "\n");
     }
-    
+
 }
 
 
 real
 gmx_sparsematrix_value(gmx_sparsematrix_t *    A,
-                       int                     row, 
+                       int                     row,
                        int                     col)
 {
     gmx_bool found  = FALSE;
-    int  i;
-    real value;
-    
-    assert(row<A->nrow);
+    int      i;
+    real     value;
+
+    assert(row < A->nrow);
 
     value = 0;
-    
+
     /* Find previous value */
-    for(i=0;i<A->ndata[row] && (found==FALSE);i++)
+    for (i = 0; i < A->ndata[row] && (found == FALSE); i++)
     {
-        if(A->data[row][i].col==col) 
+        if (A->data[row][i].col == col)
         {
             found = TRUE;
             value = A->data[row][i].value;
         }
     }
-    
+
     /* value=0 if we didn't find any match */
     return value;
 }
@@ -151,46 +160,46 @@ gmx_sparsematrix_value(gmx_sparsematrix_t *    A,
 
 void
 gmx_sparsematrix_increment_value(gmx_sparsematrix_t *    A,
-                                 int                     row, 
+                                 int                     row,
                                  int                     col,
                                  real                    difference)
 {
     gmx_bool found  = FALSE;
-    int i;
-    
-    assert(row<A->nrow);
-    
+    int      i;
+
+    assert(row < A->nrow);
+
     /* Try to find a previous entry with this row/col */
-    for(i=0;i<A->ndata[row] && !found;i++)
+    for (i = 0; i < A->ndata[row] && !found; i++)
     {
-        if(A->data[row][i].col==col) 
+        if (A->data[row][i].col == col)
         {
-            found = TRUE;
+            found                  = TRUE;
             A->data[row][i].value += difference;
         }
     }
-    
+
     /* Add a new entry if nothing was found */
-    if(!found) 
+    if (!found)
     {
         /* add the value at the end of the row */
-        if(A->ndata[row] == A->nalloc[row]) 
+        if (A->ndata[row] == A->nalloc[row])
         {
             A->nalloc[row] += 100;
-            if(A->data[row]==NULL)
+            if (A->data[row] == NULL)
             {
-                snew(A->data[row],A->nalloc[row]);
+                snew(A->data[row], A->nalloc[row]);
             }
             else
             {
-                srenew(A->data[row],A->nalloc[row]);
+                srenew(A->data[row], A->nalloc[row]);
             }
         }
         A->data[row][A->ndata[row]].col   = col;
         /* Previous value was 0.0 */
         A->data[row][A->ndata[row]].value = difference;
         A->ndata[row]++;
-    }     
+    }
 }
 
 
@@ -204,34 +213,42 @@ compare_columns(const void *v1, const void *v2)
 {
     int c1 = ((gmx_sparsematrix_entry_t *)v1)->col;
     int c2 = ((gmx_sparsematrix_entry_t *)v2)->col;
-    
-    if(c1<c2)
+
+    if (c1 < c2)
+    {
         return -1;
-    else if(c1>c2)
+    }
+    else if (c1 > c2)
+    {
         return 1;
-    else 
+    }
+    else
+    {
         return 0;
+    }
 }
 
 
 void
 gmx_sparsematrix_compress(gmx_sparsematrix_t *    A)
 {
-    int i,j;
-    
-    for (i=0;i<A->nrow;i++) 
-    {        
+    int i, j;
+
+    for (i = 0; i < A->nrow; i++)
+    {
         /* Remove last value on this row while it is zero */
-        while(A->ndata[i]>0 && A->data[i][A->ndata[i]-1].value==0)
+        while (A->ndata[i] > 0 && A->data[i][A->ndata[i]-1].value == 0)
+        {
             A->ndata[i]--;
-        
+        }
+
         /* Go through values on this row and look for more zero elements */
-        for(j=0;j<A->ndata[i];j++)
+        for (j = 0; j < A->ndata[i]; j++)
         {
             /* If this element was zero, exchange it with the last non-zero
              * element on the row (yes, this will invalidate the sort order)
              */
-            if(A->data[i][j].value==0)
+            if (A->data[i][j].value == 0)
             {
                 A->data[i][j].value = A->data[i][A->ndata[i]-1].value;
                 A->data[i][j].col   = A->data[i][A->ndata[i]-1].col;
@@ -243,7 +260,7 @@ gmx_sparsematrix_compress(gmx_sparsematrix_t *    A)
               A->ndata[i],
               sizeof(gmx_sparsematrix_entry_t),
               compare_columns);
-    }        
+    }
 }
 
 
@@ -252,51 +269,52 @@ gmx_sparsematrix_vector_multiply(gmx_sparsematrix_t *    A,
                                  real *                  x,
                                  real *                  y)
 {
-    real                        s,v,xi;
-    int                         i,j,k;
+    real                        s, v, xi;
+    int                         i, j, k;
     gmx_sparsematrix_entry_t *  data; /* pointer to simplify data access */
-    
-    for (i = 0; i < A->nrow; i ++)
+
+    for (i = 0; i < A->nrow; i++)
+    {
         y[i] = 0;
-    
-    if(A->compressed_symmetric)
+    }
+
+    if (A->compressed_symmetric)
     {
-        for (i = 0; i < A->nrow; i ++)
+        for (i = 0; i < A->nrow; i++)
         {
-            xi = x[i];
-            s = 0.0;
+            xi   = x[i];
+            s    = 0.0;
             data = A->data[i];
-            
-            for (k=0;k<A->ndata[i];k++)
+
+            for (k = 0; k < A->ndata[i]; k++)
             {
-                j = data[k].col;
-                v = data[k].value;
+                j  = data[k].col;
+                v  = data[k].value;
                 s += v * x[j];
-                if(i!=j)
-                    y[j] += v * xi; 
+                if (i != j)
+                {
+                    y[j] += v * xi;
+                }
             }
-            y[i] += s; 
-        }    
+            y[i] += s;
+        }
     }
     else
     {
         /* not compressed symmetric storage */
-        for (i = 0; i < A->nrow; i ++)
+        for (i = 0; i < A->nrow; i++)
         {
-            xi = x[i];
-            s = 0.0;
+            xi   = x[i];
+            s    = 0.0;
             data = A->data[i];
-            
-            for (k=0;k<A->ndata[i];k++) 
+
+            for (k = 0; k < A->ndata[i]; k++)
             {
-                j = data[k].col;
-                v = data[k].value;
+                j  = data[k].col;
+                v  = data[k].value;
                 s += v * x[j];
             }
-            y[i] += s; 
-        }    
+            y[i] += s;
+        }
     }
 }
-
-
-