[cig-commits] commit: fix the adaptivity.
Mercurial
hg at geodynamics.org
Thu Jun 7 13:35:39 PDT 2012
changeset: 258:39366675b6ee
user: Sylvain Barbot <sbarbot at caltech.edu>
date: Tue Apr 10 16:52:35 2012 -0700
files: input/sinker.input src/FACStokes/applyGradientDetector.C
description:
fix the adaptivity.
diff -r 47bb8c85e9aa -r 39366675b6ee input/sinker.input
--- a/input/sinker.input Tue Apr 10 13:39:52 2012 -0700
+++ b/input/sinker.input Tue Apr 10 16:52:35 2012 -0700
@@ -28,8 +28,8 @@ FACStokes {
// It owns the solver and contains the code to set up the solver.
// The inputs for FACStokes is simply the inputs for the individual
// parts owned by the FACStokes class.
- adaption_threshold = 1.0e-3
- min_full_refinement_level = 5
+ adaption_threshold = 1.0e-2
+ min_full_refinement_level = 3
lambda_ijk= 2, 2
lambda_coord_min= -0.001, -0.001
@@ -191,7 +191,7 @@ PatchHierarchy {
// [level 0 entry]
// etc....
// }
- max_levels = 7
+ max_levels = 6
proper_nesting_buffer = 2, 2, 2, 2, 2, 2
ratio_to_coarser {
level_1 = 2, 2
diff -r 47bb8c85e9aa -r 39366675b6ee src/FACStokes/applyGradientDetector.C
--- a/src/FACStokes/applyGradientDetector.C Tue Apr 10 13:39:52 2012 -0700
+++ b/src/FACStokes/applyGradientDetector.C Tue Apr 10 16:52:35 2012 -0700
@@ -34,39 +34,47 @@ void SAMRAI::FACStokes::applyGradientDet
{
TBOX_ERROR("Data index " << tag_index << " is not cell int data.\n");
}
- tbox::Pointer<hier::PatchData> soln_data = patch.getPatchData(p_id);
+ tbox::Pointer<hier::PatchData> soln_data = patch.getPatchData(v_id);
if (soln_data.isNull())
{
- TBOX_ERROR("Data index " << p_id << " does not exist for patch.\n");
+ TBOX_ERROR("Data index " << v_id << " does not exist for patch.\n");
}
- tbox::Pointer<pdat::CellData<double> > soln_cell_data_ = soln_data;
- if (soln_cell_data_.isNull())
+ tbox::Pointer<pdat::SideData<double> > soln_side_data_ = soln_data;
+ if (soln_side_data_.isNull())
{
- TBOX_ERROR("Data index " << p_id << " is not cell int data.\n");
+ TBOX_ERROR("Data index " << v_id << " is not side data.\n");
}
- pdat::CellData<double>& soln_cell_data = *soln_cell_data_;
+ pdat::SideData<double>& v = *soln_side_data_;
pdat::CellData<int>& tag_cell_data = *tag_cell_data_;
- pdat::CellData<double> estimate_data(patch.getBox(),1,
- hier::IntVector(d_dim, 0));
- computeAdaptionEstimate(estimate_data,soln_cell_data);
tag_cell_data.fill(0);
for (pdat::CellIterator ci(patch.getBox()); ci; ci++)
{
const pdat::CellIndex cell_index(*ci);
- if (maxestimate < estimate_data(cell_index))
- maxestimate=estimate_data(cell_index);
- // tbox::plog << "estimate "
- // << cell_index << " "
- // << d_adaptation_threshold << " "
- // << estimate_data(cell_index) << " "
- // << std::boolalpha
- // << (estimate_data(cell_index) > d_adaptation_threshold)
- // << " "
- // << "\n";
- if (estimate_data(cell_index) > d_adaption_threshold
- || ln<min_full_refinement_level)
+ double curve(0.);
+ for (int ix=0; ix<d_dim.getValue(); ++ix)
+ {
+ const pdat::SideIndex x(cell_index,ix,pdat::SideIndex::Lower);
+ for (int d=0; d<d_dim.getValue(); ++d){
+ const hier::Index ip(1,0), jp(0,1);
+ const hier::Index pp[]={ip,jp};
+
+ curve=std::max(curve,std::abs(v(x+pp[ix]+pp[ix])-v(x+pp[ix])-v(x)+v(x-pp[ix])));
+ }
+ }
+ tbox::plog << "estimate "
+ << cell_index << " "
+ << d_adaption_threshold << " "
+ << curve << " "
+ << std::boolalpha
+ << (curve > d_adaption_threshold)
+ << " "
+ << "\n";
+
+ if (maxestimate < curve)
+ maxestimate=curve;
+ if (curve > d_adaption_threshold || ln<min_full_refinement_level)
{
tag_cell_data(cell_index) = 1;
++ntag;
@@ -79,89 +87,3 @@ void SAMRAI::FACStokes::applyGradientDet
tbox::plog << "Max estimate is " << maxestimate << "\n";
}
-
-void SAMRAI::FACStokes::computeAdaptionEstimate
-(pdat::CellData<double>& estimate_data,
- const pdat::CellData<double>& soln_cell_data) const
-{
- const int* lower = &estimate_data.getBox().lower()[0];
- const int* upper = &estimate_data.getBox().upper()[0];
- if (d_dim == tbox::Dimension(2))
- {
- MDA_AccessConst<double, 2, MDA_OrderColMajor<2> > co =
- pdat::ArrayDataAccess::access<2, double>(soln_cell_data.getArrayData());
- MDA_Access<double, 2, MDA_OrderColMajor<2> > es =
- pdat::ArrayDataAccess::access<2, double>(estimate_data.getArrayData());
- int i, j;
- double estimate, est0, est1, est2, est3, est4, est5;
- for (j = lower[1]; j <= upper[1]; ++j)
- {
- for (i = lower[0]; i <= upper[0]; ++i)
- {
- est0=tbox::MathUtilities<double>::Abs(co(i+1,j) + co(i-1,j)
- - 2*co(i,j));
-
- est1=tbox::MathUtilities<double>::Abs(co(i,j+1) + co(i,j-1)
- - 2*co(i,j));
- est2=0.5 * tbox::MathUtilities<double>::Abs(co(i+1,j+1)
- + co(i-1,j-1)
- - 2*co(i,j));
- est3=0.5 * tbox::MathUtilities<double>::Abs(co(i+1,j-1)
- + co(i-1,j+1)
- - 2*co(i,j));
- est4=tbox::MathUtilities<double>::Max(est0,est1);
- est5=tbox::MathUtilities<double>::Max(est2,est3);
- estimate=tbox::MathUtilities<double>::Max(est4,est5);
- es(i,j)=estimate;
- }
- }
- }
- else if (d_dim == tbox::Dimension(3))
- {
- MDA_AccessConst<double, 3, MDA_OrderColMajor<3> > co =
- pdat::ArrayDataAccess::access<3, double>(soln_cell_data.getArrayData());
- MDA_Access<double, 3, MDA_OrderColMajor<3> > es =
- pdat::ArrayDataAccess::access<3, double>(estimate_data.getArrayData());
- int i, j, k;
- double estimate, est0, est1, est2, est3, est4, est5, est6, est7, est8,
- esta, estb, estc, estd, este, estf, estg;
- for (k = lower[2]; k <= upper[2]; ++k)
- for (j = lower[1]; j <= upper[1]; ++j)
- for (i = lower[0]; i <= upper[0]; ++i)
- {
- est0=tbox::MathUtilities<double>::Abs(co(i+1,j,k) + co(i-1,j,k)
- - 2*co(i,j,k));
- est1=tbox::MathUtilities<double>::Abs(co(i,j+1,k) + co(i,j-1,k)
- - 2*co(i,j,k));
- est2=tbox::MathUtilities<double>::Abs(co(i,j,k+1) + co(i,j,k-1)
- - 2*co(i,j,k));
- est3=0.5 * tbox::MathUtilities<double>::Abs(co(i,j+1,k+1)
- + co(i,j-1,k-1)
- - 2*co(i,j,k));
- est4=0.5 * tbox::MathUtilities<double>::Abs(co(i,j+1,k-1)
- + co(i,j-1,k+1)
- - 2*co(i,j,k));
- est5=0.5 * tbox::MathUtilities<double>::Abs(co(i+1,j,k+1)
- + co(i-1,j,k-1)
- - 2*co(i,j,k));
- est6=0.5 * tbox::MathUtilities<double>::Abs(co(i+1,j,k-1)
- + co(i-1,j,k+1)
- - 2*co(i,j,k));
- est7=0.5 * tbox::MathUtilities<double>::Abs(co(i+1,j+1,k)
- + co(i-1,j-1,k)
- - 2*co(i,j,k));
- est8=0.5 * tbox::MathUtilities<double>::Abs(co(i+1,j-1,k)
- + co(i-1,j+1,k)
- - 2*co(i,j,k));
- esta=tbox::MathUtilities<double>::Max(est0,est1);
- estb=tbox::MathUtilities<double>::Max(est2,est3);
- estc=tbox::MathUtilities<double>::Max(est4,est5);
- estd=tbox::MathUtilities<double>::Max(est6,est7);
- este=tbox::MathUtilities<double>::Max(esta,estb);
- estf=tbox::MathUtilities<double>::Max(estc,estd);
- estg=tbox::MathUtilities<double>::Max(este,estf);
- estimate=tbox::MathUtilities<double>::Max(estg,est8);
- es(i,j,k)=estimate;
- }
- }
-}
More information about the CIG-COMMITS
mailing list