Commit df68542d authored by Carsten Gräser's avatar Carsten Gräser
Browse files

Re-add old TNNMG implementation for transition

There is stil it lot of code related to the old implementation.
To make it work, we need to re-add the old iteration step. This
is only done for a transition period and will be removed again,
once the new implementation works in the same generality.
parent 78660533
Pipeline #32279 passed with stage
in 1 minute and 56 seconds
install(FILES
acceleratednonlineargsstep.hh
blockseparabletnnmgstep.hh
genericnonlineargs.hh
genericnonlinearjacobi.hh
linearcorrection.hh
......
#ifndef DUNE_TNNMG_ITERATIONSTEPS_BLOCK_SEPARABLE_TNNMG_STEP_HH
#define DUNE_TNNMG_ITERATIONSTEPS_BLOCK_SEPARABLE_TNNMG_STEP_HH
#include <string>
#include <sstream>
#include <vector>
#include <iomanip>
#include <dune/common/timer.hh>
#include "dune/solvers/iterationsteps/iterationstep.hh"
#include "dune/solvers/iterationsteps/lineariterationstep.hh"
/**
* \brief This is the old implementation of a TNNMG step.
*
* \deprecated
* This is only added for a transition period.
*/
template<class TNNMGProblemTypeTEMPLATE, class NonlinearSmootherType>
class TruncatedNonsmoothNewtonMultigrid : public IterationStep<typename TNNMGProblemTypeTEMPLATE::VectorType>
{
public:
typedef TNNMGProblemTypeTEMPLATE TNNMGProblemType;
typedef typename TNNMGProblemType::MatrixType MatrixType;
typedef typename TNNMGProblemType::VectorType VectorType;
typedef typename TNNMGProblemType::Linearization Linearization;
typedef typename Linearization::MatrixType CoarseMatrixType;
typedef typename Linearization::VectorType CoarseVectorType;
typedef typename Linearization::BitVectorType CoarseBitVectorType;
typedef LinearIterationStep<CoarseMatrixType, CoarseVectorType, CoarseBitVectorType> CoarseLinearIterationStep;
struct Statistics
{
int iterationCount;
double smoothingTime;
double linearizationTime;
double coarseCorrectionTime;
double postProcessingTime;
};
TruncatedNonsmoothNewtonMultigrid() :
problem_(0),
nonlinearSmoother_(0),
linearIterationStep_(0),
preSmoothingSteps_(1),
postSmoothingSteps_(0),
linearIterationSteps_(1)
{};
TruncatedNonsmoothNewtonMultigrid(CoarseLinearIterationStep& linearIterationStep, NonlinearSmootherType& nonlinearSmoother) :
problem_(0),
nonlinearSmoother_(&nonlinearSmoother),
linearIterationStep_(&linearIterationStep),
preSmoothingSteps_(1),
postSmoothingSteps_(0),
linearIterationSteps_(1)
{};
void setNonlinearSmoother(NonlinearSmootherType& nonlinearSmoother)
{
nonlinearSmoother_ = &nonlinearSmoother;
};
void setLinearIterationStep(CoarseLinearIterationStep& linearIterationStep)
{
linearIterationStep_ = &linearIterationStep;
};
void setProblem(VectorType& x, TNNMGProblemType& problem)
{
this->x_ = &x;
this->problem_ = &problem;
};
void setSmoothingSteps(int pre, int linear, int post)
{
preSmoothingSteps_ = pre;
linearIterationSteps_ = linear;
postSmoothingSteps_ = post;
};
virtual void preprocess()
{
outStream_.str("");
outStream_ << problem_->getOutput(true);
outStream_ << " step size ";
outStream_ << " energy ";
statistics_.iterationCount = 0;
statistics_.smoothingTime = 0.0;
statistics_.linearizationTime = 0.0;
statistics_.coarseCorrectionTime = 0.0;
statistics_.postProcessingTime = 0.0;
}
void iterate()
{
// clear previous output data
outStream_.str("");
// do some time measurement
Dune::Timer timer;
++statistics_.iterationCount;
VectorType& x = *x_;
TNNMGProblemType& problem = *problem_;
// apply nonlinear smoother (pre smoothing)
timer.reset();
nonlinearSmoother_->setProblem(x, problem);
nonlinearSmoother_->ignoreNodes_ = ignoreNodes_;
for(int i=0; i<preSmoothingSteps_; ++i)
nonlinearSmoother_->iterate();
x = nonlinearSmoother_->getSol();
statistics_.smoothingTime += timer.elapsed();
// assemble and truncate linear system for coarse correction
timer.reset();
Linearization linearization;
problem.assembleTruncate(x, linearization, *ignoreNodes_);
statistics_.linearizationTime += timer.elapsed();
// solve linear system for coarse correction
CoarseVectorType v(linearization.b.size());
v = 0.0;
// apply linear solver to coarse problem
timer.reset();
{
linearIterationStep_->setProblem(linearization.A, v, linearization.b);
linearIterationStep_->ignoreNodes_ = &(linearization.ignore);
linearIterationStep_->preprocess();
for(int i=0; i<linearIterationSteps_; ++i)
linearIterationStep_->iterate();
v = linearIterationStep_->getSol();
}
statistics_.coarseCorrectionTime += timer.elapsed();
// post processing
timer.reset();
// project correction
VectorType projected_v;
problem.projectCoarseCorrection(x, v, projected_v, linearization);
// compute damping parameter
double alpha = problem.computeDampingParameter(x, projected_v);
// apply damped coarse correction
x.axpy(alpha, projected_v);
statistics_.postProcessingTime += timer.elapsed();
// apply nonlinear smoother (post smoothing)
timer.reset();
nonlinearSmoother_->setProblem(x, problem);
for(int i=0; i<postSmoothingSteps_; ++i)
nonlinearSmoother_->iterate();
x = nonlinearSmoother_->getSol();
statistics_.smoothingTime += timer.elapsed();
outStream_ << problem.getOutput();
outStream_.setf(std::ios::scientific);
outStream_ << std::setw(15) << alpha;
outStream_ << std::setw(15) << problem.computeEnergy(x);
}
VectorType getSol()
{
return *x_;
}
std::string getOutput() const
{
std::string s = outStream_.str();
outStream_.str("");
return s;
}
const Statistics& getStatistics() const
{
return statistics_;
}
using IterationStep<VectorType>::ignoreNodes_;
private:
Statistics statistics_;
using IterationStep<VectorType>::x_;
TNNMGProblemType* problem_;
NonlinearSmootherType* nonlinearSmoother_;
CoarseLinearIterationStep* linearIterationStep_;
int preSmoothingSteps_;
int postSmoothingSteps_;
int linearIterationSteps_;
mutable std::ostringstream outStream_;
};
#endif // #define DUNE_TNNMG_ITERATIONSTEPS_BLOCK_SEPARABLE_TNNMG_STEP_HH
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment