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

tan2 at geodynamics.org tan2 at geodynamics.org
Mon Mar 16 15:50:06 PDT 2009


Author: tan2
Date: 2009-03-16 15:50:06 -0700 (Mon, 16 Mar 2009)
New Revision: 14352

Modified:
   mc/3D/CitcomS/trunk/CitcomS/Components/Sphere/Sphere.py
   mc/3D/CitcomS/trunk/CitcomS/Solver/Solver.py
   mc/3D/CitcomS/trunk/lib/Full_geometry_cartesian.c
   mc/3D/CitcomS/trunk/lib/Full_sphere_related.c
   mc/3D/CitcomS/trunk/lib/Instructions.c
   mc/3D/CitcomS/trunk/lib/Regional_geometry_cartesian.c
   mc/3D/CitcomS/trunk/lib/Regional_sphere_related.c
   mc/3D/CitcomS/trunk/module/bindings.c
   mc/3D/CitcomS/trunk/module/mesher.c
   mc/3D/CitcomS/trunk/module/mesher.h
   mc/3D/CitcomS/trunk/module/setProperties.c
Log:
seperate initialization and parameter input


Modified: mc/3D/CitcomS/trunk/CitcomS/Components/Sphere/Sphere.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Components/Sphere/Sphere.py	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/CitcomS/Components/Sphere/Sphere.py	2009-03-16 22:50:06 UTC (rev 14352)
@@ -31,6 +31,14 @@
 
 class Sphere(CitcomComponent):
 
+
+    def initialize(self, all_variables):
+        CitcomComponent.initialize(self, all_variables)
+        from CitcomSLib import set_3dsphere_defaults
+        set_3dsphere_defaults(self.all_variables)
+        return
+
+
     def run(self):
         start_time = CPU_time()
         self.launch()

Modified: mc/3D/CitcomS/trunk/CitcomS/Solver/Solver.py
===================================================================
--- mc/3D/CitcomS/trunk/CitcomS/Solver/Solver.py	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/CitcomS/Solver/Solver.py	2009-03-16 22:50:06 UTC (rev 14352)
@@ -78,6 +78,9 @@
         # information about clock time
         self.start_cpu_time = CPU_time()
 
+        set_signal()
+        global_default_values(self.all_variables)
+
         inv = self.inventory
 
         inv.mesher.initialize(all_variables)
@@ -93,9 +96,6 @@
         inv.tracer.initialize(all_variables)
         inv.visc.initialize(all_variables)
 
-        set_signal()
-        global_default_values(self.all_variables)
-
         from CitcomSLib import return_rank, return_pid
         rank = return_rank(self.all_variables)
         if rank == 0:

Modified: mc/3D/CitcomS/trunk/lib/Full_geometry_cartesian.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Full_geometry_cartesian.c	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/lib/Full_geometry_cartesian.c	2009-03-16 22:50:06 UTC (rev 14352)
@@ -62,65 +62,11 @@
 void full_set_3dsphere_defaults(E)
      struct All_variables *E;
 {
-  void full_set_3dsphere_defaults2(struct All_variables *);
-  int m=E->parallel.me;
-
-  input_double("radius_outer",&(E->sphere.ro),"1",m);
-  input_double("radius_inner",&(E->sphere.ri),"0.55",m);
-
-  full_set_3dsphere_defaults2(E);
-
-  return;
-}
-
-
-void full_set_3dsphere_defaults2(struct All_variables *E)
-{
-  int i,j;
-  double offset;
-
   E->mesh.nsd = 3;
   E->mesh.dof = 3;
 
   E->sphere.caps = 12;
   E->sphere.max_connections = 6;
 
-  /* adjust the corner coordinates so that the size (surface area) of
-     each cap is about the same. */
-  offset = 9.736/180.0*M_PI;
-
-  for (i=1;i<=4;i++)  {
-    E->sphere.cap[(i-1)*3+1].theta[1] = 0.0;
-    E->sphere.cap[(i-1)*3+1].theta[2] = M_PI/4.0+offset;
-    E->sphere.cap[(i-1)*3+1].theta[3] = M_PI/2.0;
-    E->sphere.cap[(i-1)*3+1].theta[4] = M_PI/4.0+offset;
-    E->sphere.cap[(i-1)*3+1].fi[1] = 0.0;
-    E->sphere.cap[(i-1)*3+1].fi[2] = (i-1)*M_PI/2.0;
-    E->sphere.cap[(i-1)*3+1].fi[3] = (i-1)*M_PI/2.0 + M_PI/4.0;
-    E->sphere.cap[(i-1)*3+1].fi[4] = i*M_PI/2.0;
-
-    E->sphere.cap[(i-1)*3+2].theta[1] = M_PI/4.0+offset;
-    E->sphere.cap[(i-1)*3+2].theta[2] = M_PI/2.0;
-    E->sphere.cap[(i-1)*3+2].theta[3] = 3*M_PI/4.0-offset;
-    E->sphere.cap[(i-1)*3+2].theta[4] = M_PI/2.0;
-    E->sphere.cap[(i-1)*3+2].fi[1] = i*M_PI/2.0;
-    E->sphere.cap[(i-1)*3+2].fi[2] = i*M_PI/2.0 - M_PI/4.0;
-    E->sphere.cap[(i-1)*3+2].fi[3] = i*M_PI/2.0;
-    E->sphere.cap[(i-1)*3+2].fi[4] = i*M_PI/2.0 + M_PI/4.0;
-    }
-
-  for (i=1;i<=4;i++)  {
-    j = (i-1)*3;
-    if (i==1) j=12;
-    E->sphere.cap[j].theta[1] = M_PI/2.0;
-    E->sphere.cap[j].theta[2] = 3*M_PI/4.0-offset;
-    E->sphere.cap[j].theta[3] = M_PI;
-    E->sphere.cap[j].theta[4] = 3*M_PI/4.0-offset;
-    E->sphere.cap[j].fi[1] = (i-1)*M_PI/2.0 + M_PI/4.0;
-    E->sphere.cap[j].fi[2] = (i-1)*M_PI/2.0;
-    E->sphere.cap[j].fi[3] = 0.0;
-    E->sphere.cap[j].fi[4] = i*M_PI/2.0;
-    }
-
   return;
 }

Modified: mc/3D/CitcomS/trunk/lib/Full_sphere_related.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Full_sphere_related.c	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/lib/Full_sphere_related.c	2009-03-16 22:50:06 UTC (rev 14352)
@@ -167,6 +167,7 @@
   double *px, *py, *qx, *qy;
   double theta, fi, cost, sint, cosf, sinf, efac2,rfac;
   double a, b;
+  double offset;
   double myatan();
 
   void even_divide_arc12();
@@ -202,8 +203,46 @@
   qx = malloc((temp+1)*sizeof(double));
   qy = malloc((temp+1)*sizeof(double));
 
+  /* define the cap corners */
 
+  /* adjust the corner coordinates so that the size (surface area) of
+     each cap is about the same. */
+  offset = 9.736/180.0*M_PI;
 
+  for (i=1;i<=4;i++)  {
+    E->sphere.cap[(i-1)*3+1].theta[1] = 0.0;
+    E->sphere.cap[(i-1)*3+1].theta[2] = M_PI/4.0+offset;
+    E->sphere.cap[(i-1)*3+1].theta[3] = M_PI/2.0;
+    E->sphere.cap[(i-1)*3+1].theta[4] = M_PI/4.0+offset;
+    E->sphere.cap[(i-1)*3+1].fi[1] = 0.0;
+    E->sphere.cap[(i-1)*3+1].fi[2] = (i-1)*M_PI/2.0;
+    E->sphere.cap[(i-1)*3+1].fi[3] = (i-1)*M_PI/2.0 + M_PI/4.0;
+    E->sphere.cap[(i-1)*3+1].fi[4] = i*M_PI/2.0;
+
+    E->sphere.cap[(i-1)*3+2].theta[1] = M_PI/4.0+offset;
+    E->sphere.cap[(i-1)*3+2].theta[2] = M_PI/2.0;
+    E->sphere.cap[(i-1)*3+2].theta[3] = 3*M_PI/4.0-offset;
+    E->sphere.cap[(i-1)*3+2].theta[4] = M_PI/2.0;
+    E->sphere.cap[(i-1)*3+2].fi[1] = i*M_PI/2.0;
+    E->sphere.cap[(i-1)*3+2].fi[2] = i*M_PI/2.0 - M_PI/4.0;
+    E->sphere.cap[(i-1)*3+2].fi[3] = i*M_PI/2.0;
+    E->sphere.cap[(i-1)*3+2].fi[4] = i*M_PI/2.0 + M_PI/4.0;
+    }
+
+  for (i=1;i<=4;i++)  {
+    j = (i-1)*3;
+    if (i==1) j=12;
+    E->sphere.cap[j].theta[1] = M_PI/2.0;
+    E->sphere.cap[j].theta[2] = 3*M_PI/4.0-offset;
+    E->sphere.cap[j].theta[3] = M_PI;
+    E->sphere.cap[j].theta[4] = 3*M_PI/4.0-offset;
+    E->sphere.cap[j].fi[1] = (i-1)*M_PI/2.0 + M_PI/4.0;
+    E->sphere.cap[j].fi[2] = (i-1)*M_PI/2.0;
+    E->sphere.cap[j].fi[3] = 0.0;
+    E->sphere.cap[j].fi[4] = i*M_PI/2.0;
+    }
+
+
   /* 4 corners of the cap in Cartesian coordinates */
   /* the order of corners is: */
   /*  1 - 4 */

Modified: mc/3D/CitcomS/trunk/lib/Instructions.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Instructions.c	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/lib/Instructions.c	2009-03-16 22:50:06 UTC (rev 14352)
@@ -372,7 +372,7 @@
   input_int("nprocy",&(E->parallel.nprocy),"1",m);
   input_int("nprocz",&(E->parallel.nprocz),"1",m);
 
-  if ( strcmp(E->control.SOLVER_TYPE,"cgrad") == 0) {
+  if (E->control.CONJ_GRAD) {
       input_int("nodex",&(E->mesh.nox),"essential",m);
       input_int("nodez",&(E->mesh.noz),"essential",m);
       input_int("nodey",&(E->mesh.noy),"essential",m);
@@ -395,6 +395,16 @@
       E->mesh.noz = E->mesh.mgunitz * (int) pow(2.0,levmax) * E->parallel.nprocz + 1;
   }
 
+  input_double("radius_outer",&(E->sphere.ro),"1",m);
+  input_double("radius_inner",&(E->sphere.ri),"0.55",m);
+
+  if(E->sphere.caps == 1) {
+      input_double("theta_min",&(E->control.theta_min),"essential",m);
+      input_double("theta_max",&(E->control.theta_max),"essential",m);
+      input_double("fi_min",&(E->control.fi_min),"essential",m);
+      input_double("fi_max",&(E->control.fi_max),"essential",m);
+  }
+
   input_int("coor",&(E->control.coor),"0",m);
   if(E->control.coor == 2){
     /*
@@ -768,7 +778,7 @@
 void check_settings_consistency(struct All_variables *E)
 {
 
-    if (strcmp(E->control.SOLVER_TYPE, "cgrad") == 0) {
+    if (E->control.CONJ_GRAD) {
         /* conjugate gradient has only one level */
         if(E->mesh.levels != 1)
             myerror(E, "Conjugate gradient solver is used. 'levels' must be 1.\n");

Modified: mc/3D/CitcomS/trunk/lib/Regional_geometry_cartesian.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Regional_geometry_cartesian.c	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/lib/Regional_geometry_cartesian.c	2009-03-16 22:50:06 UTC (rev 14352)
@@ -60,39 +60,11 @@
 void regional_set_3dsphere_defaults(E)
      struct All_variables *E;
 {
-  void regional_set_3dsphere_defaults2(struct All_variables *E);
-  int m = E->parallel.me;
-
-  input_double("radius_outer",&(E->sphere.ro),"1",m);
-  input_double("radius_inner",&(E->sphere.ri),"0.55",m);
-
-  input_double("theta_min",&(E->control.theta_min),"essential",m);
-  input_double("theta_max",&(E->control.theta_max),"essential",m);
-  input_double("fi_min",&(E->control.fi_min),"essential",m);
-  input_double("fi_max",&(E->control.fi_max),"essential",m);
-
-  regional_set_3dsphere_defaults2(E);
-
-  return;
-}
-
-
-void regional_set_3dsphere_defaults2(struct All_variables *E)
-{
   E->mesh.nsd = 3;
   E->mesh.dof = 3;
 
   E->sphere.caps = 1;
   E->sphere.max_connections = 6;
 
-  E->sphere.cap[1].theta[1] = E->control.theta_min;
-  E->sphere.cap[1].theta[2] = E->control.theta_max;
-  E->sphere.cap[1].theta[3] = E->control.theta_max;
-  E->sphere.cap[1].theta[4] = E->control.theta_min;
-  E->sphere.cap[1].fi[1] = E->control.fi_min;
-  E->sphere.cap[1].fi[2] = E->control.fi_min;
-  E->sphere.cap[1].fi[3] = E->control.fi_max;
-  E->sphere.cap[1].fi[4] = E->control.fi_max;
-
   return;
 }

Modified: mc/3D/CitcomS/trunk/lib/Regional_sphere_related.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Regional_sphere_related.c	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/lib/Regional_sphere_related.c	2009-03-16 22:50:06 UTC (rev 14352)
@@ -71,6 +71,15 @@
   nnproc=nprocyl*nprocxl*nproczl;
   temp = max(E->mesh.NOY[E->mesh.levmax],E->mesh.NOX[E->mesh.levmax]);
 
+  /* define the cap corners */
+  E->sphere.cap[1].theta[1] = E->control.theta_min;
+  E->sphere.cap[1].theta[2] = E->control.theta_max;
+  E->sphere.cap[1].theta[3] = E->control.theta_max;
+  E->sphere.cap[1].theta[4] = E->control.theta_min;
+  E->sphere.cap[1].fi[1] = E->control.fi_min;
+  E->sphere.cap[1].fi[2] = E->control.fi_min;
+  E->sphere.cap[1].fi[3] = E->control.fi_max;
+  E->sphere.cap[1].fi[4] = E->control.fi_max;
 
   if(E->control.coor==1) {
 
@@ -106,6 +115,15 @@
     
     fclose(fp);
     
+    /* redefine the cap corners */
+    E->sphere.cap[1].theta[1] = E->control.theta_min;
+    E->sphere.cap[1].theta[2] = E->control.theta_max;
+    E->sphere.cap[1].theta[3] = E->control.theta_max;
+    E->sphere.cap[1].theta[4] = E->control.theta_min;
+    E->sphere.cap[1].fi[1] = E->control.fi_min;
+    E->sphere.cap[1].fi[2] = E->control.fi_min;
+    E->sphere.cap[1].fi[3] = E->control.fi_max;
+    E->sphere.cap[1].fi[4] = E->control.fi_max;
     
     for (lev=E->mesh.gridmin;lev<=E->mesh.gridmax;lev++)  {
       

Modified: mc/3D/CitcomS/trunk/module/bindings.c
===================================================================
--- mc/3D/CitcomS/trunk/module/bindings.c	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/module/bindings.c	2009-03-16 22:50:06 UTC (rev 14352)
@@ -221,6 +221,11 @@
 
     /* from mesher.h */
 
+    {pyCitcom_set_3dsphere_defaults__name__,
+     pyCitcom_set_3dsphere_defaults,
+     METH_VARARGS,
+     pyCitcom_set_3dsphere_defaults__doc__},
+
     {pyCitcom_full_sphere_launch__name__,
      pyCitcom_full_sphere_launch,
      METH_VARARGS,

Modified: mc/3D/CitcomS/trunk/module/mesher.c
===================================================================
--- mc/3D/CitcomS/trunk/module/mesher.c	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/module/mesher.c	2009-03-16 22:50:06 UTC (rev 14352)
@@ -38,6 +38,27 @@
 
 
 
+char pyCitcom_set_3dsphere_defaults__doc__[] = "";
+char pyCitcom_set_3dsphere_defaults__name__[] = "set_3dsphere_defaults";
+
+PyObject * pyCitcom_set_3dsphere_defaults(PyObject *self, PyObject *args)
+{
+    PyObject *obj;
+    struct All_variables* E;
+
+    if (!PyArg_ParseTuple(args, "O:set_3dsphere_defaults", &obj))
+        return NULL;
+
+    E = (struct All_variables*)(PyCObject_AsVoidPtr(obj));
+
+    (E->solver.set_3dsphere_defaults)(E);
+
+    Py_INCREF(Py_None);
+    return Py_None;
+}
+
+
+
 char pyCitcom_full_sphere_launch__doc__[] = "";
 char pyCitcom_full_sphere_launch__name__[] = "full_sphere_launch";
 

Modified: mc/3D/CitcomS/trunk/module/mesher.h
===================================================================
--- mc/3D/CitcomS/trunk/module/mesher.h	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/module/mesher.h	2009-03-16 22:50:06 UTC (rev 14352)
@@ -29,6 +29,11 @@
 #define pyCitcom_mesher_h
 
 
+extern char pyCitcom_set_3dsphere_defaults__name__[];
+extern char pyCitcom_set_3dsphere_defaults__doc__[];
+PyObject * pyCitcom_set_3dsphere_defaults(PyObject *, PyObject *);
+
+
 extern char pyCitcom_full_sphere_launch__name__[];
 extern char pyCitcom_full_sphere_launch__doc__[];
 PyObject * pyCitcom_full_sphere_launch(PyObject *, PyObject *);

Modified: mc/3D/CitcomS/trunk/module/setProperties.c
===================================================================
--- mc/3D/CitcomS/trunk/module/setProperties.c	2009-03-16 22:44:42 UTC (rev 14351)
+++ mc/3D/CitcomS/trunk/module/setProperties.c	2009-03-16 22:50:06 UTC (rev 14352)
@@ -539,7 +539,7 @@
     getFloatVectorProperty(properties, "coor_refine", E->control.coor_refine, 4, fp);
     getStringProperty(properties, "coor_file", E->control.coor_file, fp);
 
-    if ( strcmp(E->control.SOLVER_TYPE,"cgrad") == 0) {
+    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);
@@ -575,16 +575,11 @@
     getDoubleProperty(properties, "radius_inner", E->sphere.ri, fp);
 
 
-    if (E->parallel.nprocxy == 12) {
-        full_set_3dsphere_defaults2(E);
-    }
-    else {
+    if (E->sphere.caps == 1) {
 	getDoubleProperty(properties, "theta_min", E->control.theta_min, fp);
 	getDoubleProperty(properties, "theta_max", E->control.theta_max, fp);
 	getDoubleProperty(properties, "fi_min", E->control.fi_min, fp);
 	getDoubleProperty(properties, "fi_max", E->control.fi_max, fp);
-
-        regional_set_3dsphere_defaults2(E);
     }
 
     E->mesh.layer[1] = 1;



More information about the CIG-COMMITS mailing list