[cig-commits] r14381 - in mc/3D/CitcomS/trunk: CitcomS/Components/Sphere module

[tan2 at geodynamics.org]

Author: tan2
Date: 2009-03-18 12:39:54 -0700 (Wed, 18 Mar 2009)
New Revision: 14381

Modified:
   mc/3D/CitcomS/trunk/CitcomS/Components/Sphere/Sphere.py
   mc/3D/CitcomS/trunk/module/setProperties.c
Log:
Remove mgunitx etc from pyre input. This restores the behavior in v3.0 and earlier version.

In pyre version, the mesh size is always specified by nodex etc. In C version, the mesh size is specified by nodex if Solver=cgrad, and by mgunitx and levels if Solver=multigrid. 


Modified: mc/3D/CitcomS/trunk/CitcomS/Components/Sphere/Sphere.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Components/Sphere/Sphere.py	2009-03-18 13:44:22 UTC (rev 14380)
+++ mc/3D/CitcomS/trunk/CitcomS/Components/Sphere/Sphere.py	2009-03-18 19:39:54 UTC (rev 14381)
@@ -85,10 +85,7 @@
         nodex = pyre.inventory.int("nodex", default=9)
         nodey = pyre.inventory.int("nodey", default=9)
         nodez = pyre.inventory.int("nodez", default=9)
-        mgunitx = pyre.inventory.int("mgunitx", default=2)
-        mgunity = pyre.inventory.int("mgunity", default=2)
-        mgunitz = pyre.inventory.int("mgunitz", default=2)
-        levels = pyre.inventory.int("levels", default=3)
+        levels = pyre.inventory.int("levels", default=1)
 
         radius_outer = pyre.inventory.float("radius_outer", default=1.0)
         radius_inner = pyre.inventory.float("radius_inner", default=0.55)

Modified: mc/3D/CitcomS/trunk/module/setProperties.c
===================================================================
--- mc/3D/CitcomS/trunk/module/setProperties.c	2009-03-18 13:44:22 UTC (rev 14380)
+++ mc/3D/CitcomS/trunk/module/setProperties.c	2009-03-18 19:39:54 UTC (rev 14381)
@@ -541,28 +541,42 @@
     getFloatVectorProperty(properties, "coor_refine", E->control.coor_refine, 4, fp);
     getStringProperty(properties, "coor_file", E->control.coor_file, fp);
 
+    getIntProperty(properties, "nodex", E->mesh.nox, fp);
+    getIntProperty(properties, "nodey", E->mesh.noy, fp);
+    getIntProperty(properties, "nodez", E->mesh.noz, fp);
+    getIntProperty(properties, "levels", E->mesh.levels, fp);
+
     if (E->control.CONJ_GRAD) {
-        getIntProperty(properties, "nodex", E->mesh.nox, fp);
-        getIntProperty(properties, "nodey", E->mesh.noy, fp);
-        getIntProperty(properties, "nodez", E->mesh.noz, fp);
-
+        if (E->mesh.levels != 1) {
+            fprintf(stderr, "!!! cgrad solver must have levels=1\n");
+            abort();
+        }
         E->mesh.mgunitx = E->mesh.nox - 1;
         E->mesh.mgunity = E->mesh.noy - 1;
         E->mesh.mgunitz = E->mesh.noz - 1;
-        E->mesh.levels = 1;
     }
     else {
         double levmax;
+        int nox, noy, noz;
 
-        getIntProperty(properties, "mgunitx", E->mesh.mgunitx, fp);
-        getIntProperty(properties, "mgunity", E->mesh.mgunity, fp);
-        getIntProperty(properties, "mgunitz", E->mesh.mgunitz, fp);
-        getIntProperty(properties, "levels", E->mesh.levels, fp);
+        if (E->mesh.levels <= 1) {
+            fprintf(stderr, "!!! multigrid solver must have levels>1\n");
+            abort();
+        }
+        levmax = E->mesh.levels - 1;
 
-        levmax = E->mesh.levels - 1;
-        E->mesh.nox = E->mesh.mgunitx * (int) pow(2.0,levmax) * E->parallel.nprocx + 1;
-        E->mesh.noy = E->mesh.mgunity * (int) pow(2.0,levmax) * E->parallel.nprocy + 1;
-        E->mesh.noz = E->mesh.mgunitz * (int) pow(2.0,levmax) * E->parallel.nprocz + 1;
+        E->mesh.mgunitx = (E->mesh.nox - 1) / E->parallel.nprocx / pow(2.0, levmax);
+        E->mesh.mgunity = (E->mesh.noy - 1) / E->parallel.nprocy / pow(2.0, levmax);
+        E->mesh.mgunitz = (E->mesh.noz - 1) / E->parallel.nprocz / pow(2.0, levmax);
+
+        nox = E->mesh.mgunitx * (int) pow(2.0,levmax) * E->parallel.nprocx + 1;
+        noy = E->mesh.mgunity * (int) pow(2.0,levmax) * E->parallel.nprocy + 1;
+        noz = E->mesh.mgunitz * (int) pow(2.0,levmax) * E->parallel.nprocz + 1;
+
+        if ((nox != E->mesh.nox) || (noy != E->mesh.noy) || (noz != E->mesh.noz)) {
+           fprintf(stderr, "!!! inconsistent mesh size and levels.\n");
+           abort();
+        }
     }
 
     if (E->parallel.nprocxy == 12) {