33 #ifndef __H__UG__LIB_DISC__OPERATOR__NON_LINEAR_OPERATOR__NI_SOLVER__NI_IMPL__
34 #define __H__UG__LIB_DISC__OPERATOR__NON_LINEAR_OPERATOR__NI_SOLVER__NI_IMPL__
43 #define PROFILE_NESTED_ITER
44 #ifdef PROFILE_NESTED_ITER
45 #define NESTED_ITER_PROFILE_FUNC() PROFILE_FUNC_GROUP("NestedIteration")
46 #define NESTED_ITER_PROFILE_BEGIN(name) PROFILE_BEGIN_GROUP(name, "NestedIteration")
47 #define NESTED_ITER_PROFILE_END() PROFILE_END()
49 #define NESTED_ITER_PROFILE_FUNC()
50 #define NESTED_ITER_PROFILE_BEGIN(name)
51 #define NESTED_ITER_PROFILE_END()
56 template <
typename TDomain,
typename TAlgebra>
65 m_bUseAdaptiveRefinement(true),
66 m_maxSteps(100), m_TOL(1e-3), m_absTOL(1e-12)
69 template <
typename TDomain,
typename TAlgebra>
72 m_spLinearSolver(NULL),
78 m_bUseAdaptiveRefinement(true),
79 m_maxSteps(100), m_TOL(1e-3), m_absTOL(1e-12)
82 template <
typename TDomain,
typename TAlgebra>
85 m_spLinearSolver(NULL),
91 m_bUseAdaptiveRefinement(true),
92 m_maxSteps(100), m_TOL(1e-3), m_absTOL(1e-12)
97 template <
typename TDomain,
typename TAlgebra>
100 m_spLinearSolver(NULL),
108 m_bUseAdaptiveRefinement(true),
109 m_maxSteps(100), m_TOL(1e-3), m_absTOL(1e-12)
115 template <
typename TDomain,
typename TAlgebra>
118 m_spLinearSolver(NULL),
122 m_spDomErr(spDomErr),
126 m_bUseAdaptiveRefinement(true)
133 template <
typename TDomain,
typename TAlgebra>
137 m_N = N.template cast_dynamic<AssembledOperator<TAlgebra> >();
139 UG_THROW(
"NestedIterationSolver: currently only works for AssembledDiscreteOperator.");
141 m_spAss = m_N->discretization();
146 template <
typename TDomain,
typename TAlgebra>
156 template <
typename TDomain,
typename TAlgebra>
159 UG_LOG(
"Computing error... "<< std::endl);
167 UG_ASSERT(spDomainEstimator.
valid(),
"Huhh: No estimatable object???")
168 UG_ASSERT(m_spRefiner.valid(),
"Huhh: Invalid refiner???");
171 if (m_spElemError.valid()) {
173 spDomainEstimator->calc_error(u, u.
dd(), m_spElemError.get());
176 spDomainEstimator->calc_error(u, u.
dd());
180 if (bClearMarks) m_spRefiner->clear_marks();
181 spDomainEstimator->mark_with_strategy(*m_spRefiner, spMarking);
182 UG_LOG(
"Estimated error: " << spMarking->global_estimated_error());
191 template <
typename TDomain,
typename TAlgebra>
207 template <
typename TDomain,
typename TAlgebra>
224 if(m_spLinearSolver.invalid())
225 UG_THROW(
"NestedIterationSolver::apply: Linear Solver not set.");
228 if(m_J.invalid() || m_J->discretization() != m_spAss) {
231 m_J->set_level(m_N->level());
233 UG_LOG(
"N level: "<<m_N->level()<< std::endl);
244 sprintf(ext,
"_call%03d", m_dgbCall);
245 std::string
name(
"NESTED_ITER_StartSolution");
247 write_debug(u,
name.c_str());
252 m_topLevel = surfLevel+m_maxSteps;
253 for (
int lev=surfLevel; lev<m_topLevel; ++lev)
265 m_spAss->assemble_linear(*m_J, *spD, surfGridLevel);
266 m_spAss->adjust_solution(u, surfGridLevel);
269 sprintf(ext,
"_call%03d_iter%03d", m_dgbCall, loopCnt);
272 std::string name0(
"NESTED_ITER_InitialSolution"); name0.append(ext);
273 write_debug(u, name0.c_str());
274 std::string
name(
"NESTED_ITER_InitialDefect");
name.append(ext);
275 write_debug(*spD,
name.c_str());
278 std::string matname(
"NESTED_ITER_Jacobian");
280 write_debug(m_J->get_matrix(), matname.c_str());
286 if(!m_spLinearSolver->init(m_J, u))
288 UG_LOG(
"ERROR in 'NestedIterationSolver::apply': Cannot init Inverse Linear "
289 "Operator for Jacobi-Operator.\n");
293 }
UG_CATCH_THROW(
"NestedIterationSolver::apply: Initialization of Linear Solver failed.");
298 if(!m_spLinearSolver->apply(u, *spD))
300 UG_LOG(
"ERROR in 'NestedIterationSolver::apply': Cannot apply Inverse Linear "
301 "Operator for Jacobi-Operator.\n");
305 }
UG_CATCH_THROW(
"NestedIterationSolver::apply: Application of Linear Solver failed.");
308 m_spAss->adjust_solution(u, surfGridLevel);
311 std::string
name(
"NESTED_ITER_FinalDefect");
name.append(ext);
312 write_debug(*spD,
name.c_str());
313 std::string name3(
"NESTED_ITER_FinalSolution"); name3.append(ext);
314 write_debug(u, name3.c_str());
321 if (use_adaptive_refinement() ==
false) {
322 UG_LOG(
"SUCCESS: Non-adaptive version!" << std::endl);
327 this->estimate_and_mark_domain(usol, m_spRefinementMarking);
330 if (m_spElemError.valid())
332 std::string
name(
"NESTED_ITER_Error");
335 outError.
print(
name.c_str(), *m_spElemError);
339 const number err = m_spRefinementMarking->global_estimated_error();
340 double desiredTOL = (m_spAssociatedSpace.valid()) ? m_TOL*m_spAssociatedSpace->norm(usol) + m_absTOL : m_TOL;
342 UG_LOG(
"NI *** Desired tolerance: " << m_TOL <<
" * "<< m_spAssociatedSpace->norm(usol) <<
" + "<< m_absTOL <<std::endl);
348 UG_LOG(
"SUCCESS: Error smaller than tolerance: " << err <<
" < "<< desiredTOL << std::endl);
352 UG_LOG(
"FAILED: Error (still) greater than tolerance: " << err <<
" > "<< desiredTOL << std::endl);
353 m_spRefiner->refine();
362 template <
typename TDomain,
typename TAlgebra>
365 std::stringstream ss;
366 ss <<
"NestedIterationSolver\n";
367 ss <<
" LinearSolver: ";
368 if(m_spLinearSolver.valid()) ss <<
ConfigShift(m_spLinearSolver->config_string()) <<
"\n";
369 else ss <<
" NOT SET!\n";
location name
Definition: checkpoint_util.lua:128
Definition: smart_pointer.h:108
bool valid() const
returns true if the pointer is valid, false if not.
Definition: smart_pointer.h:206
matrix operator based on the assembling of a problem
Definition: assembled_linear_operator.h:60
Definition: assembled_non_linear_operator.h:43
TAlgebra::vector_type vector_type
type of vector
Definition: debug_writer.h:360
domain discretization implementing the interface
Definition: domain_disc.h:638
represents numerical solutions on a grid using an algebraic vector
Definition: grid_function.h:121
const GridLevel & grid_level() const
returns the grid level
Definition: grid_function.h:247
SmartPtr< DoFDistribution > dd()
Definition: grid_function.h:237
Definition: grid_level.h:42
bool top() const
returns if top level
Definition: grid_level.h:83
bool is_surface() const
returns if type is surface
Definition: grid_level.h:89
int level() const
returns the level
Definition: grid_level.h:74
Interface providing Jacobian and Defect of a discretization.
Definition: assemble_interface.h:110
Abstract base class for element marking (in adaptive refinement)
Definition: error_elem_marking_strategy.h:125
describes an inverse linear mapping X->Y
Definition: linear_operator_inverse.h:80
describes a mapping X->Y
Definition: operator.h:86
linear solver using abstract preconditioner interface
Definition: linear_solver.h:57
SmartPtr< AssembledOperator< algebra_type > > m_N
assembling routine
Definition: nested_iteration.h:193
virtual std::string config_string() const
returns information about configuration parameters
Definition: nested_iteration_impl.h:363
virtual bool apply(vector_type &u)
apply Operator, i.e. N^{-1}(0) = u
Definition: nested_iteration_impl.h:208
SmartPtr< IAssemble< TAlgebra > > m_spAss
assembling
Definition: nested_iteration.h:197
number coarsen_domain(const grid_function_type &u)
Coarsen all elements.
Definition: nested_iteration_impl.h:192
virtual bool init(SmartPtr< IOperator< vector_type > > N)
Definition: nested_iteration_impl.h:134
void estimate_and_mark_domain(const grid_function_type &u, SmartPtr< IElementMarkingStrategy< TDomain > > spMarking, bool bClearMarks=true)
Refines domain and provides error estimate.
Definition: nested_iteration_impl.h:157
NestedIterationSolver()
default constructor
Definition: nested_iteration_impl.h:71
virtual bool prepare(vector_type &u)
prepare Operator
Definition: nested_iteration_impl.h:147
output writer to the VTK file format
Definition: vtkoutput.h:150
void print(const char *filename, TFunction &u, int step, number time, bool makeConsistent)
Definition: vtkoutput_impl.h:184
function common sprintf(s,...)
#define UG_ASSERT(expr, msg)
Definition: assert.h:70
#define UG_CATCH_THROW(msg)
Definition: error.h:64
#define UG_THROW(msg)
Definition: error.h:57
#define UG_LOG(msg)
Definition: log.h:367
double number
Definition: types.h:124
string ConfigShift(string s)
Definition: string_util.cpp:457
#define NESTED_ITER_PROFILE_END()
Definition: nested_iteration_impl.h:47
#define NESTED_ITER_PROFILE_BEGIN(name)
Definition: nested_iteration_impl.h:46
SmartPtr< T, FreePolicy > make_sp(T *inst)
returns a SmartPtr for the passed raw pointer
Definition: smart_pointer.h:836