[cig-commits] commit: Add in 3D coarsening for velocity
Mercurial
hg at geodynamics.org
Fri Apr 15 15:56:41 PDT 2011
changeset: 159:cf34865e7575
tag: tip
user: Walter Landry <wlandry at caltech.edu>
date: Fri Apr 15 15:55:38 2011 -0700
files: V_Coarsen.C V_Coarsen.h V_Coarsen/coarsen_2D.C V_Coarsen/coarsen_3D.C wscript
description:
Add in 3D coarsening for velocity
diff -r 9a8edb0b18c9 -r cf34865e7575 V_Coarsen.C
--- a/V_Coarsen.C Fri Apr 15 15:55:14 2011 -0700
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
@@ -1,180 +0,0 @@
-/*************************************************************************
- *
- * This file is part of the SAMRAI distribution. For full copyright
- * information, see COPYRIGHT and COPYING.LESSER.
- *
- * Copyright: (c) 1997-2010 Lawrence Livermore National Security, LLC
- * Description: Weighted averaging operator for side-centered double data on
- * a Cartesian mesh.
- *
- ************************************************************************/
-
-#ifndef included_geom_V_Coarsen_C
-#define included_geom_V_Coarsen_C
-
-#include "V_Coarsen.h"
-
-#include <float.h>
-#include <math.h>
-#include "SAMRAI/geom/CartesianPatchGeometry.h"
-#include "SAMRAI/hier/Index.h"
-#include "SAMRAI/pdat/SideData.h"
-#include "SAMRAI/pdat/SideVariable.h"
-#include "SAMRAI/tbox/Utilities.h"
-#include "Boundary.h"
-
-using namespace SAMRAI;
-
-inline void coarsen_v(const pdat::SideIndex &coarse,
- const hier::Index &ip, const hier::Index &jp,
- tbox::Pointer<pdat::SideData<double> > &v,
- tbox::Pointer<pdat::SideData<double> > &v_fine )
-{
- pdat::SideIndex center(coarse*2);
- (*v)(coarse)=((*v_fine)(center) + (*v_fine)(center+jp))/4
- + ((*v_fine)(center-ip) + (*v_fine)(center-ip+jp)
- + (*v_fine)(center+ip) + (*v_fine)(center+jp+ip))/8;
-}
-
-
-void SAMRAI::geom::V_Coarsen::coarsen(hier::Patch& coarse,
- const hier::Patch& fine,
- const int dst_component,
- const int src_component,
- const hier::Box& coarse_box,
- const hier::IntVector& ratio) const
-{
- const tbox::Dimension& dim(getDim());
-
- TBOX_DIM_ASSERT_CHECK_DIM_ARGS4(dim, coarse, fine, coarse_box, ratio);
-
- tbox::Pointer<pdat::SideData<double> >
- v_fine = fine.getPatchData(src_component);
- tbox::Pointer<pdat::SideData<double> >
- v = coarse.getPatchData(dst_component);
-
- TBOX_ASSERT(!v.isNull());
- TBOX_ASSERT(!v_fine.isNull());
- TBOX_ASSERT(v_fine->getDepth() == v->getDepth());
- TBOX_ASSERT(v->getDepth() == 1);
-
- const hier::IntVector& directions(v->getDirectionVector());
-
- TBOX_ASSERT(directions ==
- hier::IntVector::min(directions, v_fine->getDirectionVector()));
-
- const tbox::Pointer<CartesianPatchGeometry> cgeom =
- coarse.getPatchGeometry();
-
- /* Numbering of v nodes is
-
- x--x--x--x--x Fine
- 0 1 2 3 4
-
- x-----x-----x Coarse
- 0 1 2
-
- So the i'th coarse point is affected by the i*2-1,
- i*2, and i*2+1 fine points. So, for example, i_fine=3
- affects both i_coarse=1 and i_coarse=2.
-
-
- |---------------------------------------|
- | | | | |
- f f f f f f f f f
- | | | | |
- c---------c---------c---------c---------c
- | | | | |
- f f f f f f f f f
- | | | | |
- |---------------------------------------|
- | | | | |
- f f f f f f f f f
- | | | | |
- c---------c---------c---------c---------c
- | | | | |
- f f f f f f f f f
- | | | | |
- |---------------------------------------|
-
- In 2D, a coarse point depends on six points. In this
- case, (i*2,j*2), (i*2,j*2+1), (i*2-1,j*2),
- (i*2-1,j*2+1), (i*2+1,j*2), (i*2+1,j*2+1).
-
- The coarse/fine boundaries get fixed up later in
- V_Coarsen_Patch_Strategy::postprocessCoarsen.
- */
- hier::Index ip(1,0), jp(0,1);
- for(int j=coarse_box.lower(1); j<=coarse_box.upper(1)+1; ++j)
- for(int i=coarse_box.lower(0); i<=coarse_box.upper(0)+1; ++i)
- {
- // tbox::plog << "V Coarsen "
- // << coarse.getPatchLevelNumber() << " "
- // << i << " "
- // << j << " "
- // << fine.getBox().lower(0) << " "
- // << fine.getBox().upper(0) << " "
- // << fine.getBox().lower(1) << " "
- // << fine.getBox().upper(1) << " "
- // << v_fine->getGhostBox().lower(0) << " "
- // << v_fine->getGhostBox().upper(0) << " "
- // << v_fine->getGhostBox().lower(1) << " "
- // << v_fine->getGhostBox().upper(1) << " ";
-
- if(directions(0) && j!=coarse_box.upper(1)+1)
- {
- pdat::SideIndex coarse(hier::Index(i,j),0,
- pdat::SideIndex::Lower);
- pdat::SideIndex center(coarse*2);
- if((i==coarse_box.lower(0)
- && cgeom->getTouchesRegularBoundary(0,0))
- || (i==coarse_box.upper(0)+1
- && cgeom->getTouchesRegularBoundary(0,1)))
- {
- (*v)(coarse)=((*v_fine)(center) + (*v_fine)(center+jp))/2;
- }
- else
- {
- coarsen_v(coarse,ip,jp,v,v_fine);
- }
- // tbox::plog << "vx "
- // << (*v)(coarse) << " "
- // << (*v_fine)(center) << " "
- // << (*v_fine)(center+jp) << " "
- // << (*v_fine)(center+ip) << " "
- // << (*v_fine)(center+jp+ip) << " "
- // << (*v_fine)(center-ip) << " "
- // << (*v_fine)(center+jp-ip) << " ";
- }
- if(directions(1) && i!=coarse_box.upper(0)+1)
- {
- pdat::SideIndex coarse(hier::Index(i,j),1,
- pdat::SideIndex::Lower);
- pdat::SideIndex center(coarse*2);
- if((j==coarse_box.lower(1)
- && cgeom->getTouchesRegularBoundary(1,0))
- || (j==coarse_box.upper(1)+1
- && cgeom->getTouchesRegularBoundary(1,1)))
- {
- // tbox::plog << "vy bound ";
- (*v)(coarse)=((*v_fine)(center) + (*v_fine)(center+ip))/2;
- }
- else
- {
- // tbox::plog << "vy center ";
- coarsen_v(coarse,jp,ip,v,v_fine);
- }
- // tbox::plog << "vy "
- // << (*v)(coarse) << " "
- // << (*v_fine)(center) << " "
- // << (*v_fine)(center+ip) << " "
- // << (*v_fine)(center+jp) << " "
- // << (*v_fine)(center+jp+ip) << " "
- // << (*v_fine)(center-jp) << " "
- // << (*v_fine)(center+ip-jp) << " ";
- }
- // tbox::plog << "\n";
- }
-}
-
-#endif
diff -r 9a8edb0b18c9 -r cf34865e7575 V_Coarsen.h
--- a/V_Coarsen.h Fri Apr 15 15:55:14 2011 -0700
+++ b/V_Coarsen.h Fri Apr 15 15:55:38 2011 -0700
@@ -112,13 +112,34 @@ public:
* coarsening operator.
*/
void
- coarsen(
- hier::Patch& coarse,
- const hier::Patch& fine,
- const int dst_component,
- const int src_component,
- const hier::Box& coarse_box,
- const hier::IntVector& ratio) const;
+ coarsen(hier::Patch& coarse,
+ const hier::Patch& fine,
+ const int dst_component,
+ const int src_component,
+ const hier::Box& coarse_box,
+ const hier::IntVector& ratio) const
+ {
+ if(getDim().getValue()==2)
+ coarsen_2D(coarse,fine,dst_component,src_component,coarse_box,ratio);
+ else
+ coarsen_3D(coarse,fine,dst_component,src_component,coarse_box,ratio);
+ }
+
+ void
+ coarsen_2D(hier::Patch& coarse,
+ const hier::Patch& fine,
+ const int dst_component,
+ const int src_component,
+ const hier::Box& coarse_box,
+ const hier::IntVector& ratio) const;
+
+ void
+ coarsen_3D(hier::Patch& coarse,
+ const hier::Patch& fine,
+ const int dst_component,
+ const int src_component,
+ const hier::Box& coarse_box,
+ const hier::IntVector& ratio) const;
private:
std::string d_name_id;
diff -r 9a8edb0b18c9 -r cf34865e7575 V_Coarsen/coarsen_2D.C
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/V_Coarsen/coarsen_2D.C Fri Apr 15 15:55:38 2011 -0700
@@ -0,0 +1,147 @@
+/*************************************************************************
+ *
+ * This file is part of the SAMRAI distribution. For full copyright
+ * information, see COPYRIGHT and COPYING.LESSER.
+ *
+ * Copyright: (c) 1997-2010 Lawrence Livermore National Security, LLC
+ * Description: Weighted averaging operator for side-centered double data on
+ * a Cartesian mesh.
+ *
+ ************************************************************************/
+
+#ifndef included_geom_V_Coarsen_C
+#define included_geom_V_Coarsen_C
+
+#include "V_Coarsen.h"
+
+#include <float.h>
+#include <math.h>
+#include "SAMRAI/geom/CartesianPatchGeometry.h"
+#include "SAMRAI/hier/Index.h"
+#include "SAMRAI/pdat/SideData.h"
+#include "SAMRAI/pdat/SideVariable.h"
+#include "SAMRAI/tbox/Utilities.h"
+#include "Boundary.h"
+
+using namespace SAMRAI;
+
+inline void coarsen_point_2D(const pdat::SideIndex &coarse,
+ const hier::Index &ip, const hier::Index &jp,
+ tbox::Pointer<pdat::SideData<double> > &v,
+ tbox::Pointer<pdat::SideData<double> > &v_fine )
+{
+ pdat::SideIndex center(coarse*2);
+ (*v)(coarse)=((*v_fine)(center) + (*v_fine)(center+jp))/4
+ + ((*v_fine)(center-ip) + (*v_fine)(center-ip+jp)
+ + (*v_fine)(center+ip) + (*v_fine)(center+jp+ip))/8;
+}
+
+
+void SAMRAI::geom::V_Coarsen::coarsen_2D(hier::Patch& coarse,
+ const hier::Patch& fine,
+ const int dst_component,
+ const int src_component,
+ const hier::Box& coarse_box,
+ const hier::IntVector& ratio) const
+{
+ const tbox::Dimension& dim(getDim());
+
+ TBOX_DIM_ASSERT_CHECK_DIM_ARGS4(dim, coarse, fine, coarse_box, ratio);
+
+ tbox::Pointer<pdat::SideData<double> >
+ v_fine = fine.getPatchData(src_component);
+ tbox::Pointer<pdat::SideData<double> >
+ v = coarse.getPatchData(dst_component);
+
+ TBOX_ASSERT(!v.isNull());
+ TBOX_ASSERT(!v_fine.isNull());
+ TBOX_ASSERT(v_fine->getDepth() == v->getDepth());
+ TBOX_ASSERT(v->getDepth() == 1);
+
+ const hier::IntVector& directions(v->getDirectionVector());
+
+ TBOX_ASSERT(directions ==
+ hier::IntVector::min(directions, v_fine->getDirectionVector()));
+
+ const tbox::Pointer<CartesianPatchGeometry> cgeom =
+ coarse.getPatchGeometry();
+
+ /* Numbering of v nodes is
+
+ x--x--x--x--x Fine
+ 0 1 2 3 4
+
+ x-----x-----x Coarse
+ 0 1 2
+
+ So the i'th coarse point is affected by the i*2-1,
+ i*2, and i*2+1 fine points. So, for example, i_fine=3
+ affects both i_coarse=1 and i_coarse=2.
+
+ |---------------------------------------------------------------|
+ | | | | |
+ f f f f f f f f f
+ | | | | |
+ c c c c c
+ | | | | |
+ f f f f f f f f f
+ | | | | |
+ |---------------------------------------------------------------|
+ | | | | |
+ f f f f f f f f f
+ | | | | |
+ c c c c c
+ | | | | |
+ f f f f f f f f f
+ | | | | |
+ |---------------------------------------------------------------|
+
+ In 2D, a coarse point depends on six points. In this
+ case, (i*2,j*2), (i*2,j*2+1), (i*2-1,j*2),
+ (i*2-1,j*2+1), (i*2+1,j*2), (i*2+1,j*2+1).
+
+ The coarse/fine boundaries get fixed up later in
+ V_Coarsen_Patch_Strategy::postprocessCoarsen.
+ */
+ hier::Index ip(1,0), jp(0,1);
+ for(int j=coarse_box.lower(1); j<=coarse_box.upper(1)+1; ++j)
+ for(int i=coarse_box.lower(0); i<=coarse_box.upper(0)+1; ++i)
+ {
+ if(directions(0) && j!=coarse_box.upper(1)+1)
+ {
+ pdat::SideIndex coarse(hier::Index(i,j),0,
+ pdat::SideIndex::Lower);
+ pdat::SideIndex center(coarse*2);
+ if((i==coarse_box.lower(0)
+ && cgeom->getTouchesRegularBoundary(0,0))
+ || (i==coarse_box.upper(0)+1
+ && cgeom->getTouchesRegularBoundary(0,1)))
+ {
+ (*v)(coarse)=((*v_fine)(center) + (*v_fine)(center+jp))/2;
+ }
+ else
+ {
+ coarsen_point_2D(coarse,ip,jp,v,v_fine);
+ }
+ }
+ if(directions(1) && i!=coarse_box.upper(0)+1)
+ {
+ pdat::SideIndex coarse(hier::Index(i,j),1,
+ pdat::SideIndex::Lower);
+ pdat::SideIndex center(coarse*2);
+ if((j==coarse_box.lower(1)
+ && cgeom->getTouchesRegularBoundary(1,0))
+ || (j==coarse_box.upper(1)+1
+ && cgeom->getTouchesRegularBoundary(1,1)))
+ {
+ (*v)(coarse)=((*v_fine)(center) + (*v_fine)(center+ip))/2;
+ }
+ else
+ {
+ coarsen_point_2D(coarse,jp,ip,v,v_fine);
+ }
+ }
+ }
+}
+
+#endif
diff -r 9a8edb0b18c9 -r cf34865e7575 V_Coarsen/coarsen_3D.C
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/V_Coarsen/coarsen_3D.C Fri Apr 15 15:55:38 2011 -0700
@@ -0,0 +1,199 @@
+/*************************************************************************
+ *
+ * This file is part of the SAMRAI distribution. For full copyright
+ * information, see COPYRIGHT and COPYING.LESSER.
+ *
+ * Copyright: (c) 1997-2010 Lawrence Livermore National Security, LLC
+ * Description: Weighted averaging operator for side-centered double data on
+ * a Cartesian mesh.
+ *
+ ************************************************************************/
+
+#ifndef included_geom_V_Coarsen_C
+#define included_geom_V_Coarsen_C
+
+#include "V_Coarsen.h"
+
+#include <float.h>
+#include <math.h>
+#include "SAMRAI/geom/CartesianPatchGeometry.h"
+#include "SAMRAI/hier/Index.h"
+#include "SAMRAI/pdat/SideData.h"
+#include "SAMRAI/pdat/SideVariable.h"
+#include "SAMRAI/tbox/Utilities.h"
+#include "Boundary.h"
+
+using namespace SAMRAI;
+
+inline void coarsen_point_3D(const pdat::SideIndex &coarse,
+ const hier::Index &ip, const hier::Index &jp,
+ const hier::Index &kp,
+ tbox::Pointer<pdat::SideData<double> > &v,
+ tbox::Pointer<pdat::SideData<double> > &v_fine )
+{
+ pdat::SideIndex center(coarse*2);
+ (*v)(coarse)=((*v_fine)(center) + (*v_fine)(center+jp)
+ + (*v_fine)(center+kp) + (*v_fine)(center+jp+kp))/8
+ + ((*v_fine)(center-ip) + (*v_fine)(center-ip+jp)
+ + (*v_fine)(center-ip+kp) + (*v_fine)(center-ip+jp+kp)
+ + (*v_fine)(center+ip) + (*v_fine)(center+jp+ip)
+ + (*v_fine)(center+ip+kp) + (*v_fine)(center+jp+ip+kp))/16;
+}
+
+
+void SAMRAI::geom::V_Coarsen::coarsen_3D(hier::Patch& coarse,
+ const hier::Patch& fine,
+ const int dst_component,
+ const int src_component,
+ const hier::Box& coarse_box,
+ const hier::IntVector& ratio) const
+{
+ const tbox::Dimension& dim(getDim());
+
+ TBOX_DIM_ASSERT_CHECK_DIM_ARGS4(dim, coarse, fine, coarse_box, ratio);
+
+ tbox::Pointer<pdat::SideData<double> >
+ v_fine = fine.getPatchData(src_component);
+ tbox::Pointer<pdat::SideData<double> >
+ v = coarse.getPatchData(dst_component);
+
+ TBOX_ASSERT(!v.isNull());
+ TBOX_ASSERT(!v_fine.isNull());
+ TBOX_ASSERT(v_fine->getDepth() == v->getDepth());
+ TBOX_ASSERT(v->getDepth() == 1);
+
+ const hier::IntVector& directions(v->getDirectionVector());
+
+ TBOX_ASSERT(directions ==
+ hier::IntVector::min(directions, v_fine->getDirectionVector()));
+
+ const tbox::Pointer<CartesianPatchGeometry> cgeom =
+ coarse.getPatchGeometry();
+
+ /* Numbering of v nodes is
+
+ x--x--x--x--x Fine
+ 0 1 2 3 4
+
+ x-----x-----x Coarse
+ 0 1 2
+
+ So the i'th coarse point is affected by the i*2-1,
+ i*2, and i*2+1 fine points. So, for example, i_fine=3
+ affects both i_coarse=1 and i_coarse=2.
+
+ |---------------------------------------------------------------|
+ | | | | |
+ f f f f f f f f f
+ | | | | |
+ c c c c c
+ | | | | |
+ f f f f f f f f f
+ | | | | |
+ |---------------------------------------------------------------|
+ | | | | |
+ f f f f f f f f f
+ | | | | |
+ c c c c c
+ | | | | |
+ f f f f f f f f f
+ | | | | |
+ |---------------------------------------------------------------|
+
+ In 2D, a coarse point depends on six points. In this
+ case, (i*2,j*2), (i*2,j*2+1), (i*2-1,j*2),
+ (i*2-1,j*2+1), (i*2+1,j*2), (i*2+1,j*2+1).
+
+
+ --------------------
+ / /|
+ / / |
+ / / |
+ / / |
+ ------------------- |
+ | | |
+ | f f | |
+ | | |
+ | C | /
+ | | /
+ | f f | /
+ | |/
+ -------------------
+
+ In 3D, a coarse point depend on 12 points
+ (i*2,j*2,k*2), (i*2,j*2+1,k*2), (i*2,j*2,k*2+1), (i*2,j*2+1,k*2+1),
+ (i*2+1,j*2,k*2), (i*2+1,j*2+1,k*2), (i*2+1,j*2,k*2+1), (i*2+1,j*2+1,k*2+1),
+ (i*2-1,j*2,k*2), (i*2-1,j*2+1,k*2), (i*2-1,j*2,k*2+1), (i*2-1,j*2+1,k*2+1)
+
+ The coarse/fine boundaries get fixed up later in
+ V_Coarsen_Patch_Strategy::postprocessCoarsen.
+ */
+ hier::Index ip(1,0,0), jp(0,1,0), kp(0,0,1);
+ for(int k=coarse_box.lower(2); k<=coarse_box.upper(2)+1; ++k)
+ for(int j=coarse_box.lower(1); j<=coarse_box.upper(1)+1; ++j)
+ for(int i=coarse_box.lower(0); i<=coarse_box.upper(0)+1; ++i)
+ {
+ if(directions(0) && j!=coarse_box.upper(1)+1
+ && k!=coarse_box.upper(2)+1)
+ {
+ pdat::SideIndex coarse(hier::Index(i,j,k),0,
+ pdat::SideIndex::Lower);
+ pdat::SideIndex center(coarse*2);
+ if((i==coarse_box.lower(0)
+ && cgeom->getTouchesRegularBoundary(0,0))
+ || (i==coarse_box.upper(0)+1
+ && cgeom->getTouchesRegularBoundary(0,1)))
+ {
+ (*v)(coarse)=
+ ((*v_fine)(center) + (*v_fine)(center+jp)
+ + (*v_fine)(center+kp) + (*v_fine)(center+jp+kp))/4;
+ }
+ else
+ {
+ coarsen_point_3D(coarse,ip,jp,kp,v,v_fine);
+ }
+ }
+ if(directions(1) && i!=coarse_box.upper(0)+1
+ && k!=coarse_box.upper(2)+1)
+ {
+ pdat::SideIndex coarse(hier::Index(i,j,k),1,
+ pdat::SideIndex::Lower);
+ pdat::SideIndex center(coarse*2);
+ if((j==coarse_box.lower(1)
+ && cgeom->getTouchesRegularBoundary(1,0))
+ || (j==coarse_box.upper(1)+1
+ && cgeom->getTouchesRegularBoundary(1,1)))
+ {
+ (*v)(coarse)=
+ ((*v_fine)(center) + (*v_fine)(center+ip)
+ + (*v_fine)(center+kp) + (*v_fine)(center+ip+kp))/4;
+ }
+ else
+ {
+ coarsen_point_3D(coarse,jp,kp,ip,v,v_fine);
+ }
+ }
+ if(directions(2) && i!=coarse_box.upper(0)+1
+ && j!=coarse_box.upper(1)+1)
+ {
+ pdat::SideIndex coarse(hier::Index(i,j,k),2,
+ pdat::SideIndex::Lower);
+ pdat::SideIndex center(coarse*2);
+ if((k==coarse_box.lower(2)
+ && cgeom->getTouchesRegularBoundary(2,0))
+ || (k==coarse_box.upper(2)+1
+ && cgeom->getTouchesRegularBoundary(2,1)))
+ {
+ (*v)(coarse)=
+ ((*v_fine)(center) + (*v_fine)(center+ip)
+ + (*v_fine)(center+jp) + (*v_fine)(center+ip+jp))/4;
+ }
+ else
+ {
+ coarsen_point_3D(coarse,kp,ip,jp,v,v_fine);
+ }
+ }
+ }
+}
+
+#endif
diff -r 9a8edb0b18c9 -r cf34865e7575 wscript
--- a/wscript Fri Apr 15 15:55:14 2011 -0700
+++ b/wscript Fri Apr 15 15:55:38 2011 -0700
@@ -22,7 +22,8 @@ def build(bld):
'P_Refine.C',
'P_MDPI_Refine.C',
'V_Refine.C',
- 'V_Coarsen.C',
+ 'V_Coarsen/coarsen_2D.C',
+ 'V_Coarsen/coarsen_3D.C',
'P_Boundary_Refine.C',
'V_Boundary_Refine/refine.C',
'V_Boundary_Refine/Update_V.C',
More information about the CIG-COMMITS
mailing list