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

tan2 at geodynamics.org tan2 at geodynamics.org
Thu Mar 1 16:09:48 PST 2007


Author: tan2
Date: 2007-03-01 16:09:47 -0800 (Thu, 01 Mar 2007)
New Revision: 6152

Modified:
   mc/3D/CitcomS/trunk/lib/Full_tracer_advection.c
   mc/3D/CitcomS/trunk/lib/tracer_defs.h
Log:
* Simplified checking in the radial direction
* Removed obsolete variables in trace struct


Modified: mc/3D/CitcomS/trunk/lib/Full_tracer_advection.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Full_tracer_advection.c	2007-03-01 22:24:39 UTC (rev 6151)
+++ mc/3D/CitcomS/trunk/lib/Full_tracer_advection.c	2007-03-02 00:09:47 UTC (rev 6152)
@@ -204,7 +204,6 @@
     void write_trace_instructions();
     void viscosity_checks();
     void make_tracer_array();
-    void setup_shared_cap_information();
     void initialize_old_composition();
     void fill_composition();
     void find_tracers();
@@ -302,7 +301,7 @@
 	    parallel_process_termination();
 	}
 
-    setup_shared_cap_information(E);
+    get_neighboring_caps(E);
 
     if (E->trace.itracer_restart==0) make_tracer_array(E);
     else if (E->trace.itracer_restart==-1) read_tracer_file(E);
@@ -1814,7 +1813,6 @@
     double *REC[13];
 
     int iget_element();
-    int icheck_that_processor_shell();
     int icheck_cap();
     void expand_tracer_arrays();
 
@@ -2302,7 +2300,6 @@
 
 				    rad=REC[j][irad];
 
-
 				    ival=icheck_that_processor_shell(E,j,ithat_processor,rad);
 
 
@@ -4063,60 +4060,6 @@
 }
 
 
-/************** SETUP SHARED CAP INFORMATION ***************/
-void setup_shared_cap_information(E)
-     struct All_variables *E;
-{
-
-    int noz;
-    int j;
-    int nroot;
-    int nproc;
-
-    double diff;
-
-    static double eps=0.0001;
-
-    get_neighboring_caps(E);
-
-
-    nproc=E->parallel.nproc;
-    noz=E->lmesh.noz;
-
-    for (j=1;j<=E->sphere.caps_per_proc;j++)
-	{
-
-	    E->trace.BOTTOM_RAD[E->parallel.me][j]=E->sx[j][3][1];
-	    E->trace.TOP_RAD[E->parallel.me][j]=E->sx[j][3][noz];
-
-	    diff=fabs(E->sx[j][3][noz]-E->sphere.ro);
-
-	    if (diff<=eps) E->trace.ITOP[E->parallel.me][j]=1;
-	    else E->trace.ITOP[E->parallel.me][j]=0;
-
-
-	}
-
-    /* broadcast */
-
-    for (nroot=0;nroot<=(nproc-1);nroot++)
-	{
-	    for (j=1;j<=E->sphere.caps_per_proc;j++)
-		{
-
-
-		    MPI_Bcast(&E->trace.BOTTOM_RAD[nroot][j],1,MPI_DOUBLE,
-			      nroot,E->parallel.world);
-		    MPI_Bcast(&E->trace.TOP_RAD[nroot][j],1,MPI_DOUBLE,
-			      nroot,E->parallel.world);
-
-
-		}
-	}
-
-    return;
-}
-
 /************** RESTART TRACERS ******************************************/
 /*                                                                       */
 /* This function restarts tracers written from previous calculation      */
@@ -4548,6 +4491,7 @@
     return;
 }
 
+
 /********* ICHECK THAT PROCESSOR SHELL ********/
 /*                                            */
 /* Checks whether a given radius is within    */
@@ -4559,47 +4503,40 @@
 /* is at the surface (then both boundaries are   */
 /* included).                                    */
 
-
 int icheck_that_processor_shell(E,j,nprocessor,rad)
      struct All_variables *E;
      int j;
      int nprocessor;
      double rad;
 {
+    int icheck_processor_shell();
+    int me = E->parallel.me;
 
+    /* nprocessor is right on top of me */
+    if (nprocessor == me+1) {
+        if (icheck_processor_shell(E, j, rad) == 0) return 1;
+        else return 0;
+    }
 
-    double top_r,bottom_r;
+    /* nprocessor is right on bottom of me */
+    if (nprocessor == me-1) {
+        if (icheck_processor_shell(E, j, rad) == -99) return 1;
+        else return 0;
+    }
 
-    top_r=E->trace.TOP_RAD[nprocessor][j];
-    bottom_r=E->trace.BOTTOM_RAD[nprocessor][j];
-
-
-    /* First check bottom */
-
-    if (rad<bottom_r) return 0;
-
-    /* Check top */
-
-
-    if (E->trace.ITOP[nprocessor][j]==1)
-	{
-	    if (rad<=top_r) return 1;
-	}
-    else
-	{
-	    if (rad<top_r) return 1;
-	}
-
-    /* If here, means point is above processor */
-
-    return 0;
+    /* Shouldn't be here */
+    fprintf(E->trace.fpt, "Should not be here\n");
+    fprintf(E->trace.fpt, "Error(check_shell) nprocessor: %d, radius: %f\n",
+            nprocessor, rad);
+    fflush(E->trace.fpt);
+    exit(10);
 }
 
 
 /********** ICHECK PROCESSOR SHELL *************/
-/* returns -99 if rad is below correct shell  */
-/* returns 0 if rad is above correct shell    */
-/* returns 1 if correct shell                 */
+/* returns -99 if rad is below current shell  */
+/* returns 0 if rad is above current shell    */
+/* returns 1 if rad is within current shell   */
 /*                                            */
 /* Shell, here, refers to processor shell     */
 /*                                            */
@@ -4614,58 +4551,30 @@
 
 {
 
-    static int top=0;
-    static int bottom=0;
-    static int been_here=0;
-    double eps=0.0000001;
+    const int noz = E->lmesh.noz;
+    const int nprocz = E->parallel.nprocz;
+    double top_r, bottom_r;
 
+    if (nprocz==1) return 1;
 
-    if (E->parallel.nprocz==1) return 1;
+    top_r = E->sx[j][3][noz];
+    bottom_r = E->sx[j][3][1];
 
-    if (been_here==0)
-	{
-	    if ((E->sx[j][3][1]-eps)<E->sphere.ri) bottom=1;
-	    if ((E->sx[j][3][E->lmesh.noz]+eps)>E->sphere.ro) top=1;
-	    been_here++;
-	}
+    /* First check bottom */
 
-    if (top==1)
-	{
-	    if (rad>=E->sx[j][3][1])
-		{
-		    return 1;
-		}
-	    return -99;
-	}
+    if (rad<bottom_r) return -99;
 
-    if (bottom==1)
-	{
-	    if (rad<E->sx[j][3][E->lmesh.noz])
-		{
-		    return 1;
-		}
-	    return 0;
-	}
 
-    if ((rad>=E->sx[j][3][1])&&(rad<E->sx[j][3][E->lmesh.noz]))
-	{
-	    return 1;
-	}
-    if (rad>=E->sx[j][3][E->lmesh.noz])
-	{
-	    return 0;
-	}
-    if (rad<E->sx[j][3][1])
-	{
-	    return -99;
-	}
+    /* Check top */
 
+    if (rad<top_r) return 1;
 
-    fprintf(E->trace.fpt,"Should not be here\n");
-    fprintf(E->trace.fpt,"Error(check_shell) radius: %f\n",rad);
-    fflush(E->trace.fpt);
-    exit(10);
+    /* top processor */
 
+    if ( (rad<=top_r) && (E->parallel.me_loc[3]==nprocz-1) ) return 1;
+
+    /* If here, means point is above processor */
+    return 0;
 }
 
 /******* ICHECK ELEMENT *************************************/

Modified: mc/3D/CitcomS/trunk/lib/tracer_defs.h
===================================================================
--- mc/3D/CitcomS/trunk/lib/tracer_defs.h	2007-03-01 22:24:39 UTC (rev 6151)
+++ mc/3D/CitcomS/trunk/lib/tracer_defs.h	2007-03-02 00:09:47 UTC (rev 6152)
@@ -128,20 +128,6 @@
 double stat_stokes_time;
 
 /* Mesh information */
-double XCAP[13][5];
-double YCAP[13][5];
-double ZCAP[13][5];
-double THETA_CAP[13][5];
-double PHI_CAP[13][5];
-double RAD_CAP[13][5];
-double COS_THETA[13][5];
-double SIN_THETA[13][5];
-double COS_PHI[13][5];
-double SIN_PHI[13][5];
-double TOP_RAD[101][13];
-double BOTTOM_RAD[101][13];
-int ITOP[101][13];
-
     double xcap[13][5];
     double ycap[13][5];
     double zcap[13][5];



More information about the cig-commits mailing list