[cig-commits] commit: Make smoother_Tackley work in parallel.
Mercurial
hg at geodynamics.org
Sat Mar 12 05:37:19 PST 2011
changeset: 118:7d4eceb1a9e1
tag: tip
user: Walter Landry <wlandry at caltech.edu>
date: Fri Mar 11 20:37:48 2011 -0800
files: FACStokes/initializeLevelData.C StokesFACOps/smooth_Tackley.C input/shear_corner.input
description:
Make smoother_Tackley work in parallel.
diff -r c6b0d0a8b7f5 -r 7d4eceb1a9e1 FACStokes/initializeLevelData.C
--- a/FACStokes/initializeLevelData.C Fri Mar 11 16:56:26 2011 -0800
+++ b/FACStokes/initializeLevelData.C Fri Mar 11 20:37:48 2011 -0800
@@ -63,7 +63,7 @@ void SAMRAI::FACStokes::initializeLevelD
}
const double inclusion_radius=0.5;
- const double inclusion_viscosity=10;
+ const double inclusion_viscosity=1e2;
const double background_viscosity=1;
/*
@@ -126,6 +126,7 @@ void SAMRAI::FACStokes::initializeLevelD
double y=geom->getXLower()[1]
+ geom->getDx()[1]*(e[1]-visc_box.lower()[1]);
if(x*x + y*y < inclusion_radius*inclusion_radius)
+ // if(x<inclusion_radius && y<inclusion_radius)
edge_viscosity(e)=inclusion_viscosity;
else
edge_viscosity(e)=background_viscosity;
diff -r c6b0d0a8b7f5 -r 7d4eceb1a9e1 StokesFACOps/smooth_Tackley.C
--- a/StokesFACOps/smooth_Tackley.C Fri Mar 11 16:56:26 2011 -0800
+++ b/StokesFACOps/smooth_Tackley.C Fri Mar 11 20:37:48 2011 -0800
@@ -112,10 +112,11 @@ void SAMRAI::solv::StokesFACOps::smooth_
maxres=0;
/* vx sweep */
+
+ xeqScheduleGhostFillNoCoarse(p_id,invalid_id,ln);
for(int rb=0;rb<2;++rb)
{
- // Need to sync
- xeqScheduleGhostFillNoCoarse(p_id,v_id,ln);
+ xeqScheduleGhostFillNoCoarse(invalid_id,v_id,ln);
for (hier::PatchLevel::Iterator pi(*level); pi; pi++)
{
tbox::Pointer<hier::Patch> patch = *pi;
@@ -225,10 +226,11 @@ void SAMRAI::solv::StokesFACOps::smooth_
/* vy sweep */
+
+ xeqScheduleGhostFillNoCoarse(p_id,invalid_id,ln);
for(int rb=0;rb<2;++rb)
{
- // Need to sync
- xeqScheduleGhostFillNoCoarse(p_id,v_id,ln);
+ xeqScheduleGhostFillNoCoarse(invalid_id,v_id,ln);
for (hier::PatchLevel::Iterator pi(*level); pi; pi++)
{
tbox::Pointer<hier::Patch> patch = *pi;
@@ -338,11 +340,132 @@ void SAMRAI::solv::StokesFACOps::smooth_
- /* p sweep */
+ /* dp sweep */
+ xeqScheduleGhostFillNoCoarse(-1,v_id,ln);
for(int rb=0;rb<2;++rb)
{
- // Need to sync
- xeqScheduleGhostFillNoCoarse(-1,v_id,ln);
+ for (hier::PatchLevel::Iterator pi(*level); pi; pi++)
+ {
+ tbox::Pointer<hier::Patch> patch = *pi;
+
+ tbox::Pointer<pdat::CellData<double> > p_ptr =
+ patch->getPatchData(p_id);
+ pdat::CellData<double> &p(*p_ptr);
+ tbox::Pointer<pdat::CellData<double> > dp_ptr =
+ patch->getPatchData(dp_id);
+ pdat::CellData<double> &dp(*dp_ptr);
+
+ tbox::Pointer<pdat::SideData<double> > v_ptr =
+ patch->getPatchData(v_id);
+ pdat::SideData<double> &v(*v_ptr);
+
+ tbox::Pointer<pdat::CellData<double> > cell_visc_ptr
+ = patch->getPatchData(cell_viscosity_id);
+ pdat::CellData<double> &cell_viscosity(*cell_visc_ptr);
+ tbox::Pointer<pdat::NodeData<double> > edge_visc_ptr
+ = patch->getPatchData(edge_viscosity_id);
+ pdat::NodeData<double> &edge_viscosity(*edge_visc_ptr);
+
+ hier::Box pbox=patch->getBox();
+ tbox::Pointer<geom::CartesianPatchGeometry>
+ geom = patch->getPatchGeometry();
+ double dx = geom->getDx()[0];
+ double dy = geom->getDx()[1];
+
+ /* Set an array of bools that tells me whether a point
+ should set the pressure or just let it be. This is
+ needed at coarse/fine boundaries where the pressure
+ is fixed. */
+ hier::Box gbox=p.getGhostBox();
+ std::vector<bool> set_p(gbox.size(),true);
+
+ const tbox::Array<hier::BoundaryBox >&edges
+ =d_cf_boundary[ln]->getEdgeBoundaries(patch->getGlobalId());
+ for(int mm=0; mm<edges.size(); ++mm)
+ for(int j=edges[mm].getBox().lower(1);
+ j<=edges[mm].getBox().upper(1); ++j)
+ for(int i=edges[mm].getBox().lower(0);
+ i<=edges[mm].getBox().upper(0); ++i)
+ set_p[(i-gbox.lower(0))
+ + (gbox.upper(0)-gbox.lower(0)+1)*(j-gbox.lower(1))]=false;
+
+ const tbox::Array<hier::BoundaryBox >&nodes
+ =d_cf_boundary[ln]->getNodeBoundaries(patch->getGlobalId());
+ for(int mm=0; mm<nodes.size(); ++mm)
+ for(int j=nodes[mm].getBox().lower(1);
+ j<=nodes[mm].getBox().upper(1); ++j)
+ for(int i=nodes[mm].getBox().lower(0);
+ i<=nodes[mm].getBox().upper(0); ++i)
+ set_p[(i-gbox.lower(0))
+ + (gbox.upper(0)-gbox.lower(0)+1)*(j-gbox.lower(1))]=false;
+
+ if(geom->getTouchesRegularBoundary(0,0))
+ for(int j=gbox.lower(1); j<=gbox.upper(1); ++j)
+ set_p[(gbox.upper(0)-gbox.lower(0)+1)*(j-gbox.lower(1))]=false;
+
+ if(geom->getTouchesRegularBoundary(0,1))
+ for(int j=gbox.lower(1); j<=gbox.upper(1); ++j)
+ set_p[(gbox.upper(0)-gbox.lower(0))
+ + (gbox.upper(0)-gbox.lower(0)+1)*(j-gbox.lower(1))]=false;
+
+ if(geom->getTouchesRegularBoundary(1,0))
+ for(int i=gbox.lower(0); i<=gbox.upper(0); ++i)
+ set_p[i-gbox.lower(0)]=false;
+
+ if(geom->getTouchesRegularBoundary(1,1))
+ for(int i=gbox.lower(0); i<=gbox.upper(0); ++i)
+ set_p[(i-gbox.lower(0))
+ + (gbox.upper(0)-gbox.lower(0)+1)*(gbox.upper(1)-gbox.lower(1))]=
+ false;
+
+ for(int j=pbox.lower(1); j<=pbox.upper(1); ++j)
+ {
+ /* Do the red-black skip */
+ int i_min=pbox.lower(0) + (abs(pbox.lower(0) + j + rb))%2;
+ for(int i=i_min; i<=pbox.upper(0)+1; i+=2)
+ {
+ pdat::CellIndex center(tbox::Dimension(2));
+ center[0]=i;
+ center[1]=j;
+
+ pdat::CellIndex up(center), down(center), right(center),
+ left(center);
+
+ ++up[1];
+ --down[1];
+ ++right[0];
+ --left[0];
+
+ const pdat::SideIndex
+ center_x(center,0,pdat::SideIndex::Lower),
+ left_x(left,0,pdat::SideIndex::Lower),
+ right_x(right,0,pdat::SideIndex::Lower),
+ center_y(center,1,pdat::SideIndex::Lower),
+ up_y(up,1,pdat::SideIndex::Lower),
+ down_y(down,1,pdat::SideIndex::Lower);
+
+ /* Update p */
+ if(set_p[(i-gbox.lower(0))
+ + (gbox.upper(0)-gbox.lower(0)+1)*(j-gbox.lower(1))])
+ {
+ dp(center)=dRc_dp(pbox,center,left,right,down,up,
+ left_x,right_x,down_y,up_y,
+ cell_viscosity,edge_viscosity,
+ v,dx,dy);
+ }
+ }
+ }
+ }
+ }
+
+
+
+
+ /* p sweep */
+
+ xeqScheduleGhostFillNoCoarse(dp_id,invalid_id,ln);
+ for(int rb=0;rb<2;++rb)
+ {
for (hier::PatchLevel::Iterator pi(*level); pi; pi++)
{
tbox::Pointer<hier::Patch> patch = *pi;
@@ -361,13 +484,6 @@ void SAMRAI::solv::StokesFACOps::smooth_
patch->getPatchData(v_id);
pdat::SideData<double> &v(*v_ptr);
- tbox::Pointer<pdat::CellData<double> > cell_visc_ptr
- = patch->getPatchData(cell_viscosity_id);
- pdat::CellData<double> &cell_viscosity(*cell_visc_ptr);
- tbox::Pointer<pdat::NodeData<double> > edge_visc_ptr
- = patch->getPatchData(edge_viscosity_id);
- pdat::NodeData<double> &edge_viscosity(*edge_visc_ptr);
-
hier::Box pbox=patch->getBox();
tbox::Pointer<geom::CartesianPatchGeometry>
geom = patch->getPatchGeometry();
@@ -460,33 +576,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
maxres=std::max(maxres,std::fabs(delta_R_continuity));
dp(center)=delta_R_continuity*theta_continuity
- /dRc_dp(pbox,center,left,right,down,up,
- left_x,right_x,down_y,up_y,
- cell_viscosity,edge_viscosity,v,dx,dy);
-
- // tbox::plog << "smooth p "
- // << i << " "
- // << j << " "
- // << p(center) << " "
- // << dRc_dp(pbox,center,left,right,down,up,
- // left_x,right_x,down_y,up_y,
- // cell_viscosity,edge_viscosity,v,dx,dy) << " "
- // << dp(center) << " "
- // << delta_R_continuity << " "
- // << dvx_dx << " "
- // << dvy_dy << " "
- // << p_rhs(center) << " "
- // // << v(pdat::SideIndex(center,pdat::SideIndex::X,
- // // pdat::SideIndex::Upper)) << " "
- // // << v(pdat::SideIndex(center,pdat::SideIndex::X,
- // // pdat::SideIndex::Lower)) << " "
- // // << v(pdat::SideIndex(center,pdat::SideIndex::Y,
- // // pdat::SideIndex::Upper)) << " "
- // // << v(pdat::SideIndex(center,pdat::SideIndex::Y,
- // // pdat::SideIndex::Lower)) << " "
-
- // << "\n";
-
+ /dp(center);
p(center)+=dp(center);
}
}
@@ -496,10 +586,9 @@ void SAMRAI::solv::StokesFACOps::smooth_
/* fix v sweep */
+ xeqScheduleGhostFillNoCoarse(dp_id,invalid_id,ln);
for(int rb=0;rb<2;++rb)
{
- // Need to sync
- xeqScheduleGhostFillNoCoarse(dp_id,invalid_id,ln);
for (hier::PatchLevel::Iterator pi(*level); pi; pi++)
{
tbox::Pointer<hier::Patch> patch = *pi;
@@ -604,13 +693,6 @@ void SAMRAI::solv::StokesFACOps::smooth_
up_e(up,pdat::NodeIndex::LowerLeft),
right_e(right,pdat::NodeIndex::LowerLeft);
- // tbox::plog << "Smooth "
- // << i << " "
- // << j << " ";
- // if(set_p[(i-gbox.lower(0))
- // + (gbox.upper(0)-gbox.lower(0)+1)*(j-gbox.lower(1))])
- // tbox::plog << p(center) << " ";
-
/* Update v */
if(set_p[(i-gbox.lower(0))
+ (gbox.upper(0)-gbox.lower(0)+1)*(j-gbox.lower(1))]
@@ -623,9 +705,6 @@ void SAMRAI::solv::StokesFACOps::smooth_
v(center_x)+=(dp(center) - dp(left))
/(dx*dRm_dv(cell_viscosity,edge_viscosity,center,
left,up_e,center_e,dx,dy));
- // tbox::plog << "vx "
- // << v(center_x) << " "
- // << v_rhs(center_x) << " ";
}
if(set_p[(i-gbox.lower(0))
+ (gbox.upper(0)-gbox.lower(0)+1)*(j-gbox.lower(1))]
@@ -638,11 +717,7 @@ void SAMRAI::solv::StokesFACOps::smooth_
v(center_y)+=(dp(center) - dp(down))
/(dy*dRm_dv(cell_viscosity,edge_viscosity,center,
down,right_e,center_e,dy,dx));
- // tbox::plog << "vy "
- // << v(center_y) << " "
- // << v_rhs(center_y) << " ";
}
- // tbox::plog << "\n";
}
}
}
diff -r c6b0d0a8b7f5 -r 7d4eceb1a9e1 input/shear_corner.input
--- a/input/shear_corner.input Fri Mar 11 16:56:26 2011 -0800
+++ b/input/shear_corner.input Fri Mar 11 20:37:48 2011 -0800
@@ -34,16 +34,16 @@ FACStokes {
enable_logging = TRUE // Bool flag to switch logging on/off
max_cycles = 1000 // Max number of FAC cycles to use
residual_tol = 1e-8 // Residual tolerance to solve for
- num_pre_sweeps = 5 // Number of presmoothing sweeps to use
- num_post_sweeps = 5 // Number of postsmoothing sweeps to use
+ num_pre_sweeps = 15 // Number of presmoothing sweeps to use
+ num_post_sweeps = 15 // Number of postsmoothing sweeps to use
smoothing_choice = "Tackley"
coarse_solver_choice = "Tackley"
// smoothing_choice = "Gerya"
// coarse_solver_choice = "Gerya"
coarse_solver_max_iterations = 25
coarse_solver_tolerance = 1e-8
- p_prolongation_method = "P_MDPI_REFINE"
- // p_prolongation_method = "P_REFINE"
+ // p_prolongation_method = "P_MDPI_REFINE"
+ p_prolongation_method = "P_REFINE"
v_prolongation_method = "V_REFINE"
}
bc_coefs {
@@ -98,9 +98,9 @@ StandardTagAndInitialize {
// level_3 = [(0,0),(31,31)]
// level_4 = [(0,0),(63,63)]
// level_5 = [(0,0),(127,127)]
- level_1 = [(2,2),(5,5)]
- level_2 = [(6,6),(9,9)]
- level_3 = [(14,14),(15,15)]
+ level_1 = [(0,0),(5,5)]
+ level_2 = [(0,0),(9,9)]
+ level_3 = [(0,0),(16,16)]
//etc.
}
}
More information about the CIG-COMMITS
mailing list