Peter Bastian committed Mar 17, 2006 1 2 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- // vi: set et ts=4 sw=2 sts=2:  Christoph Grüninger committed Oct 24, 2014 3 #include // know what grids are present  Peter Bastian committed Mar 17, 2006 4 5 6 #include // for input/output to shell #include // for input/output to files #include // STL vector class  Robert Klöfkorn committed Jul 29, 2007 7   Christian Engwer committed May 16, 2006 8 #include // mapper class  Christian Engwer committed Nov 06, 2012 9 #include // include mpi helper class  Peter Bastian committed Mar 17, 2006 10   Robert Klöfkorn committed Oct 25, 2006 11 #include "vtkout.hh"  Peter Bastian committed Mar 21, 2006 12 #include "transportproblem2.hh"  Peter Bastian committed Mar 17, 2006 13 14 15 16 17 18 19 20 21 22 23 24 #include "initialize.hh" #include "evolve.hh" #include "finitevolumeadapt.hh" //=============================================================== // the time loop function working for all types of grids //=============================================================== template void timeloop (G& grid, double tend, int lmin, int lmax) { // make a mapper for codim 0 entities in the leaf grid  Ansgar Burchardt committed Aug 09, 2017 25 26  Dune::LeafMultipleCodimMultipleGeomTypeMapper mapper(grid, Dune::mcmgElementLayout());  Peter Bastian committed Mar 17, 2006 27 28 29 30 31 32 33 34  // allocate a vector for the concentration std::vector c(mapper.size()); // initialize concentration with initial values initialize(grid,mapper,c); for (int i=grid.maxLevel(); i=lmax) break;  Christian Engwer committed Mar 09, 2009 36  finitevolumeadapt(grid,mapper,c,lmin,lmax,0); /*@\label{afv:in}@*/  Peter Bastian committed Mar 17, 2006 37 38  initialize(grid,mapper,c); }  Robert Klöfkorn committed Jul 29, 2007 39 40  // write initial data  Christian Engwer committed Mar 09, 2009 41  vtkout(grid,c,"concentration",0,0);  Peter Bastian committed Mar 17, 2006 42   Robert Klöfkorn committed Jul 29, 2007 43  // variables for time, timestep etc.  Peter Bastian committed Mar 17, 2006 44  double dt, t=0;  Robert Klöfkorn committed Jul 29, 2007 45  double saveStep = 0.1;  Andreas Dedner committed Nov 26, 2011 46 47  const double saveInterval = t + 0.1; int counter = 1;  Robert Klöfkorn committed Jul 29, 2007 48 49 50  int k = 0; std::cout << "s=" << grid.size(0) << " k=" << k << " t=" << t << std::endl;  Peter Bastian committed Mar 17, 2006 51 52  while (t= saveStep) { // write data  Christian Engwer committed Mar 09, 2009 66  vtkout(grid,c,"concentration",counter,t);  Robert Klöfkorn committed Jul 29, 2007 67 68 69 70 71 72 73  // increase counter and saveStep for next interval saveStep += saveInterval; ++counter; } // print info about time, timestep size and counter  Christian Engwer committed Mar 09, 2009 74 75  std::cout << "s=" << grid.size(0) << " k=" << k << " t=" << t << " dt=" << dt << std::endl;  Robert Klöfkorn committed Jul 29, 2007 76 77  // for unstructured grids call adaptation algorithm  Oliver Sander committed Aug 26, 2010 78  finitevolumeadapt(grid,mapper,c,lmin,lmax,k); /*@\label{afv:ad}@*/  Peter Bastian committed Mar 17, 2006 79  }  Robert Klöfkorn committed Jul 29, 2007 80 81  // write last time step  Christian Engwer committed Mar 09, 2009 82 83 84  vtkout(grid,c,"concentration",counter,tend); // write  Peter Bastian committed Mar 17, 2006 85 86 87 88 89 90 91 92 } //=============================================================== // The main function creates objects and does the time loop //=============================================================== int main (int argc , char ** argv) {  Robert Klöfkorn committed Jul 29, 2007 93 94  // initialize MPI, finalize is done automatically on exit Dune::MPIHelper::instance(argc,argv);  Peter Bastian committed Mar 17, 2006 95   Christian Engwer committed Dec 05, 2006 96 97  // start try/catch block to get error messages from dune try {  Robert Klöfkorn committed Jul 29, 2007 98 99  using namespace Dune;  Andreas Dedner committed Jan 14, 2010 100 101 102 103  // the GridSelector :: GridType is defined in gridtype.hh and is // set during compilation typedef GridSelector :: GridType Grid;  Robert Klöfkorn committed Jul 29, 2007 104 105  // use unitcube from grids std::stringstream dgfFileName;  Christoph Grüninger committed May 10, 2013 106 107  dgfFileName << DUNE_GRID_HOWTO_EXAMPLE_GRIDS_PATH << "unitcube" << Grid::dimension << ".dgf";  Robert Klöfkorn committed Jul 29, 2007 108   Andreas Dedner committed Jan 14, 2010 109 110  // create grid pointer GridPtr gridPtr( dgfFileName.str() );  Robert Klöfkorn committed Jul 29, 2007 111 112  // grid reference  Andreas Dedner committed Jan 14, 2010 113  Grid& grid = *gridPtr;  Robert Klöfkorn committed Jul 29, 2007 114 115  // minimal allowed level during refinement  Andreas Dedner committed Jan 14, 2010 116  int minLevel = 2 * DGFGridInfo::refineStepsForHalf();  Robert Klöfkorn committed Jul 29, 2007 117 118 119 120 121  // refine grid until upper limit of level grid.globalRefine(minLevel); // maximal allowed level during refinement  Andreas Dedner committed Jan 14, 2010 122  int maxLevel = minLevel + 3 * DGFGridInfo::refineStepsForHalf();  Robert Klöfkorn committed Jul 29, 2007 123 124 125  // do time loop until end time 0.5 timeloop(grid, 0.5, minLevel, maxLevel);  Christian Engwer committed Dec 05, 2006 126 127  } catch (std::exception & e) {  Christoph Grüninger committed Mar 08, 2016 128  std::cout << "ERROR: " << e.what() << std::endl;  Christian Engwer committed Dec 05, 2006 129 130 131 132 133 134  return 1; } catch (...) { std::cout << "Unknown ERROR" << std::endl; return 1; }  Peter Bastian committed Mar 17, 2006 135 136 137 138  // done return 0; }