Plugins
time_integrator.hpp
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2014-2020: G-CSC, Goethe University Frankfurt
3  * Author: Arne Naegel, Andreas Kreienbuehl
4  *
5  * This file is part of UG4.
6  *
7  * UG4 is free software: you can redistribute it and/or modify it under the
8  * terms of the GNU Lesser General Public License version 3 (as published by the
9  * Free Software Foundation) with the following additional attribution
10  * requirements (according to LGPL/GPL v3 §7):
11  *
12  * (1) The following notice must be displayed in the Appropriate Legal Notices
13  * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14  *
15  * (2) The following notice must be displayed at a prominent place in the
16  * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17  *
18  * (3) The following bibliography is recommended for citation and must be
19  * preserved in all covered files:
20  * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21  * parallel geometric multigrid solver on hierarchically distributed grids.
22  * Computing and visualization in science 16, 4 (2013), 151-164"
23  * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24  * flexible software system for simulating pde based models on high performance
25  * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26  *
27  * This program is distributed in the hope that it will be useful,
28  * but WITHOUT ANY WARRANTY; without even the implied warranty of
29  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30  * GNU Lesser General Public License for more details.
31  */
32 
33 #ifndef TIME_INTEGRATOR_HPP_
34 #define TIME_INTEGRATOR_HPP_
35 
36 #if __cplusplus >= 201103L
37 #define OVERRIDE override
38 #else
39 #define OVERRIDE
40 #endif
41 
42 // std headers.
43 #include <string>
44 
45 // UG4 headers.
50 #include "lib_disc/assemble_interface.h" // TODO: missing IAssemble in following file:
59 #include "lib_disc/function_spaces/grid_function_util.h" // SaveVectorForConnectionViewer
60 #include "lib_disc/function_spaces/interpolate.h" //Interpolate
61 #include "lib_disc/function_spaces/integrate.h" //Integral
62 #include "lib_disc/function_spaces/grid_function.h" //GridFunction
63 #include "lib_disc/io/vtkoutput.h"
64 
65 
66 
67 // My headers.
68 #include "time_extrapolation.h"
69 #include "../limex_tools.h"
70 
71 namespace ug {
72 
73 
75 template<class TDomain, class TAlgebra>
77 : public ITimeIntegratorObserver<TDomain, TAlgebra>,
78  public ITimeIntegratorStageObserver_start<TDomain, TAlgebra>,
79  public ITimeIntegratorStageObserver_finalize<TDomain, TAlgebra>,
80  public ITimeIntegratorStageObserver_end<TDomain, TAlgebra>
81 {
82 public:
86 
88  : m_sp_vtk(SPNULL), m_filename("0000"), m_startTime(0.0), m_plotStep(0.0) {}
89 
90  VTKOutputObserver(const char *filename, SmartPtr<vtk_type> vtk)
91  : m_sp_vtk(vtk), m_filename(filename), m_startTime(0.0), m_plotStep(0.0) {}
92 
93  VTKOutputObserver(const char *filename, SmartPtr<vtk_type> vtk, number plotStep)
94  : m_sp_vtk(vtk), m_filename(filename), m_startTime(0.0), m_plotStep(plotStep) {}
95 
97  { m_sp_vtk = SPNULL; }
98 
99 
100  void set_output_scales(const std::vector<number>& vScales)
101  {
102  m_vOutputScales = vScales;
103  }
104 
105 
107  {
108  return true;
109  }
110 
111 
113  {
114  if (!m_sp_vtk.valid())
115  return false;
116 
117  writeToFile(u, step, time);
118 
119  m_startTime = time;
120  m_uOld = u->clone();
121 
122  return true;
123  }
124 
125 
127  {
128  if (!m_sp_vtk.valid())
129  return false;
130 
131  if (m_plotStep == 0.0)
132  {
133  writeToFile(uNew, step, time);
134  return true;
135  }
136 
137  // otherwise, only plot data at multiples of given time step (interpolate if necessary)
138  number rem = fmod(time - dt - m_startTime, m_plotStep);
139  number nextPlotTimePt = time - dt - rem + m_plotStep;
140  int nextStep = (int) ((nextPlotTimePt - m_startTime + 0.5 * m_plotStep) / m_plotStep);
141 
142  if (nextPlotTimePt > time)
143  {
144  m_uOld = uNew->clone();
145  return true;
146  }
147 
148  SmartPtr<grid_function_type> uCur = uNew->clone_without_values();
149  while (nextPlotTimePt <= time)
150  {
151  number alpha = (time - nextPlotTimePt) / dt;
152  VecScaleAdd(static_cast<typename TAlgebra::vector_type&>(*uCur),
153  alpha, static_cast<typename TAlgebra::vector_type&>(*m_uOld),
154  1.0 - alpha, static_cast<typename TAlgebra::vector_type&>(*uNew));
155 
156  writeToFile(uCur, nextStep, nextPlotTimePt);
157 
158  nextPlotTimePt = (++nextStep) * m_plotStep;
159  }
160 
161  m_uOld = uNew->clone();
162  return true;
163  }
164 
165 
167  {
168  if (!m_sp_vtk.valid())
169  return false;
170 
171  m_sp_vtk->write_time_pvd(m_filename.c_str(), *u);
172  return true;
173  }
174 
175 
176 private:
178  {
179  if (m_vOutputScales.size())
180  {
181  SmartPtr<grid_function_type> uTmp = u->clone_without_values();
182  ScaleGF<grid_function_type>(uTmp, u, m_vOutputScales);
183  m_sp_vtk->print(m_filename.c_str(), *uTmp, step, time);
184  }
185  else
186  m_sp_vtk->print(m_filename.c_str(), *u, step, time);
187  }
188 
189 protected:
192  std::string m_filename;
195  std::vector<number> m_vOutputScales;
196 };
197 
198 
199 
201 template<class TDomain, class TAlgebra>
203 : public ITimeIntegratorObserver<TDomain, TAlgebra>
204 {
205 public:
208 
209  ConnectionViewerOutputObserver(const char *filename)
210  : m_filename(filename), m_outputTime(-1.0) {}
211 
212  ConnectionViewerOutputObserver(const char *filename, number t_out)
213  : m_filename(filename), m_outputTime(t_out) {}
214 
216  {}
217 
218  bool step_process(SmartPtr<grid_function_type> uNew, /*SmartPtr<grid_function_type> uOld, */int step, number time, number dt) OVERRIDE
219  {
220  // quit, if time does not match
221  if (m_outputTime >=0.0 && time != m_outputTime) return true;
222 
223  SaveVectorForConnectionViewer<grid_function_type>(*uNew, m_filename.c_str());
224  return true;
225  }
226 
227 protected:
228  std::string m_filename;
230 
231 };
232 
233 
234 #if 0
235 template <typename TData, typename TDataIn1, typename TDataIn2>
236 class LuaFunction2 // : public IFunction<TData, TDataIn1, typename TDataIn2>
237 {
238  public:
240  LuaFunction2();
241  virtual ~LuaFunction2() {};
242 
244 
249  void set_lua_callback(const char* luaCallback, size_t numArgs);
250 
252  virtual void operator() (TData& out, int numArgs1, int numArgs2,...);
253 
254  protected:
256  std::string m_cbValueName;
257 
259  int m_cbValueRef;
260 
262  lua_State* m_L;
263 
265  size_t m_numArgs;
266 };
267 
268 
269 template <typename TData, typename TDataIn1, typename TDataIn2>
270 LuaFunction2<TData,TDataIn1,TDataIn2>::LuaFunction2() : m_numArgs(0)
271 {
273  m_cbValueRef = LUA_NOREF;
274 }
275 
276 template <typename TData, typename TDataIn1, typename TDataIn2>
277 void LuaFunction2<TData,TDataIn1,TDataIn2>::set_lua_callback(const char* luaCallback, size_t numArgs)
278 {
279 // store name (string) of callback
280  m_cbValueName = luaCallback;
281 
282 // obtain a reference
283  lua_getglobal(m_L, m_cbValueName.c_str());
284 
285 // make sure that the reference is valid
286  if(lua_isnil(m_L, -1)){
287  UG_THROW("LuaFunction::set_lua_callback(...):"
288  "Specified lua callback does not exist: " << m_cbValueName);
289  }
290 
291 // store reference to lua function
292  m_cbValueRef = luaL_ref(m_L, LUA_REGISTRYINDEX);
293 
294 // remember number of arguments to be used
295  m_numArgs = numArgs;
296 }
297 #endif
298 
299 /*
300 SmartUserDataWrapper* CreateNewUserData(lua_State* L, const SmartPtr<void>& ptr,
301  const char* metatableName)
302 {
303 // create the userdata
304  SmartUserDataWrapper* udata = (SmartUserDataWrapper*)lua_newuserdata(L,
305  sizeof(SmartUserDataWrapper));
306  new(udata) SmartUserDataWrapper;
307 
308 // associate the object with the userdata.
309  udata->smartPtr = ptr;
310  udata->type = SMART_POINTER;
311 
312 // associate the metatable (userdata is already on the stack)
313  luaL_getmetatable(L, metatableName);
314  lua_setmetatable(L, -2);
315 
316  return udata;
317 }
318 */
319 /*
320 template <typename TData, typename TDataIn1, typename TDataIn2>
321 void LuaFunction2<TData,TDataIn1,TDataIn2>::operator() (TData& out, int numArgs1, SmartPtr<TDataIn1> valsArgs1[],
322  int numArgs2, ...)
323 {
324  PROFILE_CALLBACK_BEGIN(operatorBracket);
325 
326  UG_ASSERT((numArgs1+numArgs2) == (int)m_numArgs, "Number of arguments mismatched.");
327 
328  // push the callback function on the stack
329  lua_rawgeti(m_L, LUA_REGISTRYINDEX, m_cbValueRef);
330 
331  // get list of arguments
332  va_list ap;
333  va_start(ap, numArgs2);
334 
335  // read all arguments and push them to the lua stack
336  for(int i = 0; i < numArgs1; ++i)
337  {
338 
339  CreateNewUserData(m_L, &valArgs1[i], "");
340 
341  }
342  for(int i = 0; i < numArgs2; ++i)
343  {
344  TDataIn2 val = va_arg(ap, TDataIn2);
345  lua_traits<TDataIn2>::push(m_L, val);
346  }
347 
348 
349  // end read in of parameters
350  va_end(ap);
351 
352  // compute total args size
353  size_t argSize = lua_traits<TDataIn1>::size * numArgs1;
354  argSize += lua_traits<TDataIn2>::size * numArgs2;
355 
356  // compute total return size
357  size_t retSize = lua_traits<TData>::size;
358 
359  // call lua function
360  if(lua_pcall(m_L, argSize, retSize, 0) != 0)
361  UG_THROW("LuaFunction::operator(...): Error while "
362  "running callback '" << m_cbValueName << "',"
363  " lua message: "<< lua_tostring(m_L, -1));
364 
365  try{
366  // read return value
367  lua_traits<TData>::read(m_L, out);
368  UG_COND_THROW(IsFiniteAndNotTooBig(out)==false, out);
369  }
370  UG_CATCH_THROW("LuaFunction::operator(...): Error while running "
371  "callback '" << m_cbValueName << "'");
372 
373  // pop values
374  lua_pop(m_L, retSize);
375 
376  PROFILE_CALLBACK_END();
377 }
378 */
379 
380 template<class TDomain, class TAlgebra>
382 : public ITimeIntegratorObserver<TDomain, TAlgebra>
383 {
384 public:
388 
390  { m_spReference = spExactSol; }
391 
392 #ifdef UG_FOR_LUA
393  PlotRefOutputObserver(const char *ExactSol)
394  : m_sp_vtk(SPNULL)
396 
397  PlotRefOutputObserver(const char *ExactSol, SmartPtr<vtk_type> vtk)
398  : m_sp_vtk(vtk)
399  { m_spReference = make_sp(new LuaUserData<number, grid_function_type::dim>(ExactSol)); }
400 
401 #endif
402 
404  {}
405 
406  // TODO: replace by call 'func (SmartPtr<G> u, int step, number dt, number t)'
407  bool step_process(SmartPtr<grid_function_type> uNew, /* SmartPtr<grid_function_type> uOld, */ int step, number time, number dt) OVERRIDE
408  {
409  UG_LOG("L2Error(\t"<< time << "\t) = \t" << L2Error(m_spReference, uNew, "c", time, 4) << std::endl);
410  if (m_sp_vtk.valid())
411  {
412  SmartPtr<grid_function_type> ref = uNew->clone();
413  Interpolate<grid_function_type> (m_spReference, ref, "c", time);
414  m_sp_vtk->print("MyReference", *ref, step, time);
415  return true;
416  }
417 
418  return false;
419 
420  }
421 
422 protected:
423  // TODO: replace by appropriate call-back
426 };
427 
429 
432 template<class TDomain, class TAlgebra>
434 : public ITimeIntegratorObserver<TDomain, TAlgebra>
435 {
436 public:
440 protected:
442  {
443  IntegralSpecs(const char* cmp, const char* subsets, int quadOrder, const char *idString) :
444  m_cmp(cmp), m_subsets(subsets), m_quadOrder(quadOrder), m_idString(idString)
445  {};
446  std::string m_cmp;
447  std::string m_subsets;
449  std::string m_idString;
450  };
451 
452 public:
454  {}
455 
457  {}
458 
459  // TODO: replace by call 'func (SmartPtr<G> u, int step, number dt, number t)'
460  bool step_process(SmartPtr<grid_function_type> uNew, /*SmartPtr<grid_function_type> uOld,*/ int step, number time, number dt) OVERRIDE
461  {
462 
463  for (typename std::vector<IntegralSpecs>::iterator it = m_vIntegralData.begin();
464  it!=m_vIntegralData.end(); ++it)
465  {
466  number value=Integral(uNew, it->m_cmp.c_str(), it->m_subsets.c_str(), it->m_quadOrder);
467  UG_LOG("Integral(\t"<< it->m_idString << "\t"<< time << "\t)=\t" << value << std::endl);
468  }
469 
470  return true;
471 
472 
473  }
474 
475 
476  void add_integral_specs(const char* cmp, const char* subsets, int quadOrder, const char* idString)
477  {
478  m_vIntegralData.push_back(IntegralSpecs(cmp, subsets, quadOrder, idString));
479  }
480 
481 protected:
482 
483  std::vector<IntegralSpecs> m_vIntegralData;
484  //const char* cmp,
485  // const char* subsets,
486  // int quadOrder
487 };
488 
489 
490 
492 template<class TDomain, class TAlgebra>
494  : public IOperator< GridFunction<TDomain, TAlgebra> >,
495  public TimeIntegratorSubject<TDomain, TAlgebra>
496 {
497  public:
498 
499  typedef TAlgebra algebra_type;
502 
503 
504  typedef TDomain domain_type;
506 
507  protected:
508  double m_dt;
509  double m_lower_tim;
510  double m_upper_tim;
511 
513 
515 
516 
517  public:
518  // constructor
520  : m_dt(1.0), m_lower_tim(0.0), m_upper_tim(0.0), m_precisionBound(1e-10), m_bNoLogOut(false)
521  {}
522 
524  virtual ~ITimeIntegrator() {};
525 
527 
538  virtual void init(grid_function_type const& u)
539  {
540  // UG_ASSERT(m_spDomainDisc.valid(), "TimeIntegrator<TDomain, TAlgebra>::init: m_spDomainDisc invalid.");
541  // m_spTimeDisc = make_sp(new time_disc_type(m_spDomainDisc, m_theta));
542  }
543 
544 
546 
551  void init()
552  { UG_THROW("Please init with grid function!"); }
553 
555 
556  /* This method is used to prepare the in- and output vector used later in apply.
557  * It can be called, after the init method has been called at least once.
558  * The prepare method is e.g. used to set dirichlet values.
559  */
561 
563 
565  { UG_THROW("Fix interfaces!"); }
566 
568 
570  void set_time_step(double dt)
571  { m_dt = dt; }
572 
573  double get_time_step()
574  { return m_dt; }
575 
576  void set_precision_bound(double precisionBound)
577  { m_precisionBound = precisionBound; return; }
578 
579  void set_no_log_out(bool bNoLogOut)
580  { m_bNoLogOut = bNoLogOut; return; }
581 
582 
583 
584 };
585 
586 
587 
589 template<class TDomain, class TAlgebra>
590 class ILinearTimeIntegrator : public ITimeIntegrator<TDomain, TAlgebra>
591 {
592 
593 public:
598 
599  // forward constructor
601  : base_type() {}
602 
604  : base_type(), m_spLinearSolver(lSolver)
605  {}
606 
607 
609  { m_spLinearSolver=lSolver;}
610 
611 protected:
613 
614 };
615 
616 
618 template<class TAlgebra>
620 {
621 public:
623 
625  m_spTimeDisc(spTimeDisc)
626  {}
627 
629 protected:
631 };
632 
633 
635 template<class TDomain, class TAlgebra>
637  public ILinearTimeIntegrator<TDomain, TAlgebra>,
638  public ITimeDiscDependentObject<TAlgebra>
639 
640 {
641 private:
642 
643 protected:
645 public:
650 
651  // constructor
653  : base_type(), ITimeDiscDependentObject<TAlgebra>(tDisc) {}
654 
656  : base_type(lSolver), ITimeDiscDependentObject<TAlgebra>(tDisc) {}
657 
659 };
660 
661 
662 
663 template<typename TDomain, typename TAlgebra>
665 {
666 
668  UG_COND_THROW(!base_type::m_spLinearSolver.valid(), "Linear solver invalid");
669 
670  // short-cuts
671  GridLevel const &gl = u0->grid_level();
672  time_disc_type &tdisc = *tdisc_dep_type::m_spTimeDisc;
673 
674  // create solution vector & right hand side
676  SmartPtr<typename base_type::vector_type> b= u0->clone_without_values();
677 
678  // solution time series
679  SmartPtr<vector_time_series_type> m_spSolTimeSeries;
680  m_spSolTimeSeries=make_sp(new vector_time_series_type());
681  m_spSolTimeSeries->clear();
682  m_spSolTimeSeries->push(uold, t0);
683 
684  // create matrix operator
685  SmartPtr<typename base_type::assembled_operator_type> spAssOp=make_sp(new typename base_type::assembled_operator_type(tdisc_dep_type::m_spTimeDisc, gl));
686 
687  // integrate
688  if(!base_type::m_bNoLogOut)
689  UG_LOG("+++ Integrating: ["<< t0 <<", "<< t1 <<"]\n");
690 
691  double t = t0;
692  number dt_assembled = -1.0; // invalid
693  int step = 1;
694 
695  number currdt = base_type::m_dt;
696 
697  while((t < t1) && (t1-t > base_type::m_precisionBound))
698  {
699 
700  if(!base_type::m_bNoLogOut)
701  { UG_LOG("+++ Timestep +++" << step << "\n"); }
702 
703  // determine step size
704  number dt = std::min(currdt, t1-t);
705 
706  // prepare step
707  tdisc.prepare_step(m_spSolTimeSeries, dt);
708  if (fabs(dt-dt_assembled) > base_type::m_precisionBound)
709  {
710  // re-assemble operator
711  if(!base_type::m_bNoLogOut)
712  UG_LOG("+++ Reassemble (t=" << t << ", dt=" << dt <<")\n");
713 
714  tdisc.assemble_linear(*spAssOp, *b, gl);
715  (base_type::m_spLinearSolver)->init(spAssOp, *u1);
716  dt_assembled = dt;
717  }
718  else
719  {
720  // keep old operator
721  tdisc.assemble_rhs(*b, gl);
722  }
723 
724  // execute step
725  if (base_type::m_spLinearSolver->apply(*u1, *b))
726  {
727  // ACCEPTING:
728  // push updated solution into time series
729  t += dt;
730  SmartPtr<typename base_type::vector_type> tmp = m_spSolTimeSeries->oldest();
731  VecAssign(*tmp, *u1.template cast_dynamic<typename base_type::vector_type>());
732  m_spSolTimeSeries->push_discard_oldest(tmp, t);
733 
734  this->notify_finalize_step(u1, step++, t+dt, dt);
735  }
736  else
737  {
738  // DISCARDING
739  currdt *= 0.5;
740  }
741 
742  }
743 
744  this->notify_end(u1, step, t1, currdt);
745 
746  return true;
747 };
748 
749 
751 template<class TDomain, class TAlgebra>
753  public ILinearTimeIntegrator<TDomain, TAlgebra>,
754  public ITimeDiscDependentObject<TAlgebra>
755 {
756 protected:
758 public:
764 
765 private:
766 
767 protected:
768 
770 public:
771 
772  // constructor
774  : base_type(), ITimeDiscDependentObject<TAlgebra>(tDisc), m_numSteps(1) {}
775 
777  : base_type(lSolver), ITimeDiscDependentObject<TAlgebra>(tDisc), m_numSteps(1) {}
778 
779  void set_num_steps(int steps) {m_numSteps = steps;}
780 
782 };
783 
784 
785 
786 template<typename TDomain, typename TAlgebra>
788 {
790  UG_COND_THROW(!base_type::m_spLinearSolver.valid(), "Linear solver invalid");
791 
792  // short-cuts
793  GridLevel const &gl = u0->grid_level();
794  time_disc_type &tdisc = *tdisc_dep_type::m_spTimeDisc;
795 
796  // create solution vector & right hand side
798  SmartPtr<typename base_type::vector_type> b= u0->clone_without_values();
799 
800  // solution time series
801  SmartPtr<vector_time_series_type> m_spSolTimeSeries;
802  m_spSolTimeSeries=make_sp(new vector_time_series_type());
803  m_spSolTimeSeries->push(uold, t0);
804 
806 
807  // select number of steps
808  double t = t0;
809  int numSteps = round((t1-t0) / base_type::m_dt);
810  number currdt = (t1-t0) / numSteps;
811 
812 
813 
814  //std::cerr << "+++ Integrating: ["<< t0 <<", "<< t1 <<"] with dt=" << currdt << "("<< numSteps<< " iters)\n";
815  if(!base_type::m_bNoLogOut)
816  {
817  UG_LOG("+++ Integrating: [\t"<< t0 <<"\t, \t"<< t1 <<"\t] with dt=\t" << currdt << "("<< numSteps<< " iters)" << std::endl);
818  }
819 
820  // integrate
821  for(int step = 1; step<=numSteps; ++step)
822  {
823  // determine step size
824  // number dt = std::min(currdt, t1-t);
825  const number dt = currdt;
826 
827  if(!base_type::m_bNoLogOut)
828  {
829  UG_LOG("+++ Const timestep +++" << step<< "(t=" << t << ", dt=" << dt << ")"<< std::endl);
830  }
831  this->notify_init_step(u1, step, t, dt);
832 
833  // prepare step
834  tdisc.prepare_step(m_spSolTimeSeries, dt);
835  if (spAssOp==SPNULL)
836  {
837  // Assemble operator.
838  if(!base_type::m_bNoLogOut) UG_LOG("+++ Assemble (t=" << t << ", dt=" << dt <<")" << std::endl);
839 
840  spAssOp=make_sp(new typename base_type::assembled_operator_type(tdisc_dep_type::m_spTimeDisc, gl));
841  tdisc.assemble_linear(*spAssOp, *b, gl);
842  (base_type::m_spLinearSolver)->init(spAssOp, *u1);
843  }
844  else
845  {
846  // Recycle existing operator.
847  // std::cerr << "Recycling timestep " << step << "\n";
848  tdisc.assemble_rhs(*b, gl);
849  }
850 
851  // execute step
852  if (base_type::m_spLinearSolver->apply(*u1, *b))
853  {
854  // ACCEPTING: push updated solution into time series
855  t += dt;
856  SmartPtr<typename base_type::vector_type> tmp = m_spSolTimeSeries->oldest();
857  VecAssign(*tmp, *u1.template cast_dynamic<typename base_type::vector_type>());
858  m_spSolTimeSeries->push_discard_oldest(tmp, t);
859  this->notify_finalize_step(u1, step, t, dt); // t updated already!
860  }
861  else
862  {
863  UG_THROW("ConstStepLinearTimeIntegrator::apply failed!!!");
864  this->notify_rewind_step(u1, step, t+dt, dt);
865  }
866 
867  }
868 
869  this->notify_end(u1, numSteps, t1, currdt);
870 
871  return true;
872 };
873 
874 
875 
877 template<class TDomain, class TAlgebra>
879  public ILinearTimeIntegrator<TDomain, TAlgebra>,
880  public ITimeDiscDependentObject<TAlgebra>
881 {
882 protected:
884 public:
889 
890 private:
891 
892 
893 protected:
894 
896 
898 
899 
900 public:
902  : base_type(), ITimeDiscDependentObject<TAlgebra>(tDisc1), m_tol(1e-2), m_dtmin(-1.0), m_dtmax(-1.0)
903  {
904  m_spTimeDisc2 = tDisc2;
905  }
906 
907  void init(grid_function_type const& u);
909 
910  void set_tol(double tol) {m_tol = tol;}
913 };
914 
915 
916 template<typename TDomain, typename TAlgebra>
918 {
919  // call base
920  base_type::init(u);
921  //m_spTimeDisc2 = make_sp(new typename base_type::time_disc_type(base_type::m_spDomainDisc, base_type::m_theta));
922 }
923 
924 template<typename TDomain, typename TAlgebra>
926 {
927  // short-cuts
928  GridLevel const &gl = u0->grid_level();
929  time_disc_type &tdisc = *tdisc_dep_type::m_spTimeDisc;
930  time_disc_type &tdisc2 = *m_spTimeDisc2;
931 
932  // create solution vector & right hand side
934  SmartPtr<typename base_type::vector_type> b= u0->clone_without_values();
935 
936  // create matrix operator
937  SmartPtr<typename base_type::assembled_operator_type> spAssOp=make_sp(new typename base_type::assembled_operator_type(tdisc_dep_type::m_spTimeDisc, gl));
938 
939  // create additional solutions
942 
943  // solution time series
944  SmartPtr<vector_time_series_type> m_spSolTimeSeries;
945  m_spSolTimeSeries=make_sp(new vector_time_series_type());
946  m_spSolTimeSeries->push(uold, t0);
947 
948  SmartPtr<vector_time_series_type> m_spSolTimeSeries2;
949  m_spSolTimeSeries2=make_sp(new vector_time_series_type());
950  m_spSolTimeSeries2->push(u2old, t0);
951 
952  // automatic selection of min/max time step
953  if (m_dtmin <= 0.0) { m_dtmin = (t1 - t0)/1.0e+5; }
954  if (m_dtmax <= 0.0) { m_dtmax = (t1 - t0)/10.0; }
955 
956  // Aitken Neville extrapolation
957  const size_t tsteps[2] = {1,2};
958  std::vector<size_t> vsteps (tsteps, tsteps+2);
960 
961  // integrate
962  if(!base_type::m_bNoLogOut)
963  UG_LOG("+++ Integrating: ["<< t0 <<", "<< t1 <<"]\n");
964 
965  double t = t0;
966  int step = 0;
967 
968  number dt = base_type::m_dt;
969  while((t < t1) && (t1-t > base_type::m_precisionBound))
970  {
971  // step: t -> t+dt
972  bool bSuccess = false;
973  while (!bSuccess){
974 
975  // determine step size
976  if (dt<m_dtmin){
977  if(!base_type::m_bNoLogOut)
978  UG_LOG("Step size below minimum")
979  }
980  dt = std::min(dt, t1-t);
981 
982  // basic step
983  if(!base_type::m_bNoLogOut)
984  UG_LOG("+++ Timestep: " << ++step << "\n");
985 
986  tdisc.prepare_step(m_spSolTimeSeries, dt);
987  tdisc.assemble_linear(*spAssOp, *b, gl);
988  base_type::m_spLinearSolver->init(spAssOp, *u1);
989  base_type::m_spLinearSolver->apply(*u1, *b);
990 
991 
992  // control 1/2
993  if(!base_type::m_bNoLogOut)
994  UG_LOG("+++ Control: " << step << "\n");
995 
996  tdisc2.prepare_step(m_spSolTimeSeries, 0.5*dt);
997  tdisc2.assemble_linear(*spAssOp, *b, gl);
998  base_type::m_spLinearSolver->init(spAssOp, *u2);
999  base_type::m_spLinearSolver->apply(*u2, *b);
1000 
1001  // push back solution
1002  SmartPtr<typename base_type::vector_type> tmp2 = m_spSolTimeSeries2->oldest();
1003  (*tmp2)=*(u2.template cast_static<typename base_type::vector_type>());
1004  m_spSolTimeSeries2->push_discard_oldest(tmp2, t + 0.5*dt);
1005 
1006  // control 2/2
1007  tdisc2.prepare_step(m_spSolTimeSeries2, 0.5*dt);
1008  tdisc2.assemble_linear(*spAssOp, *b, gl);
1009  base_type::m_spLinearSolver->init(spAssOp, *u2);
1010  base_type::m_spLinearSolver->apply(*u2, *b);
1011 
1012  // obtain extrapolated solution
1013  timex.set_solution(u1, 0);
1014  timex.set_solution(u2, 1);
1015  timex.apply();
1016 
1017  // predict (subsequent) time step
1018  number eps = timex.get_error_estimates()[0];
1019  number lambda = std::pow(0.8* m_tol/eps, 0.5);
1020 
1021  number dtEst= dt*lambda;
1022  dtEst = std::min(dtEst, 1.5*dt);
1023  dtEst = std::min(dtEst, m_dtmax);
1024 
1025  dt = dtEst;
1026  if (eps <= m_tol)
1027  {
1028  // ACCEPT STEP (and thus solution u2)
1029  if(!base_type::m_bNoLogOut)
1030  UG_LOG("ACCEPTING solution, dtnew=" << dt);
1031 
1032  bSuccess = true;
1033  }
1034  else
1035  {
1036  // DISCARD step
1037  if(!base_type::m_bNoLogOut)
1038  UG_LOG("DISCARDING solution, dtnew=" << dt);
1039 
1040  // => reset solutions
1041  VecAssign(*u1.template cast_dynamic<typename base_type::vector_type>(), *uold);
1042 
1043  // => timeSeries2 has been updated...
1044  SmartPtr<typename base_type::vector_type> tmp = m_spSolTimeSeries2->oldest();
1045  VecAssign(*tmp, *uold);
1046  m_spSolTimeSeries2->push_discard_oldest(tmp, t);
1047 
1048  }
1049 
1050  }
1051 
1052 
1053  // prepare next loop
1054  t += dt;
1055 
1056  // push updated solution into time series (and continue)
1057  SmartPtr<typename base_type::vector_type> tmp = m_spSolTimeSeries->oldest();
1058  VecAssign(*tmp, static_cast<typename base_type::vector_type> (*u2));
1059  m_spSolTimeSeries->push_discard_oldest(tmp, t);
1060 
1061  }
1062 
1063  return true;
1064 };
1065 
1067 {
1068 public:
1070  : m_dtMin(0.0), m_dtMax(std::numeric_limits<double>::max()), m_redFac(1.0), m_incFac(1.0)
1071  {}
1072 
1073  void set_dt_min(double min) { m_dtMin = min; }
1074  double get_dt_min() { return m_dtMin; }
1075 
1076  void set_dt_max(double max) { m_dtMax = max; }
1077  double get_dt_max() { return m_dtMax; }
1078 
1079  void set_reduction_factor(double dec) { m_redFac = dec; }
1080  double get_reduction_factor() { return m_redFac; }
1081 
1082  void set_increase_factor(double inc) { m_incFac = inc; }
1083  double get_increase_factor() { return m_incFac; }
1084 
1085  void rescale(double alpha)
1086  { m_dtMin*= alpha; m_dtMax*= alpha;}
1087 
1088 protected:
1089  double m_dtMin, m_dtMax;
1091 };
1092 
1094 template<class TDomain, class TAlgebra>
1096 : public ITimeIntegrator<TDomain, TAlgebra>
1097 {
1098 public:
1103 
1105  : m_dtBounds() {}
1106 
1108  { m_spSolver=solver;}
1109 
1111  { return m_spSolver;}
1112 
1114  { return m_spSolver;}
1115 
1116  void set_dt_min(double min) { m_dtBounds.set_dt_min(min); }
1117  double get_dt_min() { return m_dtBounds. get_dt_min(); }
1118 
1119  void set_dt_max(double max) { m_dtBounds.set_dt_max(max); }
1120  double get_dt_max() { return m_dtBounds.get_dt_max(); }
1121 
1124 
1127 
1128 protected:
1131 
1132 
1133 };
1134 
1135 
1136 } // ug
1137 
1138 #include "simple_integrator.hpp"
1139 
1140 namespace ug {
1141 
1143 template<class TDomain, class TAlgebra>
1145  public INonlinearTimeIntegrator<TDomain, TAlgebra>
1146 {
1147 public:
1148 
1151 
1153  base_type(), m_wrappedIntegrator(baseIntegrator), m_timePoints() {};
1154 
1156  {
1157  int dstep = 0;
1158  auto tpoint = m_timePoints.begin();
1159  double tcurr = (tpoint == m_timePoints.end()) ? t1 : (*tpoint);
1160  double eps = 1e-8;
1161 
1162  // Perform first step.
1163  //this->notify_init_step(u0, dstep, t0, tcurr-t0);
1164  bool status = m_wrappedIntegrator->apply(u1, tcurr*(1.0-eps), u0, t0);
1165  this->notify_finalize_step(u1, dstep++, tcurr, tcurr-t0);
1166 
1167  // Repeat for all intermediate points.
1168  while (tpoint != m_timePoints.end())
1169  {
1170  tpoint++;
1171  double tnext = (tpoint == m_timePoints.end()) ? t1 : (*tpoint);
1172 
1173  // Perform step.
1174  //this->notify_init_step(u1, dstep, tcurr, tnext-tcurr);
1175  status = status && m_wrappedIntegrator->apply(u1, tnext*(1.0-eps), u1, tcurr);
1176  this->notify_finalize_step(u1, dstep++, tnext, tnext-tcurr);
1177 
1178  tcurr = tnext;
1179  }
1180  this->notify_end(u1, dstep, t1, 0.0);
1181  return status;
1182  }
1183 
1184  void insert_points (std::vector<double> points)
1185  {
1186  m_timePoints = points;
1187  }
1188 protected:
1190  std::vector<double> m_timePoints;
1191 };
1192 
1193 } // namespace ug
1194 
1195 #endif /* TIME_INTEGRATOR_HPP_ */
Definition: time_extrapolation.h:984
void apply(size_t nstages, bool with_error=true)
best error estimate
Definition: time_extrapolation.h:1063
void set_solution(SmartPtr< vector_type > soli, int i)
set solution (for stage i)
Definition: time_extrapolation.h:1014
const std::vector< number > & get_error_estimates() const
Definition: time_extrapolation.h:1028
Sample class for integration observer: Output to VTK.
Definition: time_integrator.hpp:204
ConnectionViewerOutputObserver(const char *filename)
Definition: time_integrator.hpp:209
std::string m_filename
Definition: time_integrator.hpp:228
ITimeIntegratorObserver< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:206
virtual ~ConnectionViewerOutputObserver()
Definition: time_integrator.hpp:215
number m_outputTime
Definition: time_integrator.hpp:229
ConnectionViewerOutputObserver(const char *filename, number t_out)
Definition: time_integrator.hpp:212
GridFunction< TDomain, TAlgebra > grid_function_type
Definition: time_integrator.hpp:207
bool step_process(SmartPtr< grid_function_type > uNew, int step, number time, number dt) OVERRIDE
Definition: time_integrator.hpp:218
Integrate over a given time interval (for a linear problem)
Definition: time_integrator.hpp:755
ConstStepLinearTimeIntegrator(SmartPtr< time_disc_type > tDisc)
Definition: time_integrator.hpp:773
base_type::linear_solver_type linear_solver_type
Definition: time_integrator.hpp:761
base_type::grid_function_type grid_function_type
Definition: time_integrator.hpp:762
void set_num_steps(int steps)
Definition: time_integrator.hpp:779
VectorTimeSeries< typename base_type::vector_type > vector_time_series_type
Definition: time_integrator.hpp:763
bool apply(SmartPtr< grid_function_type > u1, number t1, ConstSmartPtr< grid_function_type > u0, number t0)
Definition: time_integrator.hpp:787
int m_numSteps
Definition: time_integrator.hpp:769
ConstStepLinearTimeIntegrator(SmartPtr< time_disc_type > tDisc, SmartPtr< typename base_type::linear_solver_type > lSolver)
Definition: time_integrator.hpp:776
ILinearTimeIntegrator< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:759
ITimeDiscDependentObject< TAlgebra > tdisc_dep_type
Definition: time_integrator.hpp:757
ITimeDiscretization< TAlgebra > time_disc_type
Definition: time_integrator.hpp:760
This class integrates (t0, t1] with stops at intermediate points tk.
Definition: time_integrator.hpp:1146
DiscontinuityIntegrator(SmartPtr< base_type > baseIntegrator)
Definition: time_integrator.hpp:1152
void insert_points(std::vector< double > points)
Definition: time_integrator.hpp:1184
SmartPtr< base_type > m_wrappedIntegrator
Definition: time_integrator.hpp:1189
base_type::grid_function_type grid_function_type
Definition: time_integrator.hpp:1150
INonlinearTimeIntegrator< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:1149
bool apply(SmartPtr< grid_function_type > u1, number t1, ConstSmartPtr< grid_function_type > u0, number t0)
Definition: time_integrator.hpp:1155
std::vector< double > m_timePoints
Definition: time_integrator.hpp:1190
void assemble_linear(matrix_type &A, vector_type &b)
void assemble_rhs(vector_type &b)
integration of linear systems
Definition: time_integrator.hpp:591
ILinearOperatorInverse< vector_type > linear_solver_type
Definition: time_integrator.hpp:596
void set_linear_solver(SmartPtr< linear_solver_type > lSolver)
Definition: time_integrator.hpp:608
AssembledLinearOperator< TAlgebra > assembled_operator_type
Definition: time_integrator.hpp:597
ILinearTimeIntegrator()
Definition: time_integrator.hpp:600
SmartPtr< linear_solver_type > m_spLinearSolver
Definition: time_integrator.hpp:612
ITimeIntegrator< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:594
ILinearTimeIntegrator(SmartPtr< linear_solver_type > lSolver)
Definition: time_integrator.hpp:603
base_type::vector_type vector_type
Definition: time_integrator.hpp:595
integration of non-linear systems (with bounds on dt)
Definition: time_integrator.hpp:1097
double get_reduction_factor()
Definition: time_integrator.hpp:1123
AssembledOperator< TAlgebra > assembled_operator_type
Definition: time_integrator.hpp:1102
double get_dt_max()
Definition: time_integrator.hpp:1120
void set_dt_min(double min)
Definition: time_integrator.hpp:1116
SmartPtr< solver_type > get_solver()
Definition: time_integrator.hpp:1113
base_type::vector_type vector_type
Definition: time_integrator.hpp:1100
INonlinearTimeIntegrator()
Definition: time_integrator.hpp:1104
ITimeIntegrator< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:1099
TimeStepBounds m_dtBounds
Definition: time_integrator.hpp:1130
IOperatorInverse< vector_type > solver_type
Definition: time_integrator.hpp:1101
SmartPtr< solver_type > m_spSolver
Definition: time_integrator.hpp:1129
ConstSmartPtr< solver_type > get_solver() const
Definition: time_integrator.hpp:1110
double get_increase_factor()
Definition: time_integrator.hpp:1126
void set_increase_factor(double inc)
Definition: time_integrator.hpp:1125
void set_solver(SmartPtr< solver_type > solver)
Definition: time_integrator.hpp:1107
void set_reduction_factor(double dec)
Definition: time_integrator.hpp:1122
double get_dt_min()
Definition: time_integrator.hpp:1117
void set_dt_max(double max)
Definition: time_integrator.hpp:1119
ITimeDiscDependentObject.
Definition: time_integrator.hpp:620
SmartPtr< time_disc_type > m_spTimeDisc
Definition: time_integrator.hpp:630
ITimeDiscretization< TAlgebra > time_disc_type
Definition: time_integrator.hpp:622
ITimeDiscDependentObject(SmartPtr< time_disc_type > spTimeDisc)
Definition: time_integrator.hpp:624
SmartPtr< time_disc_type > get_time_disc()
Definition: time_integrator.hpp:628
virtual void prepare_step(SmartPtr< VectorTimeSeries< vector_type > > prevSol, number dt)=0
Integrates over a given time interval [a,b] with step size dt.
Definition: time_integrator.hpp:496
TDomain domain_type
Definition: time_integrator.hpp:504
void set_precision_bound(double precisionBound)
Definition: time_integrator.hpp:576
void prepare(grid_function_type &u)
prepares functions for application
Definition: time_integrator.hpp:560
void init()
init operator
Definition: time_integrator.hpp:551
double m_lower_tim
Definition: time_integrator.hpp:509
virtual void init(grid_function_type const &u)
init operator depending on a function u
Definition: time_integrator.hpp:538
GridFunction< TDomain, TAlgebra > grid_function_type
Definition: time_integrator.hpp:505
void set_no_log_out(bool bNoLogOut)
Definition: time_integrator.hpp:579
TAlgebra::matrix_type matrix_type
Definition: time_integrator.hpp:501
double m_dt
Definition: time_integrator.hpp:508
virtual ~ITimeIntegrator()
virtual destructor
Definition: time_integrator.hpp:524
bool m_bNoLogOut
Definition: time_integrator.hpp:514
void set_time_step(double dt)
Set initial time step.
Definition: time_integrator.hpp:570
double m_upper_tim
Definition: time_integrator.hpp:510
virtual bool apply(SmartPtr< grid_function_type > u1, number t1, ConstSmartPtr< grid_function_type > u0, number t0)=0
TAlgebra::vector_type vector_type
Definition: time_integrator.hpp:500
double get_time_step()
Definition: time_integrator.hpp:573
TAlgebra algebra_type
Definition: time_integrator.hpp:499
void apply(grid_function_type &u1, const grid_function_type &u0)
Apply operator.
Definition: time_integrator.hpp:564
double m_precisionBound
Definition: time_integrator.hpp:512
ITimeIntegrator()
Definition: time_integrator.hpp:519
Integration observer: Output using Lua callback.
Definition: time_integrator.hpp:435
GridFunction< TDomain, TAlgebra > grid_function_type
Definition: time_integrator.hpp:438
VTKOutput< TDomain::dim > vtk_type
Definition: time_integrator.hpp:439
ITimeIntegratorObserver< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:437
void add_integral_specs(const char *cmp, const char *subsets, int quadOrder, const char *idString)
Definition: time_integrator.hpp:476
virtual ~IntegrationOutputObserver()
Definition: time_integrator.hpp:456
bool step_process(SmartPtr< grid_function_type > uNew, int step, number time, number dt) OVERRIDE
Definition: time_integrator.hpp:460
IntegrationOutputObserver()
Definition: time_integrator.hpp:453
std::vector< IntegralSpecs > m_vIntegralData
Definition: time_integrator.hpp:483
Integrate over a given time interval (for a linear problem)
Definition: time_integrator.hpp:640
bool apply(SmartPtr< grid_function_type > u1, number t1, ConstSmartPtr< grid_function_type > u0, number t0)
Definition: time_integrator.hpp:664
ITimeDiscretization< TAlgebra > time_disc_type
Definition: time_integrator.hpp:647
ITimeDiscDependentObject< TAlgebra > tdisc_dep_type
Definition: time_integrator.hpp:644
base_type::grid_function_type grid_function_type
Definition: time_integrator.hpp:648
LinearTimeIntegrator(SmartPtr< time_disc_type > tDisc, SmartPtr< typename base_type::linear_solver_type > lSolver)
Definition: time_integrator.hpp:655
VectorTimeSeries< typename base_type::vector_type > vector_time_series_type
Definition: time_integrator.hpp:649
ILinearTimeIntegrator< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:646
LinearTimeIntegrator(SmartPtr< time_disc_type > tDisc)
Definition: time_integrator.hpp:652
Definition: time_integrator.hpp:383
SmartPtr< UserData< number, grid_function_type::dim > > m_spReference
Definition: time_integrator.hpp:424
PlotRefOutputObserver(SmartPtr< UserData< number, grid_function_type::dim > > spExactSol)
Definition: time_integrator.hpp:389
virtual ~PlotRefOutputObserver()
Definition: time_integrator.hpp:403
GridFunction< TDomain, TAlgebra > grid_function_type
Definition: time_integrator.hpp:386
VTKOutput< TDomain::dim > vtk_type
Definition: time_integrator.hpp:387
ITimeIntegratorObserver< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:385
SmartPtr< vtk_type > m_sp_vtk
Definition: time_integrator.hpp:425
bool step_process(SmartPtr< grid_function_type > uNew, int step, number time, number dt) OVERRIDE
Definition: time_integrator.hpp:407
Integrate over a given time interval (for a linear problem)
Definition: time_integrator.hpp:881
ITimeDiscDependentObject< TAlgebra > tdisc_dep_type
Definition: time_integrator.hpp:883
VectorTimeSeries< typename base_type::vector_type > vector_time_series_type
Definition: time_integrator.hpp:888
TimeIntegratorLinearAdaptive(SmartPtr< time_disc_type > tDisc1, SmartPtr< time_disc_type > tDisc2)
Definition: time_integrator.hpp:901
double m_tol
Definition: time_integrator.hpp:897
void set_time_step_max(number dt)
Definition: time_integrator.hpp:912
ITimeDiscretization< TAlgebra > time_disc_type
Definition: time_integrator.hpp:886
double m_dtmax
Definition: time_integrator.hpp:897
bool apply(SmartPtr< grid_function_type > u1, number t1, ConstSmartPtr< grid_function_type > u0, number t0)
Definition: time_integrator.hpp:925
void set_tol(double tol)
Definition: time_integrator.hpp:910
void set_time_step_min(number dt)
Definition: time_integrator.hpp:911
double m_dtmin
Definition: time_integrator.hpp:897
SmartPtr< time_disc_type > m_spTimeDisc2
Definition: time_integrator.hpp:895
ILinearTimeIntegrator< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:885
base_type::grid_function_type grid_function_type
Definition: time_integrator.hpp:887
Definition: time_integrator.hpp:1067
void set_increase_factor(double inc)
Definition: time_integrator.hpp:1082
double get_reduction_factor()
Definition: time_integrator.hpp:1080
TimeStepBounds()
Definition: time_integrator.hpp:1069
void set_dt_min(double min)
Definition: time_integrator.hpp:1073
void set_dt_max(double max)
Definition: time_integrator.hpp:1076
double m_incFac
Definition: time_integrator.hpp:1090
double m_redFac
Definition: time_integrator.hpp:1090
double m_dtMin
Definition: time_integrator.hpp:1089
double get_dt_max()
Definition: time_integrator.hpp:1077
void rescale(double alpha)
Definition: time_integrator.hpp:1085
void set_reduction_factor(double dec)
Definition: time_integrator.hpp:1079
double m_dtMax
Definition: time_integrator.hpp:1089
double get_increase_factor()
Definition: time_integrator.hpp:1083
double get_dt_min()
Definition: time_integrator.hpp:1074
Sample class for integration observer: Output to VTK.
Definition: time_integrator.hpp:81
bool step_process(SmartPtr< grid_function_type > uNew, int step, number time, number dt) OVERRIDE
Definition: time_integrator.hpp:106
virtual ~VTKOutputObserver()
Definition: time_integrator.hpp:96
SmartPtr< vtk_type > m_sp_vtk
Definition: time_integrator.hpp:190
GridFunction< TDomain, TAlgebra > grid_function_type
Definition: time_integrator.hpp:84
void set_output_scales(const std::vector< number > &vScales)
Definition: time_integrator.hpp:100
SmartPtr< grid_function_type > m_uOld
Definition: time_integrator.hpp:191
number m_startTime
Definition: time_integrator.hpp:193
VTKOutputObserver(const char *filename, SmartPtr< vtk_type > vtk, number plotStep)
Definition: time_integrator.hpp:93
VTKOutputObserver(const char *filename, SmartPtr< vtk_type > vtk)
Definition: time_integrator.hpp:90
VTKOutputObserver()
Definition: time_integrator.hpp:87
bool end_action(SmartPtr< grid_function_type > u, int step, number time, number dt) OVERRIDE
Definition: time_integrator.hpp:166
std::vector< number > m_vOutputScales
Definition: time_integrator.hpp:195
bool start_action(SmartPtr< grid_function_type > u, int step, number time, number dt) OVERRIDE
Definition: time_integrator.hpp:112
number m_plotStep
Definition: time_integrator.hpp:194
bool finalize_action(SmartPtr< grid_function_type > uNew, int step, number time, number dt) OVERRIDE
Definition: time_integrator.hpp:126
std::string m_filename
Definition: time_integrator.hpp:192
void writeToFile(SmartPtr< grid_function_type > u, int step, number time)
Definition: time_integrator.hpp:177
ITimeIntegratorObserver< TDomain, TAlgebra > base_type
Definition: time_integrator.hpp:83
VTKOutput< TDomain::dim > vtk_type
Definition: time_integrator.hpp:85
ParallelMatrix< SparseMatrix< double > > matrix_type
ParallelVector< Vector< double > > vector_type
void VecAssign(GPUVector< T > &dest, const GPUVector< T > &v1)
value_type & operator()(std::size_t row, std::size_t col)
const NullSmartPtr SPNULL
#define UG_THROW(msg)
#define UG_LOG(msg)
#define UG_COND_THROW(cond, msg)
double number
#define LIMEX_PROFILE_FUNC()
Definition: limex_tools.h:38
struct lua_State lua_State
lua_State * GetDefaultLuaState()
number L2Error(SmartPtr< TGridFunction > spGridFct1, const char *cmp1, SmartPtr< TGridFunction > spGridFct2, const char *cmp2, int quadOrder)
number Integral(number val, SmartPtr< TGridFunction > spGridFct)
void VecScaleAdd(double &dest, double alpha1, const double &v1, double alpha2, const double &v2)
SmartPtr< T, FreePolicy > make_sp(T *inst)
Definition: time_integrator.hpp:442
int m_quadOrder
Definition: time_integrator.hpp:448
std::string m_subsets
Definition: time_integrator.hpp:447
std::string m_cmp
Definition: time_integrator.hpp:445
std::string m_idString
Definition: time_integrator.hpp:449
IntegralSpecs(const char *cmp, const char *subsets, int quadOrder, const char *idString)
Definition: time_integrator.hpp:443
#define OVERRIDE
Definition: time_integrator.hpp:39