[cig-commits] commit: Implement box benchmark

[Mercurial]

changeset:   3:82190c9edc01
tag:         tip
user:        Walter Landry <wlandry at caltech.edu>
date:        Sun Feb 20 23:33:15 2011 -0800
files:       madds-1.cc
description:
Implement box benchmark


diff -r bcf3f0757289 -r 82190c9edc01 madds-1.cc
--- a/madds-1.cc	Sun Feb 20 23:02:33 2011 -0800
+++ b/madds-1.cc	Sun Feb 20 23:33:15 2011 -0800
@@ -280,6 +280,19 @@ vmult (BlockVector<double>       &dst,
 				 // those. The definition of the depth
 				 // variable $z$ is discussed in the
 				 // introduction.
+
+template <int dim>
+bool inside_box(const Point<dim> &p)
+{
+  for (unsigned int d=0; d<dim; ++d)
+    if ((p[d] < 1./3) || (p[d] > 2./3))
+      {
+        return false;
+        break;
+      }
+  return true;
+}
+
 template <int dim>
 class TopBoundaryValues : public Function<dim>
 {
@@ -306,8 +319,8 @@ TopBoundaryValues<dim>::value (const Poi
   Assert (component < this->n_components,
           ExcIndexRange (component, 0, this->n_components));
   
-  if (component == 1 && p[1]>0.999999)
-    return 1;
+  // if (component == 1 && p[1]>0.999999)
+  //   return 1;
   return 0;
 }
 
@@ -344,8 +357,8 @@ RightBoundaryValues<dim>::value (const P
   Assert (component < this->n_components,
           ExcIndexRange (component, 0, this->n_components));
 
-  if (component == 0 && p[0]>0.999999)
-    return -1;
+  // if (component == 0 && p[0]>0.999999)
+  //   return -1;
   return 0;
 }
 
@@ -382,8 +395,13 @@ RightHandSide<dim>::value (const Point<d
 RightHandSide<dim>::value (const Point<dim>  &p,
                            const unsigned int component) const
 {
-  // if(component==1)
-  //   return 10;
+  if(component==1)
+    {
+      if(inside_box(p))
+        return 1.03;
+      else
+        return 1;
+    }
   return 0;
 }
 
@@ -542,24 +560,10 @@ void StokesProblem<dim>::setup_dofs ()
   system_rhs.collect_sizes ();
 }
 
-
-
 template <int dim>
 double StokesProblem<dim>::get_viscosity (const Point<dim> &p) const
 {
-  bool inside_high_viscosity_region = true;
-  if(p[0]*p[0] + p[1]*p[1] < 0.2*0.2)
-    inside_high_viscosity_region = false;
-    
-  // for (unsigned int d=0; d<dim; ++d)
-  //   if ((p[d] < 1./3)
-  //       ||
-  //       (p[d] > 2./3))
-  //     {
-  //       inside_high_viscosity_region = false;
-  //       break;
-  //     }
-  if (inside_high_viscosity_region)
+  if (inside_box(p))
     return parameters.eta_jump * parameters.eta;
   else
     return parameters.eta;