[cig-commits] commit: Generalize set_boundary to handle stress BC. It is hard-coded for the
Mercurial
hg at geodynamics.org
Mon Apr 11 12:14:38 PDT 2011
changeset: 147:a43af5d15c5a
user: Walter Landry <wlandry at caltech.edu>
date: Mon Apr 11 01:10:32 2011 -0700
files: FACStokes/solveStokes.C StokesFACOps.h StokesFACOps/computeCompositeResidualOnLevel.C StokesFACOps/prolongErrorAndCorrect.C StokesFACOps/restrictSolution.C StokesFACOps/set_boundaries.C StokesFACOps/smooth_Gerya.C StokesFACOps/smooth_Tackley_2D.C StokesFACOps/smooth_Tackley_3D.C StokesFACOps/xeqScheduleGhostFill.C StokesFACSolver.h V_Coarsen_Patch_Strategy.h V_Refine_Patch_Strategy.C set_V_boundary.C set_V_boundary.h set_boundary.C set_boundary.h wscript
description:
Generalize set_boundary to handle stress BC. It is hard-coded for the
Geomod extension benchmark right now.
diff -r 4fa91a89c29d -r a43af5d15c5a FACStokes/solveStokes.C
--- a/FACStokes/solveStokes.C Thu Mar 24 02:42:29 2011 -0700
+++ b/FACStokes/solveStokes.C Mon Apr 11 01:10:32 2011 -0700
@@ -52,7 +52,7 @@ int SAMRAI::FACStokes::solveStokes()
v = patch->getPatchData(v_id);
v->fill(0.0);
}
- d_stokes_fac_solver.set_boundaries(v_id,level,false);
+ d_stokes_fac_solver.set_boundaries(p_id,v_id,level,false);
}
fix_viscosity();
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps.h
--- a/StokesFACOps.h Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps.h Mon Apr 11 01:10:32 2011 -0700
@@ -513,20 +513,23 @@ public:
const SAMRAIVectorReal<double>& current_soln,
const SAMRAIVectorReal<double>& residual);
- void set_boundaries(const int &v_id, const int &l)
+ void set_boundaries(const int &p_id, const int &v_id, const int &l)
{
- set_boundaries(v_id,l,true);
+ set_boundaries(p_id,v_id,l,true);
}
- void set_boundaries(const int &v_id, const int &l, const bool &rhs)
+ void set_boundaries(const int &p_id, const int &v_id, const int &l,
+ const bool &rhs)
{
tbox::Pointer<hier::PatchLevel> level = d_hierarchy->getPatchLevel(l);
- set_boundaries(v_id,level,rhs);
+ set_boundaries(p_id,v_id,level,rhs);
}
- void set_boundaries(const int &v_id, tbox::Pointer<hier::PatchLevel> &level)
+ void set_boundaries(const int &p_id, const int &v_id,
+ tbox::Pointer<hier::PatchLevel> &level)
{
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
}
- void set_boundaries(const int &v_id, tbox::Pointer<hier::PatchLevel> &level,
+ void set_boundaries(const int &p_id, const int &v_id,
+ tbox::Pointer<hier::PatchLevel> &level,
const bool &rhs);
//@}
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps/computeCompositeResidualOnLevel.C
--- a/StokesFACOps/computeCompositeResidualOnLevel.C Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps/computeCompositeResidualOnLevel.C Mon Apr 11 01:10:32 2011 -0700
@@ -53,7 +53,7 @@ void SAMRAI::solv::StokesFACOps::compute
/* S1. Fill solution ghost data. */
- set_boundaries(v_id,ln,error_equation_indicator);
+ set_boundaries(p_id,v_id,ln,error_equation_indicator);
if (ln > d_ln_min) {
/* Fill from current, next coarser level and physical boundary */
xeqScheduleGhostFill(p_id, v_id, ln);
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps/prolongErrorAndCorrect.C
--- a/StokesFACOps/prolongErrorAndCorrect.C Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps/prolongErrorAndCorrect.C Mon Apr 11 01:10:32 2011 -0700
@@ -98,7 +98,8 @@ void SAMRAI::solv::StokesFACOps::prolong
d_side_scratch_id,
dest_ln);
- set_boundaries(s.getComponentDescriptorIndex(1),fine_level,true);
+ set_boundaries(s.getComponentDescriptorIndex(0),
+ s.getComponentDescriptorIndex(1),fine_level,true);
/*
* Add the refined error in the scratch space to the error currently
* residing in the destination level.
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps/restrictSolution.C
--- a/StokesFACOps/restrictSolution.C Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps/restrictSolution.C Mon Apr 11 01:10:32 2011 -0700
@@ -67,7 +67,7 @@ void SAMRAI::solv::StokesFACOps::restric
xeqScheduleURestriction(p_dst,p_src,v_dst,v_src,dest_ln);
tbox::Pointer<hier::PatchLevel> level = d_hierarchy->getPatchLevel(dest_ln);
- set_boundaries(v_dst,level,false);
+ set_boundaries(p_dst,v_dst,level,false);
// v_refine_patch_strategy.setHomogeneousBc(false);
p_refine_patch_strategy.setTargetDataId(d.getComponentDescriptorIndex(0));
v_refine_patch_strategy.setTargetDataId(d.getComponentDescriptorIndex(1));
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps/set_boundaries.C
--- a/StokesFACOps/set_boundaries.C Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps/set_boundaries.C Mon Apr 11 01:10:32 2011 -0700
@@ -41,16 +41,17 @@
#include "SAMRAI/xfer/PatchLevelFullFillPattern.h"
#include "Boundary.h"
-#include "set_V_boundary.h"
+#include "set_boundary.h"
/* Set the physical boundaries for the velocity. */
void SAMRAI::solv::StokesFACOps::set_boundaries
-(const int &v_id, tbox::Pointer<hier::PatchLevel> &level, const bool &rhs)
+(const int &p_id, const int &v_id,
+ tbox::Pointer<hier::PatchLevel> &level, const bool &rhs)
{
for (hier::PatchLevel::Iterator pi(*level); pi; pi++)
{
tbox::Pointer<hier::Patch> patch = *pi;
- set_V_boundary(*patch,v_id,rhs);
+ set_boundary(*patch,p_id,v_id,rhs);
}
}
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps/smooth_Gerya.C
--- a/StokesFACOps/smooth_Gerya.C Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps/smooth_Gerya.C Mon Apr 11 01:10:32 2011 -0700
@@ -39,7 +39,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
correct. */
p_refine_patch_strategy.setTargetDataId(p_id);
v_refine_patch_strategy.setTargetDataId(v_id);
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
xeqScheduleGhostFillNoCoarse(p_rhs_id,v_rhs_id,ln);
if (ln > d_ln_min) {
@@ -159,7 +159,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
}
}
}
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
}
// if (residual_tolerance >= 0.0) {
/*
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps/smooth_Tackley_2D.C
--- a/StokesFACOps/smooth_Tackley_2D.C Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps/smooth_Tackley_2D.C Mon Apr 11 01:10:32 2011 -0700
@@ -38,7 +38,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
correct. */
p_refine_patch_strategy.setTargetDataId(p_id);
v_refine_patch_strategy.setTargetDataId(v_id);
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
xeqScheduleGhostFillNoCoarse(p_rhs_id,v_rhs_id,ln);
if (ln > d_ln_min) {
@@ -121,7 +121,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
}
}
}
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
}
@@ -175,7 +175,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
}
}
}
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
}
@@ -238,6 +238,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
p(center)+=dp(center);
}
}
+ set_boundaries(p_id,v_id,level,true);
/* fix v sweep */
@@ -299,7 +300,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
}
}
}
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
// if (residual_tolerance >= 0.0) {
/*
@@ -316,10 +317,10 @@ void SAMRAI::solv::StokesFACOps::smooth_
mpi.AllReduce(&tmp, 1, MPI_MIN);
}
converged=(tmp==1);
- if (d_enable_logging)
- tbox::plog
- // << d_object_name << "\n"
- << "Tackley " << ln << " " << sweep << " : " << maxres << "\n";
+ // if (d_enable_logging)
+ // tbox::plog
+ // // << d_object_name << "\n"
+ // << "Tackley " << ln << " " << sweep << " : " << maxres << "\n";
// }
}
}
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps/smooth_Tackley_3D.C
--- a/StokesFACOps/smooth_Tackley_3D.C Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps/smooth_Tackley_3D.C Mon Apr 11 01:10:32 2011 -0700
@@ -36,7 +36,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
calculating a new pressure update requires computing in the ghost
region so that the update for the velocity inside the box will be
correct. */
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
xeqScheduleGhostFillNoCoarse(p_rhs_id,v_rhs_id,ln);
p_refine_patch_strategy.setTargetDataId(p_id);
@@ -123,7 +123,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
}
}
}
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
}
/* p sweep
@@ -256,7 +256,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
}
}
}
- set_boundaries(v_id,level,true);
+ set_boundaries(p_id,v_id,level,true);
}
// if (residual_tolerance >= 0.0) {
/*
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACOps/xeqScheduleGhostFill.C
--- a/StokesFACOps/xeqScheduleGhostFill.C Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACOps/xeqScheduleGhostFill.C Mon Apr 11 01:10:32 2011 -0700
@@ -61,7 +61,7 @@ void SAMRAI::solv::StokesFACOps::xeqSche
if (!v_ghostfill_refine_schedules[dest_ln]) {
TBOX_ERROR("Expected schedule not found.");
}
- set_boundaries(v_id,dest_ln-1);
+ set_boundaries(p_id,v_id,dest_ln-1);
xfer::RefineAlgorithm refiner(d_dim);
refiner.registerRefine(v_id,v_id,v_id,v_ghostfill_refine_operator);
refiner.resetSchedule(v_ghostfill_refine_schedules[dest_ln]);
diff -r 4fa91a89c29d -r a43af5d15c5a StokesFACSolver.h
--- a/StokesFACSolver.h Thu Mar 24 02:42:29 2011 -0700
+++ b/StokesFACSolver.h Mon Apr 11 01:10:32 2011 -0700
@@ -566,11 +566,11 @@ public:
double
getResidualNorm() const;
- void set_boundaries(const int &v_id,
+ void set_boundaries(const int &p_id, const int &v_id,
tbox::Pointer<hier::PatchLevel> &level,
const bool &homogeneous)
{
- d_fac_ops.set_boundaries(v_id,level,homogeneous);
+ d_fac_ops.set_boundaries(p_id,v_id,level,homogeneous);
}
diff -r 4fa91a89c29d -r a43af5d15c5a V_Coarsen_Patch_Strategy.h
--- a/V_Coarsen_Patch_Strategy.h Thu Mar 24 02:42:29 2011 -0700
+++ b/V_Coarsen_Patch_Strategy.h Mon Apr 11 01:10:32 2011 -0700
@@ -18,7 +18,7 @@
#include "SAMRAI/pdat/SideData.h"
#include "SAMRAI/pdat/CellIndex.h"
#include "Boundary.h"
-#include "set_V_boundary.h"
+#include "set_boundary.h"
namespace SAMRAI {
namespace solv {
@@ -69,7 +69,7 @@ public:
const hier::Box& ,
const hier::IntVector& )
{
- set_V_boundary(fine,v_id,true);
+ set_boundary(fine,-1,v_id,true);
}
virtual void
diff -r 4fa91a89c29d -r a43af5d15c5a V_Refine_Patch_Strategy.C
--- a/V_Refine_Patch_Strategy.C Thu Mar 24 02:42:29 2011 -0700
+++ b/V_Refine_Patch_Strategy.C Mon Apr 11 01:10:32 2011 -0700
@@ -1,5 +1,5 @@
#include "V_Refine_Patch_Strategy.h"
-#include "set_V_boundary.h"
+#include "set_boundary.h"
void
SAMRAI::solv::V_Refine_Patch_Strategy::preprocessRefine
@@ -8,5 +8,5 @@ SAMRAI::solv::V_Refine_Patch_Strategy::p
const hier::Box& ,
const hier::IntVector& )
{
- set_V_boundary(coarse,v_id,true);
+ set_boundary(coarse,-1,v_id,true);
}
diff -r 4fa91a89c29d -r a43af5d15c5a set_V_boundary.C
--- a/set_V_boundary.C Thu Mar 24 02:42:29 2011 -0700
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
@@ -1,73 +0,0 @@
-#include "set_V_boundary.h"
-#include "Boundary.h"
-#include "SAMRAI/tbox/Pointer.h"
-#include "SAMRAI/pdat/SideData.h"
-#include "SAMRAI/pdat/CellIndex.h"
-#include "SAMRAI/geom/CartesianPatchGeometry.h"
-
-using namespace SAMRAI;
-
-void set_V_boundary(const SAMRAI::hier::Patch& patch, const int &v_id,
- const bool &rhs)
-{
- tbox::Pointer<pdat::SideData<double> > v_ptr = patch.getPatchData(v_id);
- pdat::SideData<double> &v(*v_ptr);
-
- hier::Box pbox=patch.getBox();
- hier::Box gbox=v.getGhostBox();
-
- tbox::Pointer<geom::CartesianPatchGeometry> geom = patch.getPatchGeometry();
- const tbox::Dimension Dim(patch.getDim());
- const int dim(patch.getDim().getValue());
-
- const hier::Index zero(hier::Index::getZeroIndex(Dim));
- hier::Index pp[]={zero,zero,zero};
- for(int i=0;i<dim;++i)
- pp[i][i]=1;
- /* This should really get read from the input file. */
- double lower_boundary[]={0,0,0};
- double upper_boundary[]={-1,1,0};
-
- for(int ix=0; ix<dim; ++ix)
- {
- for(pdat::SideIterator si(gbox,ix); si; si++)
- {
- pdat::SideIndex x(*si);
-
- /* Set a sentinel value for normal components */
- if((x[ix]<pbox.lower(ix) && geom->getTouchesRegularBoundary(ix,0))
- || (x[ix]>pbox.upper(ix)+1 && geom->getTouchesRegularBoundary(ix,1)))
- {
- v(x)=boundary_value;
- }
- /* Set values for normal components */
- else if(x[ix]==pbox.lower(ix) && geom->getTouchesRegularBoundary(ix,0)
- && !rhs)
- {
- v(x)=lower_boundary[ix];
- }
- else if(x[ix]==pbox.upper(ix)+1 && geom->getTouchesRegularBoundary(ix,1)
- && !rhs)
- {
- v(x)=upper_boundary[ix];
- }
- /* Set derivatives for tangential component */
- else
- {
- for(int iy=(ix+1)%dim; iy!=ix; iy=(iy+1)%dim)
- {
- if(x[iy]<pbox.lower(iy)
- && geom->getTouchesRegularBoundary(iy,0))
- {
- v(x)=v(x+pp[iy]);
- }
- else if(x[iy]>pbox.upper(iy)
- && geom->getTouchesRegularBoundary(iy,1))
- {
- v(x)=v(x-pp[iy]);
- }
- }
- }
- }
- }
-}
diff -r 4fa91a89c29d -r a43af5d15c5a set_V_boundary.h
--- a/set_V_boundary.h Thu Mar 24 02:42:29 2011 -0700
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
@@ -1,9 +0,0 @@
-#ifndef SET_V_BOUNDARY_H
-#define SET_V_BOUNDARY_H
-
-#include "SAMRAI/hier/Patch.h"
-
-void set_V_boundary(const SAMRAI::hier::Patch& patch, const int &id,
- const bool &rhs);
-
-#endif
diff -r 4fa91a89c29d -r a43af5d15c5a set_boundary.C
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/set_boundary.C Mon Apr 11 01:10:32 2011 -0700
@@ -0,0 +1,140 @@
+#include "set_boundary.h"
+#include "Boundary.h"
+#include "SAMRAI/tbox/Pointer.h"
+#include "SAMRAI/pdat/SideData.h"
+#include "SAMRAI/pdat/CellData.h"
+#include "SAMRAI/geom/CartesianPatchGeometry.h"
+
+using namespace SAMRAI;
+
+void set_boundary(const SAMRAI::hier::Patch& patch, const int &p_id,
+ const int &v_id, const bool &rhs)
+{
+ hier::Box pbox=patch.getBox();
+
+ tbox::Pointer<geom::CartesianPatchGeometry> geom = patch.getPatchGeometry();
+ const tbox::Dimension Dim(patch.getDim());
+ const int dim(patch.getDim().getValue());
+
+ const hier::Index zero(hier::Index::getZeroIndex(Dim));
+ hier::Index pp[]={zero,zero,zero};
+ for(int i=0;i<dim;++i)
+ pp[i][i]=1;
+ /* This should really get read from the input file. */
+ double lower_boundary[]={0,0,0};
+ bool lower_dirichlet[]={true,true,true};
+
+ double upper_boundary[]={-6.94444444444e4,0,0};
+ bool upper_dirichlet[]={true,false,true};
+
+ // bool upper_dirichlet[]={true,true,true};
+ // double upper_boundary[]={-1,1,0};
+
+ if(p_id!=-1)
+ {
+ tbox::Pointer<pdat::CellData<double> > p_ptr = patch.getPatchData(p_id);
+ pdat::CellData<double> &p(*p_ptr);
+
+ hier::Box gbox=p.getGhostBox();
+
+ for(pdat::CellIterator ci(gbox); ci; ci++)
+ {
+ pdat::CellIndex center(*ci);
+ for(int ix=0; ix<dim; ++ix)
+ {
+ if(center[ix]<pbox.lower(ix)
+ && geom->getTouchesRegularBoundary(ix,0))
+ {
+ if(!lower_dirichlet[ix])
+ p(center)=p(center+pp[ix]);
+ }
+ else if(center[ix]>pbox.upper(ix)
+ && geom->getTouchesRegularBoundary(ix,1))
+ {
+ if(!upper_dirichlet[ix])
+ {
+ // p(center)=p(center-pp[ix]);
+ p(center)=0;
+ }
+ }
+ }
+ }
+ }
+
+
+ tbox::Pointer<pdat::SideData<double> > v_ptr = patch.getPatchData(v_id);
+ pdat::SideData<double> &v(*v_ptr);
+
+ hier::Box gbox=v.getGhostBox();
+ for(int ix=0; ix<dim; ++ix)
+ {
+ for(pdat::SideIterator si(gbox,ix); si; si++)
+ {
+ pdat::SideIndex x(*si);
+
+ double pos_x=geom->getXLower()[0]
+ + geom->getDx()[0]*(x[0]-pbox.lower()[0]);
+
+ /* Set a sentinel value for normal components */
+ if(x[ix]<pbox.lower(ix) && geom->getTouchesRegularBoundary(ix,0))
+ {
+ if(lower_dirichlet[ix])
+ v(x)=boundary_value;
+ else
+ v(x)=v(x+pp[ix]*2);
+ }
+ else if(x[ix]>pbox.upper(ix)+1 && geom->getTouchesRegularBoundary(ix,1))
+ {
+ if(upper_dirichlet[ix])
+ v(x)=boundary_value;
+ else
+ v(x)=v(x-pp[ix]*2);
+ }
+ /* Set values for normal components */
+ else if(x[ix]==pbox.lower(ix) && geom->getTouchesRegularBoundary(ix,0)
+ && !rhs && lower_dirichlet[ix])
+ {
+ v(x)=lower_boundary[ix];
+ }
+ else if(x[ix]==pbox.upper(ix)+1 && geom->getTouchesRegularBoundary(ix,1)
+ && !rhs && upper_dirichlet[ix])
+ {
+ v(x)=upper_boundary[ix];
+ }
+ /* Set derivatives for tangential component */
+ else
+ {
+ for(int iy=(ix+1)%dim; iy!=ix; iy=(iy+1)%dim)
+ {
+ if(x[iy]<pbox.lower(iy)
+ && geom->getTouchesRegularBoundary(iy,0))
+ {
+ v(x)=v(x+pp[iy]);
+
+ if(ix==0 && iy==1)
+ {
+ if(pos_x<0.1 || rhs)
+ {
+ v(x)=-v(x+pp[iy]);
+ }
+ else
+ {
+ v(x)=-v(x+pp[iy]) + 2*upper_boundary[0];
+ }
+ }
+ }
+ else if(x[iy]>pbox.upper(iy)
+ && geom->getTouchesRegularBoundary(iy,1))
+ {
+ v(x)=v(x-pp[iy]);
+
+ if(ix==1 && iy==0)
+ {
+ v(x)=-v(x-pp[iy]);
+ }
+ }
+ }
+ }
+ }
+ }
+}
diff -r 4fa91a89c29d -r a43af5d15c5a set_boundary.h
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/set_boundary.h Mon Apr 11 01:10:32 2011 -0700
@@ -0,0 +1,9 @@
+#ifndef SET_V_BOUNDARY_H
+#define SET_V_BOUNDARY_H
+
+#include "SAMRAI/hier/Patch.h"
+
+void set_boundary(const SAMRAI::hier::Patch& patch,
+ const int &p_id, const int &v_id, const bool &rhs);
+
+#endif
diff -r 4fa91a89c29d -r a43af5d15c5a wscript
--- a/wscript Thu Mar 24 02:42:29 2011 -0700
+++ b/wscript Mon Apr 11 01:10:32 2011 -0700
@@ -30,7 +30,7 @@ def build(bld):
'V_Coarsen_Patch_Strategy.C',
'Cell_Viscosity_Coarsen.C',
'Edge_Viscosity_Coarsen.C',
- 'set_V_boundary.C',
+ 'set_boundary.C',
'StokesFACOps/StokesFACOps.C',
'StokesFACOps/checkInputPatchDataIndices.C',
'StokesFACOps/computeCompositeResidualOnLevel.C',
More information about the CIG-COMMITS
mailing list