[cig-commits] r4584 - mc/3D/CitcomS/trunk/lib

[tan2 at geodynamics.org]

Author: tan2
Date: 2006-09-21 14:17:46 -0700 (Thu, 21 Sep 2006)
New Revision: 4584

Modified:
   mc/3D/CitcomS/trunk/lib/Output.c
   mc/3D/CitcomS/trunk/lib/Output_h5.c
   mc/3D/CitcomS/trunk/lib/Process_buoyancy.c
Log:
Refactored out the code computing horizontal average of temperature and rms
velocity into a function. This fixes the bug in issue46.


Modified: mc/3D/CitcomS/trunk/lib/Output.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Output.c	2006-09-21 20:08:43 UTC (rev 4583)
+++ mc/3D/CitcomS/trunk/lib/Output.c	2006-09-21 21:17:46 UTC (rev 4584)
@@ -255,45 +255,15 @@
 void output_average(struct All_variables *E, int cycles)
 {
   /* horizontal average output of temperature and rms velocity*/
-  void return_horiz_ave_f();
+  void compute_horiz_avg();
 
-  int m, i, j;
-  float *S1[NCS],*S2[NCS],*S3[NCS];
+  int j;
   char output_file[255];
   FILE *fp1;
 
   // compute horizontal average here....
+  compute_horiz_avg(E);
 
-  for(m=1;m<=E->sphere.caps_per_proc;m++)      {
-    S1[m] = (float *)malloc((E->lmesh.nno+1)*sizeof(float));
-    S2[m] = (float *)malloc((E->lmesh.nno+1)*sizeof(float));
-    S3[m] = (float *)malloc((E->lmesh.nno+1)*sizeof(float));
-  }
-
-  for(m=1;m<=E->sphere.caps_per_proc;m++) {
-    for(i=1;i<=E->lmesh.nno;i++) {
-      S1[m][i] = E->T[m][i];
-      S2[m][i] = E->sphere.cap[m].V[1][i]*E->sphere.cap[m].V[1][i]
-          	+ E->sphere.cap[m].V[2][i]*E->sphere.cap[m].V[2][i];
-      S3[m][i] = E->sphere.cap[m].V[3][i]*E->sphere.cap[m].V[3][i];
-    }
-  }
-
-  return_horiz_ave_f(E,S1,E->Have.T);
-  return_horiz_ave_f(E,S2,E->Have.V[1]);
-  return_horiz_ave_f(E,S3,E->Have.V[2]);
-
-  for(m=1;m<=E->sphere.caps_per_proc;m++) {
-    free((void *)S1[m]);
-    free((void *)S2[m]);
-    free((void *)S3[m]);
-  }
-
-  for (i=1;i<=E->lmesh.noz;i++) {
-      E->Have.V[1][i] = sqrt(E->Have.V[1][i]);
-      E->Have.V[2][i] = sqrt(E->Have.V[2][i]);
-  }
-
   // only the first nprocz processors need to output
 
   if (E->parallel.me<E->parallel.nprocz)  {

Modified: mc/3D/CitcomS/trunk/lib/Output_h5.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Output_h5.c	2006-09-21 20:08:43 UTC (rev 4583)
+++ mc/3D/CitcomS/trunk/lib/Output_h5.c	2006-09-21 21:17:46 UTC (rev 4584)
@@ -1761,63 +1761,26 @@
 void h5output_average(struct All_variables *E, int cycles)
 {
     /* horizontal average output of temperature and rms velocity */
-    void return_horiz_ave_f();
+    void compute_horiz_avg();
 
     int cap;
     hid_t cap_group;
     hid_t dataset;
     herr_t status;
 
-    float vx, vy, vz;
-    float *S1[NCS], *S2[NCS], *S3[NCS];
-
     field_t *field;
 
     int k;
-    int n, nz;
-    int m, mz;
+    int mz;
 
 
     field = E->hdf5.scalar1d;
 
-    nz = E->lmesh.noz;
     mz = field->block[1];
 
+    /* calculate horizontal averages */
+    compute_horiz_avg(E);
 
-    /*
-     * calculate horizontal averages
-     */
-
-    S1[1] = (float *)malloc((E->lmesh.nno + 1) * sizeof(float));
-    S2[1] = (float *)malloc((E->lmesh.nno + 1) * sizeof(float));
-    S3[1] = (float *)malloc((E->lmesh.nno + 1) * sizeof(float));
-
-    for(n = 1; n < E->lmesh.nno; n++)
-    {
-        vx = E->sphere.cap[1].V[1][n];
-        vy = E->sphere.cap[1].V[2][n];
-        vz = E->sphere.cap[1].V[3][n];
-
-        S1[1][n] = E->T[1][n];
-        S2[1][n] = vx*vx + vy*vy;
-        S3[1][n] = vz*vz;
-    }
-
-    return_horiz_ave_f(E, S1, E->Have.T);
-    return_horiz_ave_f(E, S2, E->Have.V[1]);
-    return_horiz_ave_f(E, S3, E->Have.V[2]);
-
-    free((void *)S1[1]);
-    free((void *)S2[1]);
-    free((void *)S3[1]);
-
-    for(k = 1; k <= E->lmesh.noz; k++)
-    {
-        E->Have.V[1][k] = sqrt(E->Have.V[1][k]);
-        E->Have.V[2][k] = sqrt(E->Have.V[2][k]);
-    }
-
-
     /* extend all datasets -- collective I/O call */
     for(cap = 0; cap < E->sphere.caps; cap++)
     {

Modified: mc/3D/CitcomS/trunk/lib/Process_buoyancy.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Process_buoyancy.c	2006-09-21 20:08:43 UTC (rev 4583)
+++ mc/3D/CitcomS/trunk/lib/Process_buoyancy.c	2006-09-21 21:17:46 UTC (rev 4584)
@@ -185,5 +185,48 @@
   free((void *)sum_h);
 
   return;
-  }
+}
 
+
+
+/*
+  compute horizontal average of temperature and rms velocity
+*/
+void compute_horiz_avg(struct All_variables *E)
+{
+    void return_horiz_ave_f();
+
+    int m, i;
+    float *S1[NCS],*S2[NCS],*S3[NCS];
+
+    for(m=1;m<=E->sphere.caps_per_proc;m++)      {
+	S1[m] = (float *)malloc((E->lmesh.nno+1)*sizeof(float));
+	S2[m] = (float *)malloc((E->lmesh.nno+1)*sizeof(float));
+	S3[m] = (float *)malloc((E->lmesh.nno+1)*sizeof(float));
+    }
+
+    for(m=1;m<=E->sphere.caps_per_proc;m++) {
+	for(i=1;i<=E->lmesh.nno;i++) {
+	    S1[m][i] = E->T[m][i];
+	    S2[m][i] = E->sphere.cap[m].V[1][i]*E->sphere.cap[m].V[1][i]
+          	+ E->sphere.cap[m].V[2][i]*E->sphere.cap[m].V[2][i];
+	    S3[m][i] = E->sphere.cap[m].V[3][i]*E->sphere.cap[m].V[3][i];
+	}
+    }
+
+    return_horiz_ave_f(E,S1,E->Have.T);
+    return_horiz_ave_f(E,S2,E->Have.V[1]);
+    return_horiz_ave_f(E,S3,E->Have.V[2]);
+
+    for(m=1;m<=E->sphere.caps_per_proc;m++) {
+	free((void *)S1[m]);
+	free((void *)S2[m]);
+	free((void *)S3[m]);
+    }
+
+    for (i=1;i<=E->lmesh.noz;i++) {
+	E->Have.V[1][i] = sqrt(E->Have.V[1][i]);
+	E->Have.V[2][i] = sqrt(E->Have.V[2][i]);
+    }
+
+}