[macro-dumux] Cleanup time loop

Author: mathiskelm

two-scale-heat-conduction/macro-dumux/appl/main.cc

@@ -262,28 +262,29 @@ int main(int argc, char **argv)
 
   // time loop
   int n = 0; // counts timesteps for the output interval
+  const auto runWithCoupling = getParam<bool>("Precice.RunWithCoupling");
   std::cout << "Time Loop starts" << std::endl;
   timeLoop->start();
   do {
-    if (getParam<bool>("Precice.RunWithCoupling") == true) {
+    if (runWithCoupling) {
       if (couplingParticipant.isCouplingOngoing() == false)
         break;
 
       // write checkpoint
       if (couplingParticipant.requiresToWriteCheckpoint()) {
-        // TODO
-        //xOld = x;
         xCheckpoint = x;
         timeCheckpoint = timeLoop->time();
         timeStepCheckpoint = timeLoop->timeStepIndex();
       }
 
       preciceDt = couplingParticipant.getMaxTimeStepSize();
-      solverDt  = std::min(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize()), timeLoop->maxTimeStepSize());
+      solverDt  = std::min(nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize()),
+                           timeLoop->maxTimeStepSize());
       dt        = std::min(preciceDt, solverDt);
-      timeLoop->setTimeStepSize(dt);
 
       // read porosity and conductivity data from other solver
+      // TODO: data needs to be updated if Newton solver adapts time-step size
+      // and coupling data is interpolated in time
       couplingParticipant.readQuantityFromOtherSolver(meshName, readDatak00,
                                                       dt);
       couplingParticipant.readQuantityFromOtherSolver(meshName, readDatak01,
@@ -294,95 +295,71 @@ int main(int argc, char **argv)
                                                       dt);
       couplingParticipant.readQuantityFromOtherSolver(meshName,
                                                       readDataPorosity, dt);
-      // store coupling data in problem
+      // store coupling data in spatial params, exchange with MPI
       problem->spatialParams().updateCouplingData();
     } else {
       dt = std::min(
           nonLinearSolver.suggestTimeStepSize(timeLoop->timeStepSize()),
-          getParam<Scalar>("TimeLoop.MaxDt"));
-      timeLoop->setTimeStepSize(dt);
+          timeLoop->maxTimeStepSize());
     }
+    // set new dt as suggested by the newton solver or by preCICE
+    timeLoop->setTimeStepSize(dt);
 
-    std::cout << "Solver starts" << std::endl;
+    std::cout << "Solver starts with Target dt: " << dt  << std::endl;
 
-    std::cout << "Using dt: " << dt << std::endl;
     // linearize & solve
     nonLinearSolver.solve(x, *timeLoop);
     // save actual time-step size
     dt = timeLoop->timeStepSize();
+
     // DuMux advance + report
     gridVariables->advanceTimeStep();
     timeLoop->advanceTimeStep();
     timeLoop->reportTimeStep();
+    xOld = x;
+
+    // Vtk output
+    // TODO: output interval doesn't work easily with subcycling
+    n += 1;
+    if (n == vtkOutputInterval) {
+      problem->updateVtkOutput(x);
+      vtkWriter.write(timeLoop->time());
+      n = 0;
+    }
 
-    for (int solIdx = 0; solIdx < numberOfElements; ++solIdx)
-      temperatures[solIdx] = x[solIdx][problem->returnTemperatureIdx()];
+    if (runWithCoupling) {
+      // write coupling data to preCICE
+      for (int solIdx = 0; solIdx < numberOfElements; ++solIdx)
+        temperatures[solIdx] = x[solIdx][problem->returnTemperatureIdx()];
 
-    if (getParam<bool>("Precice.RunWithCoupling") == true) {
       couplingParticipant.writeQuantityVector(meshName,
                                               writeDataConcentration, temperatures);
       couplingParticipant.writeQuantityToOtherSolver(meshName,
                                                      writeDataConcentration);
-    }
-
-    // advance precice
-    if (getParam<bool>("Precice.RunWithCoupling") == true) {
 
+      // advance precice
       if ((!fabs(preciceDt - dt)) < 1e-14) {
         std::cout << "dt from preCICE is different than dt from DuMuX."
                   << " preCICE dt = " << preciceDt
                   << " and DuMuX dt =" << solverDt
                   << " resulted in dt = " << dt
                   << std::endl;
       }
+      std::flush(std::cout);
       couplingParticipant.advance(dt);
-    }
 
-    if (getParam<bool>("Precice.RunWithCoupling") == true) {
+      // reset to checkpoint if not converged
       if (couplingParticipant.requiresToReadCheckpoint()) {
-        // make the new solution the old solution
-        // TODO
-        //x = xOld;
         x = xCheckpoint;
         xOld = x;
         timeLoop->setTime(timeCheckpoint, timeStepCheckpoint);
+        // TODO: previousTimeStep might be more appropriate, last one could be small
         timeLoop->setTimeStepSize(dt);
         gridVariables->update(x);
         gridVariables->advanceTimeStep();
-      } else // coupling successful OR not at time window!!
-      {
-        // TODO
-        xOld = x;
-        n += 1;
-        if (n == vtkOutputInterval) {
-          problem->updateVtkOutput(x);
-          vtkWriter.write(timeLoop->time());
-          n = 0;
-        }
-        //gridVariables->advanceTimeStep();
-        //// advance the time loop to the next step
-        //timeLoop->advanceTimeStep();
-        //// report statistics of this time step
-        //timeLoop->reportTimeStep();
-      }
-    } else {
-      xOld = x;
-      //gridVariables->advanceTimeStep();
-      //// advance the time loop to the next step
-      //timeLoop->advanceTimeStep();
-      //// report statistics of this time step
-      //timeLoop->reportTimeStep();
-
-      // output every outputinterval steps
-      n += 1;
-      if (n == vtkOutputInterval) {
-        problem->updateVtkOutput(x);
-        vtkWriter.write(timeLoop->time());
-        n = 0;
+        continue;
       }
     }
-    // set new dt as suggested by the newton solver or by preCICE
-    timeLoop->setTimeStepSize(dt);
 
     std::cout << "Time: " << timeLoop->time() << std::endl;
 
@@ -393,7 +370,7 @@ int main(int argc, char **argv)
   ////////////////////////////////////////////////////////////
   // finalize, print dumux message to say goodbye
   ////////////////////////////////////////////////////////////
-  if (getParam<bool>("Precice.RunWithCoupling") == true) {
+  if (runWithCoupling) {
     couplingParticipant.finalize();
   }
   // print dumux end message