[cig-commits] r5981 - in mc/3D/CitcomS/branches/compressible: CitcomS/Components CitcomS/Components/Advection_diffusion lib module visual

tan2 at geodynamics.org tan2 at geodynamics.org
Wed Feb 7 16:19:06 PST 2007


Author: tan2
Date: 2007-02-07 16:19:05 -0800 (Wed, 07 Feb 2007)
New Revision: 5981

Modified:
   mc/3D/CitcomS/branches/compressible/CitcomS/Components/Advection_diffusion/Advection_diffusion.py
   mc/3D/CitcomS/branches/compressible/CitcomS/Components/IC.py
   mc/3D/CitcomS/branches/compressible/lib/Advection_diffusion.c
   mc/3D/CitcomS/branches/compressible/lib/Convection.c
   mc/3D/CitcomS/branches/compressible/lib/Full_version_dependent.c
   mc/3D/CitcomS/branches/compressible/lib/Initial_temperature.c
   mc/3D/CitcomS/branches/compressible/lib/Regional_version_dependent.c
   mc/3D/CitcomS/branches/compressible/lib/Solver_multigrid.c
   mc/3D/CitcomS/branches/compressible/module/setProperties.c
   mc/3D/CitcomS/branches/compressible/visual/execpaste.py
Log:
Porting r5956 to r5980 from trunk

Modified: mc/3D/CitcomS/branches/compressible/CitcomS/Components/Advection_diffusion/Advection_diffusion.py
===================================================================
--- mc/3D/CitcomS/branches/compressible/CitcomS/Components/Advection_diffusion/Advection_diffusion.py	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/CitcomS/Components/Advection_diffusion/Advection_diffusion.py	2007-02-08 00:19:05 UTC (rev 5981)
@@ -90,9 +90,11 @@
 
         fixed_timestep = prop.float("fixed_timestep", default=0.0)
         finetunedt = prop.float("finetunedt", default=0.9)
+        adv_gamma = prop.float("adv_gamma", default=0.5)
+        adv_sub_iterations = prop.int("adv_sub_iterations", default=2)
+
         inputdiffusivity = prop.float("inputdiffusivity", default=1)
 
-        adv_sub_iterations = prop.int("adv_sub_iterations", default=2)
 
 
 

Modified: mc/3D/CitcomS/branches/compressible/CitcomS/Components/IC.py
===================================================================
--- mc/3D/CitcomS/branches/compressible/CitcomS/Components/IC.py	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/CitcomS/Components/IC.py	2007-02-08 00:19:05 UTC (rev 5981)
@@ -106,9 +106,9 @@
         zero_elapsed_time = pyre.inventory.bool("zero_elapsed_time", default=True)
 
         tic_method = pyre.inventory.int("tic_method", default=0,
-                            validator=pyre.inventory.choice([0, 1, 2]))
+                            validator=pyre.inventory.choice([0, 1, 2, 3]))
 
-        # for tic_method=0
+        # for tic_method=0 or 3
         num_perturbations = pyre.inventory.int("num_perturbations", default=1,
                             validator=pyre.inventory.less(255))
         perturbl = pyre.inventory.list("perturbl", default=[1])

Modified: mc/3D/CitcomS/branches/compressible/lib/Advection_diffusion.c
===================================================================
--- mc/3D/CitcomS/branches/compressible/lib/Advection_diffusion.c	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/lib/Advection_diffusion.c	2007-02-08 00:19:05 UTC (rev 5981)
@@ -195,7 +195,8 @@
     input_int("maxtotstep",&(E->advection.max_total_timesteps),"1000000",m);
     input_float("finetunedt",&(E->advection.fine_tune_dt),"0.9",m);
     input_float("fixed_timestep",&(E->advection.fixed_timestep),"0.0",m);
-    input_int("adv_sub_iterations",&(E->advection.temp_iterations),"2,2,nomax",m);
+    input_float("adv_gamma",&(E->advection.gamma),"0.5",m);
+    input_int("adv_sub_iterations",&(E->advection.temp_iterations),"2,1,nomax",m);
 
     input_float("inputdiffusivity",&(E->control.inputdiff),"1.0",m);
 
@@ -311,7 +312,7 @@
     const int dims=E->mesh.nsd;
     const int dofs=E->mesh.dof;
     const int ends=enodes[dims];
-    const int sphere_key = 0;
+    const int sphere_key = 1;
 
     /* adiabatic and dissipative heating*/
     process_heating(E);
@@ -325,10 +326,14 @@
 
             velo_from_element(E,VV,m,el,sphere_key);
 
-            get_global_shape_fn(E,el,&GN,&GNx,&dOmega,0,sphere_key,rtf,E->mesh.levmax,m);
+            get_global_shape_fn(E, el, &GN, &GNx, &dOmega, 0,
+                                sphere_key, rtf, E->mesh.levmax, m);
 
-            pg_shape_fn(E,el,&PG,&GNx,VV,diff,m);
-            element_residual(E,el,PG,GNx,dOmega,VV,T,Tdot,Q0,Eres,diff,E->sphere.cap[m].TB,FLAGS,m);
+            pg_shape_fn(E, el, &PG, &GNx, VV,
+                        rtf, diff, m);
+            element_residual(E, el, PG, GNx, dOmega, VV, T, Tdot,
+                             Q0, Eres, rtf, diff, E->sphere.cap[m].TB,
+                             FLAGS, m);
 
             for(a=1;a<=ends;a++) {
                 a1 = E->ien[m][el].node[a];
@@ -356,31 +361,25 @@
    Petrov-Galerkin shape functions for a given element
    =================================================== */
 
-void pg_shape_fn(E,el,PG,GNx,VV,diffusion,m)
+void pg_shape_fn(E,el,PG,GNx,VV,rtf,diffusion,m)
      struct All_variables *E;
      int el,m;
      struct Shape_function *PG;
      struct Shape_function_dx *GNx;
      float VV[4][9];
+     double rtf[4][9];
      double diffusion;
 
 {
-    int i,j,node;
+    int i,j;
     int *ienm;
 
     double uc1,uc2,uc3;
-    double u1,u2,u3;
-    double aa,bb,cc,uxse,ueta,ufai,xse,eta,fai,dx1,dx2,dx3,adiff,rr1;
+    double u1,u2,u3,sint[9];
+    double uxse,ueta,ufai,xse,eta,fai,adiff;
 
     double prod1,unorm,twodiff;
 
-    const int dims=E->mesh.nsd;
-    const int dofs=E->mesh.dof;
-    const int lev=E->mesh.levmax;
-    const int nno=E->lmesh.nno;
-    const int ends=enodes[E->mesh.nsd];
-    const int vpts=vpoints[E->mesh.nsd];
-
     ienm=E->ien[m][el].node;
 
     twodiff = 2.0*diffusion;
@@ -393,64 +392,22 @@
       uc3 +=  E->N.ppt[GNPINDEX(i,1)]*VV[3][i];
       }
 
-    dx1=0.25*(E->x[m][1][ienm[3]]+E->x[m][1][ienm[4]]
-             +E->x[m][1][ienm[7]]+E->x[m][1][ienm[8]]
-             -E->x[m][1][ienm[1]]-E->x[m][1][ienm[2]]
-             -E->x[m][1][ienm[5]]-E->x[m][1][ienm[6]]);
-    dx2=0.25*(E->x[m][2][ienm[3]]+E->x[m][2][ienm[4]]
-             +E->x[m][2][ienm[7]]+E->x[m][2][ienm[8]]
-             -E->x[m][2][ienm[1]]-E->x[m][2][ienm[2]]
-             -E->x[m][2][ienm[5]]-E->x[m][2][ienm[6]]);
-    dx3=0.25*(E->x[m][3][ienm[3]]+E->x[m][3][ienm[4]]
-             +E->x[m][3][ienm[7]]+E->x[m][3][ienm[8]]
-             -E->x[m][3][ienm[1]]-E->x[m][3][ienm[2]]
-             -E->x[m][3][ienm[5]]-E->x[m][3][ienm[6]]);
-    uxse = fabs(uc1*dx1+uc2*dx2+uc3*dx3);
+    uxse = fabs(uc1*E->eco[m][el].size[1]);
+    ueta = fabs(uc2*E->eco[m][el].size[2]);
+    ufai = fabs(uc3*E->eco[m][el].size[3]);
 
-    dx1=0.25*(E->x[m][1][ienm[2]]+E->x[m][1][ienm[3]]
-             +E->x[m][1][ienm[6]]+E->x[m][1][ienm[7]]
-             -E->x[m][1][ienm[1]]-E->x[m][1][ienm[4]]
-             -E->x[m][1][ienm[5]]-E->x[m][1][ienm[8]]);
-    dx2=0.25*(E->x[m][2][ienm[2]]+E->x[m][2][ienm[3]]
-             +E->x[m][2][ienm[6]]+E->x[m][2][ienm[7]]
-             -E->x[m][2][ienm[1]]-E->x[m][2][ienm[4]]
-             -E->x[m][2][ienm[5]]-E->x[m][2][ienm[8]]);
-    dx3=0.25*(E->x[m][3][ienm[2]]+E->x[m][3][ienm[3]]
-             +E->x[m][3][ienm[6]]+E->x[m][3][ienm[7]]
-             -E->x[m][3][ienm[1]]-E->x[m][3][ienm[4]]
-             -E->x[m][3][ienm[5]]-E->x[m][3][ienm[8]]);
-    ueta = fabs(uc1*dx1+uc2*dx2+uc3*dx3);
-
-    dx1=0.25*(E->x[m][1][ienm[5]]+E->x[m][1][ienm[6]]
-             +E->x[m][1][ienm[7]]+E->x[m][1][ienm[8]]
-             -E->x[m][1][ienm[1]]-E->x[m][1][ienm[2]]
-             -E->x[m][1][ienm[3]]-E->x[m][1][ienm[4]]);
-    dx2=0.25*(E->x[m][2][ienm[5]]+E->x[m][2][ienm[6]]
-             +E->x[m][2][ienm[7]]+E->x[m][2][ienm[8]]
-             -E->x[m][2][ienm[1]]-E->x[m][2][ienm[2]]
-             -E->x[m][2][ienm[3]]-E->x[m][2][ienm[4]]);
-    dx3=0.25*(E->x[m][3][ienm[5]]+E->x[m][3][ienm[6]]
-             +E->x[m][3][ienm[7]]+E->x[m][3][ienm[8]]
-             -E->x[m][3][ienm[1]]-E->x[m][3][ienm[2]]
-             -E->x[m][3][ienm[3]]-E->x[m][3][ienm[4]]);
-    ufai = fabs(uc1*dx1+uc2*dx2+uc3*dx3);
-
-/*    xse = (uxse>twodiff)? (1.0-twodiff/uxse):0.0;
+    xse = (uxse>twodiff)? (1.0-twodiff/uxse):0.0;
     eta = (ueta>twodiff)? (1.0-twodiff/ueta):0.0;
     fai = (ufai>twodiff)? (1.0-twodiff/ufai):0.0;
-*/
 
-    aa = 2.0*uxse/twodiff;
-    bb = 2.0*ueta/twodiff;
-    cc = 2.0*ufai/twodiff;
-    xse = (1.0+exp(-aa))/(1.0-exp(-aa))-twodiff/uxse;
-    eta = (1.0+exp(-bb))/(1.0-exp(-bb))-twodiff/ueta;
-    fai = (1.0+exp(-cc))/(1.0-exp(-cc))-twodiff/ufai;
 
     unorm = uc1*uc1 + uc2*uc2 + uc3*uc3;
 
     adiff = (unorm>0.000001)?( (uxse*xse+ueta*eta+ufai*fai)/(2.0*unorm) ):0.0;
 
+    for(i=1;i<=VPOINTS3D;i++)
+       sint[i] = rtf[3][i]/sin(rtf[1][i]);
+
     for(i=1;i<=VPOINTS3D;i++) {
        u1 = u2 = u3 = 0.0;
        for(j=1;j<=ENODES3D;j++)  /* this line heavily used */ {
@@ -460,9 +417,9 @@
 	    }
 
        for(j=1;j<=ENODES3D;j++) {
-            prod1 = (u1 * GNx->vpt[GNVXINDEX(0,j,i)] +
-                     u2 * GNx->vpt[GNVXINDEX(1,j,i)] +
-                    u3 * GNx->vpt[GNVXINDEX(2,j,i)] ) ;
+            prod1 = (u1 * GNx->vpt[GNVXINDEX(0,j,i)]*rtf[3][i] +
+                     u2 * GNx->vpt[GNVXINDEX(1,j,i)]*sint[i] +
+                     u3 * GNx->vpt[GNVXINDEX(2,j,i)] ) ;
 
 	    PG->vpt[GNVINDEX(j,i)] = E->N.vpt[GNVINDEX(j,i)] + adiff * prod1;
 	    }
@@ -478,7 +435,7 @@
    Used to correct the Tdot term.
    =========================================  */
 
-void element_residual(E,el,PG,GNx,dOmega,VV,field,fielddot,Q0,Eres,diff,BC,FLAGS,m)
+void element_residual(E,el,PG,GNx,dOmega,VV,field,fielddot,Q0,Eres,rtf,diff,BC,FLAGS,m)
      struct All_variables *E;
      int el,m;
      struct Shape_function PG;
@@ -488,6 +445,7 @@
      double **field,**fielddot;
      struct SOURCES Q0;
      double Eres[9];
+     double rtf[4][9];
      double diff;
      float **BC;
      unsigned int **FLAGS;
@@ -496,7 +454,7 @@
     int i,j,a,k,node,nodes[5],d,aid,back_front,onedfns;
     double Q;
     double dT[9];
-    double tx1[9],tx2[9],tx3[9];
+    double tx1[9],tx2[9],tx3[9],sint[9];
     double v1[9],v2[9],v3[9];
     double adv_dT,t2[4];
     double T,DT;
@@ -523,6 +481,9 @@
       v3[i] = tx3[i]=  0.0;
       }
 
+    for(i=1;i<=vpts;i++)
+        sint[i] = rtf[3][i]/sin(rtf[1][i]);
+
     for(j=1;j<=ends;j++)       {
       node = E->ien[m][el].node[j];
       T = field[m][node];
@@ -533,8 +494,8 @@
 
       for(i=1;i<=vpts;i++)  {
 		  dT[i] += DT * E->N.vpt[GNVINDEX(j,i)];
-		  tx1[i] += GNx.vpt[GNVXINDEX(0,j,i)] * T;
-		  tx2[i] += GNx.vpt[GNVXINDEX(1,j,i)] * T;
+		  tx1[i] += GNx.vpt[GNVXINDEX(0,j,i)] * T * rtf[3][i];
+		  tx2[i] += GNx.vpt[GNVXINDEX(1,j,i)] * T * sint[i];
 	 	  tx3[i] += GNx.vpt[GNVXINDEX(2,j,i)] * T;
 		  sfn = E->N.vpt[GNVINDEX(j,i)];
 		  v1[i] += VV[1][j] * sfn;
@@ -560,8 +521,8 @@
 	  Eres[j] -=
 	    PG.vpt[GNVINDEX(j,i)] * dOmega.vpt[i]
 	    * (dT[i] - Q + v1[i]*tx1[i] + v2[i]*tx2[i] + v3[i]*tx3[i])
- 	    + diff*dOmega.vpt[i] * (GNx.vpt[GNVXINDEX(0,j,i)]*tx1[i] +
-				    GNx.vpt[GNVXINDEX(1,j,i)]*tx2[i] +
+ 	    + diff*dOmega.vpt[i] * (GNx.vpt[GNVXINDEX(0,j,i)]*tx1[i]*rtf[3][i] +
+				    GNx.vpt[GNVXINDEX(1,j,i)]*tx2[i]*sint[i] +
 				    GNx.vpt[GNVXINDEX(2,j,i)]*tx3[i] );
       }
     }

Modified: mc/3D/CitcomS/branches/compressible/lib/Convection.c
===================================================================
--- mc/3D/CitcomS/branches/compressible/lib/Convection.c	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/lib/Convection.c	2007-02-08 00:19:05 UTC (rev 5981)
@@ -48,11 +48,9 @@
     void convection_initial_fields();
     void twiddle_thumbs();
 
-    E->advection.temp_iterations = 2; /* petrov-galerkin iterations: minimum value. */
     E->advection.total_timesteps = 1;
     E->advection.sub_iterations = 1;
     E->advection.last_sub_iterations = 1;
-    E->advection.gamma = 0.5;
     E->advection.dt_reduced = 1.0;
 
     E->monitor.T_maxvaried = 1.05;

Modified: mc/3D/CitcomS/branches/compressible/lib/Full_version_dependent.c
===================================================================
--- mc/3D/CitcomS/branches/compressible/lib/Full_version_dependent.c	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/lib/Full_version_dependent.c	2007-02-08 00:19:05 UTC (rev 5981)
@@ -149,11 +149,7 @@
           fprintf(E->fp,"(Nodal_mesh.c #1) Cannot open %s\n",output_file);
           exit(8);
 	}
-      fscanf(fp1,"%s%d",a,&i);
-      if (i != E->mesh.noz ) {
-          fprintf(E->fp,"(Nodal_mesh.c #2) inconsistent file length: %s\n",output_file);
-          exit(8);
-      }
+      fscanf(fp1,"%s %d",a,&i);
       for (k=1;k<=E->mesh.noz;k++)  {
 	fscanf(fp1,"%d %f",&nn,&tt1);
 	rr[k]=tt1;
@@ -303,6 +299,45 @@
 	  }
     }
   }
+  else if (E->convection.tic_method == 3) {
+    /* set up a linear temperature profile first */
+    for(m=1;m<=E->sphere.caps_per_proc;m++)
+      for(i=1;i<=noy;i++)
+	for(j=1;j<=nox;j++)
+	  for(k=1;k<=noz;k++) {
+	    node=k+(j-1)*noz+(i-1)*nox*noz;
+	    r1=E->sx[m][3][node];
+	    E->T[m][node] = E->control.TBCbotval - (E->control.TBCtopval + E->control.TBCbotval)*(r1 - E->sphere.ri)/(E->sphere.ro - E->sphere.ri);
+	  }
+
+    /* This part put a temperature anomaly for whole mantle. The horizontal
+       pattern of the anomaly is given by spherical harmonic ll & mm. */
+
+    for (p=0; p<E->convection.number_of_perturbations; p++) {
+      mm = E->convection.perturb_mm[p];
+      ll = E->convection.perturb_ll[p];
+      con = E->convection.perturb_mag[p];
+      kk = E->convection.load_depth[p];
+
+      if ( (kk < 1) || (kk >= gnoz) ) continue;
+
+      if (E->parallel.me == 0)
+	fprintf(stderr,"Initial temperature perturbation:  layer=%d  mag=%g  l=%d  m=%d\n", kk, con, ll, mm);
+
+      for(m=1;m<=E->sphere.caps_per_proc;m++)
+	for(i=1;i<=noy;i++)
+	  for(j=1;j<=nox;j++)
+            for(k=1;k<=noz;k++) {
+	      node=k+(j-1)*noz+(i-1)*nox*noz;
+	      t1=E->sx[m][1][node];
+	      f1=E->sx[m][2][node];
+	      r1=E->sx[m][3][node];
+              E->T[m][node] += con*modified_plgndr_a(ll,mm,t1)
+                  *(cos(mm*f1)+sin(mm*f1))
+                  *sin(M_PI*(r1-E->sphere.ri)/(E->sphere.ro-E->sphere.ri));
+	  }
+    }
+  }
   else if (E->convection.tic_method == 1) {
 
   }

Modified: mc/3D/CitcomS/branches/compressible/lib/Initial_temperature.c
===================================================================
--- mc/3D/CitcomS/branches/compressible/lib/Initial_temperature.c	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/lib/Initial_temperature.c	2007-02-08 00:19:05 UTC (rev 5981)
@@ -1,6 +1,6 @@
 /*
  *~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
- * 
+ *
  *<LicenseText>
  *
  * CitcomS by Louis Moresi, Shijie Zhong, Lijie Han, Eh Tan,
@@ -22,7 +22,7 @@
  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  *
  *</LicenseText>
- * 
+ *
  *~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  */
 
@@ -64,10 +64,14 @@
         - blob_center: A comma-separated list of three float numbers.
         - blob_radius: A dmensionless length, typically a fraction
                        of the Earth's radius.
-        - blob_dT    : Dimensionless temperature.            */
+        - blob_dT    : Dimensionless temperature.
 
+     When tic_method is 3, the temperature is a linear profile + perturbation
+     for whole mantle.
+  */
 
-  if (E->convection.tic_method == 0) {
+
+  if (E->convection.tic_method == 0 || E->convection.tic_method == 3 ) {
     /* This part put a temperature anomaly at depth where the global
        node number is equal to load_depth. The horizontal pattern of
        the anomaly is given by spherical harmonic ll & mm. */

Modified: mc/3D/CitcomS/branches/compressible/lib/Regional_version_dependent.c
===================================================================
--- mc/3D/CitcomS/branches/compressible/lib/Regional_version_dependent.c	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/lib/Regional_version_dependent.c	2007-02-08 00:19:05 UTC (rev 5981)
@@ -437,7 +437,45 @@
                       E->T[m][node] += amp * exp(-1.0*distance/radius);
                 }
   }
+  else if (E->convection.tic_method == 3) {
+    /* set up a linear temperature profile first */
+    for(m=1;m<=E->sphere.caps_per_proc;m++)
+      for(i=1;i<=noy;i++)
+	for(j=1;j<=nox;j++)
+	  for(k=1;k<=noz;k++) {
+	    node=k+(j-1)*noz+(i-1)*nox*noz;
+	    r1=E->sx[m][3][node];
+	    E->T[m][node] = E->control.TBCbotval - (E->control.TBCtopval + E->control.TBCbotval)*(r1 - E->sphere.ri)/(E->sphere.ro - E->sphere.ri);
+	  }
 
+    /* This part put a temperature anomaly for whole mantle. The horizontal
+       pattern of the anomaly is given by spherical harmonic ll & mm. */
+
+    for (p=0; p<E->convection.number_of_perturbations; p++) {
+      mm = E->convection.perturb_mm[p];
+      ll = E->convection.perturb_ll[p];
+      con = E->convection.perturb_mag[p];
+      kk = E->convection.load_depth[p];
+
+      if ( (kk < 1) || (kk >= gnoz) ) continue;
+
+      if (E->parallel.me == 0)
+	fprintf(stderr,"Initial temperature perturbation:  layer=%d  mag=%g  l=%d  m=%d\n", kk, con, ll, mm);
+
+      for(m=1;m<=E->sphere.caps_per_proc;m++)
+	for(i=1;i<=noy;i++)
+	  for(j=1;j<=nox;j++)
+            for(k=1;k<=noz;k++) {
+	      node=k+(j-1)*noz+(i-1)*nox*noz;
+	      t1=E->sx[m][1][node];
+	      f1=E->sx[m][2][node];
+	      r1=E->sx[m][3][node];
+              E->T[m][node] += con*(cos(mm*f1)+sin(mm*f1))
+                  *sin(M_PI*(r1-E->sphere.ri)/(E->sphere.ro-E->sphere.ri));
+	  }
+    }
+  }
+
   temperatures_conform_bcs(E);
 
   return;

Modified: mc/3D/CitcomS/branches/compressible/lib/Solver_multigrid.c
===================================================================
--- mc/3D/CitcomS/branches/compressible/lib/Solver_multigrid.c	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/lib/Solver_multigrid.c	2007-02-08 00:19:05 UTC (rev 5981)
@@ -256,12 +256,9 @@
     sl_minus = lv -1;
 
     if (E->viscosity.smooth_cycles==1)  {
-/*      visc_from_gint_to_nodes(E,E->EVI[lv],viscU,lv);
+      visc_from_gint_to_nodes(E,E->EVI[lv],viscU,lv);
       project_scalar(E,lv,viscU,viscD);
-      visc_from_nodes_to_gint(E,viscD,E->EVI[sl_minus],sl_minus);  */
-
-      project_scalar(E,lv,E->VI[lv],E->VI[sl_minus]);
-      visc_from_nodes_to_gint(E,E->VI[sl_minus],E->EVI[sl_minus],sl_minus);
+      visc_from_nodes_to_gint(E,viscD,E->EVI[sl_minus],sl_minus);
       }
     else if (E->viscosity.smooth_cycles==2)   {
       visc_from_gint_to_ele(E,E->EVI[lv],viscU,lv);

Modified: mc/3D/CitcomS/branches/compressible/module/setProperties.c
===================================================================
--- mc/3D/CitcomS/branches/compressible/module/setProperties.c	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/module/setProperties.c	2007-02-08 00:19:05 UTC (rev 5981)
@@ -98,11 +98,12 @@
 
     getFloatProperty(properties, "finetunedt", E->advection.fine_tune_dt, fp);
     getFloatProperty(properties, "fixed_timestep", E->advection.fixed_timestep, fp);
+    getFloatProperty(properties, "adv_gamma", E->advection.gamma, fp);
+    getIntProperty(properties, "adv_sub_iterations", E->advection.temp_iterations, fp);
+
     getFloatProperty(properties, "inputdiffusivity", E->control.inputdiff, fp);
 
-    getIntProperty(properties, "adv_sub_iterations", E->advection.temp_iterations, fp);
 
-
     PUTS(("\n"));
 
     Py_INCREF(Py_None);
@@ -237,7 +238,7 @@
 
     getIntProperty(properties, "tic_method", E->convection.tic_method, fp);
 
-    if (E->convection.tic_method == 0) {
+    if (E->convection.tic_method == 0 || E->convection.tic_method == 3) {
 	int num_perturb;
 
 	getIntProperty(properties, "num_perturbations", num_perturb, fp);

Modified: mc/3D/CitcomS/branches/compressible/visual/execpaste.py
===================================================================
--- mc/3D/CitcomS/branches/compressible/visual/execpaste.py	2007-02-08 00:15:03 UTC (rev 5980)
+++ mc/3D/CitcomS/branches/compressible/visual/execpaste.py	2007-02-08 00:19:05 UTC (rev 5981)
@@ -55,12 +55,18 @@
 
 ## expand datadir
 s = "%HOSTNAME"
-if s in datadir:
+try:
+    datadir.index(s)
+except: pass
+else:
     from socket import gethostname
     datadir = datadir.replace(s, gethostname())
 
 s = "%RANK"
-if s in datadir:
+try:
+    datadir.index(s)
+except: pass
+else:
     datadir = datadir.replace(s, rank)
 
 if datadir == "%DATADIR":



More information about the cig-commits mailing list