[cig-commits] commit: Sync the edge viscosity so that it is correct in parallel

[Mercurial]

changeset:   184:f883f3778688
tag:         tip
user:        Walter Landry <wlandry at caltech.edu>
date:        Fri Apr 22 12:21:08 2011 -0700
files:       FACStokes/fix_viscosity.C
description:
Sync the edge viscosity so that it is correct in parallel


diff -r d31e8b419587 -r f883f3778688 FACStokes/fix_viscosity.C
--- a/FACStokes/fix_viscosity.C	Thu Apr 21 02:38:19 2011 -0700
+++ b/FACStokes/fix_viscosity.C	Fri Apr 22 12:21:08 2011 -0700
@@ -102,64 +102,78 @@ void SAMRAI::FACStokes::fix_viscosity()
   edge_viscosity_coarsen_algorithm.setNull();
   edge_viscosity_coarsen_schedules.setNull();
 
+
+  /* Compute edge_viscosity by averaging the cell viscosities. */
+
   hier::Index ip(hier::Index::getZeroIndex(d_dim)), jp(ip), kp(ip);
   ip[0]=1;
   jp[1]=1;
   if(d_dim.getValue()>2)
     kp[2]=1;
 
-  for (int ln = 0; ln <= d_hierarchy->getFinestLevelNumber(); ++ln) {
-    tbox::Pointer<hier::PatchLevel> level = d_hierarchy->getPatchLevel(ln);
-    hier::PatchLevel::Iterator i_p(*level);
-    for ( ; i_p; i_p++) {
-      tbox::Pointer<hier::Patch> patch = *i_p;
-      tbox::Pointer<pdat::CellData<double> >
-        cell_viscosity_ptr = patch->getPatchData(cell_viscosity_id);
-      pdat::CellData<double> &cell_viscosity(*cell_viscosity_ptr);
-      if(d_dim.getValue()==2)
+  for (int ln = 0; ln <= d_hierarchy->getFinestLevelNumber(); ++ln)
+    {
+      tbox::Pointer<hier::PatchLevel> level = d_hierarchy->getPatchLevel(ln);
+      hier::PatchLevel::Iterator i_p(*level);
+      for ( ; i_p; i_p++)
         {
-          tbox::Pointer<pdat::NodeData<double> >
-            edge_viscosity_ptr = patch->getPatchData(edge_viscosity_id);
-          pdat::NodeData<double> &edge_viscosity(*edge_viscosity_ptr);
+          tbox::Pointer<hier::Patch> patch = *i_p;
+          tbox::Pointer<pdat::CellData<double> >
+            cell_viscosity_ptr = patch->getPatchData(cell_viscosity_id);
+          pdat::CellData<double> &cell_viscosity(*cell_viscosity_ptr);
+          if(d_dim.getValue()==2)
+            {
+              tbox::Pointer<pdat::NodeData<double> >
+                edge_viscosity_ptr = patch->getPatchData(edge_viscosity_id);
+              pdat::NodeData<double> &edge_viscosity(*edge_viscosity_ptr);
 
-          for(pdat::NodeIterator ni(edge_viscosity.getBox()); ni; ni++)
+              for(pdat::NodeIterator ni(edge_viscosity.getBox()); ni; ni++)
+                {
+                  pdat::NodeIndex e=ni();
+                  pdat::CellIndex c(e);
+                  cell_viscosity(c);
+                  cell_viscosity(c-ip);
+                  cell_viscosity(c-jp);
+                  cell_viscosity(c-ip-jp);
+                  edge_viscosity(e)=
+                    pow(cell_viscosity(c)*cell_viscosity(c-ip)
+                        *cell_viscosity(c-jp)*cell_viscosity(c-ip-jp),0.25);
+                }
+            }
+          else
             {
-              pdat::NodeIndex e=ni();
-              pdat::CellIndex c(e);
-              cell_viscosity(c);
-              cell_viscosity(c-ip);
-              cell_viscosity(c-jp);
-              cell_viscosity(c-ip-jp);
-              edge_viscosity(e)=
-                pow(cell_viscosity(c)*cell_viscosity(c-ip)
-                    *cell_viscosity(c-jp)*cell_viscosity(c-ip-jp),0.25);
-            }
-        }
-      else
-        {
-          tbox::Pointer<pdat::EdgeData<double> >
-            edge_viscosity_ptr = patch->getPatchData(edge_viscosity_id);
-          pdat::EdgeData<double> &edge_viscosity(*edge_viscosity_ptr);
-          for(int axis=0;axis<3;++axis)
-            {
-              for(pdat::EdgeIterator ni(edge_viscosity.getBox(),axis); ni; ni++)
+              tbox::Pointer<pdat::EdgeData<double> >
+                edge_viscosity_ptr = patch->getPatchData(edge_viscosity_id);
+              pdat::EdgeData<double> &edge_viscosity(*edge_viscosity_ptr);
+              for(int axis=0;axis<3;++axis)
                 {
-                  pdat::EdgeIndex e=ni();
-                  pdat::CellIndex c(e);
-                  edge_viscosity(e)=
-                    pow(cell_viscosity(c)
-                        *cell_viscosity(c-ip)
-                        *cell_viscosity(c-jp)
-                        *cell_viscosity(c-ip-jp)
-                        *cell_viscosity(c-kp)
-                        *cell_viscosity(c-ip-kp)
-                        *cell_viscosity(c-jp-kp)
-                        *cell_viscosity(c-ip-jp-kp),0.125);
+                  for(pdat::EdgeIterator ni(edge_viscosity.getBox(),axis); ni; ni++)
+                    {
+                      pdat::EdgeIndex e=ni();
+                      pdat::CellIndex c(e);
+                      edge_viscosity(e)=
+                        pow(cell_viscosity(c)
+                            *cell_viscosity(c-ip)
+                            *cell_viscosity(c-jp)
+                            *cell_viscosity(c-ip-jp)
+                            *cell_viscosity(c-kp)
+                            *cell_viscosity(c-ip-kp)
+                            *cell_viscosity(c-jp-kp)
+                            *cell_viscosity(c-ip-jp-kp),0.125);
+                    }
                 }
             }
         }
+
+      /* Ghost fill */
+      xfer::RefineAlgorithm refiner(d_dim);
+      refiner.registerRefine(edge_viscosity_id,edge_viscosity_id,
+                             edge_viscosity_id,
+                             tbox::Pointer<xfer::RefineOperator>(0));
+
+      tbox::Pointer<xfer::RefineSchedule> schedule=
+        refiner.createSchedule(d_hierarchy->getPatchLevel(ln));
+        
+      schedule->fillData(0.0,false);
     }
-  }
-
-
 }