diff --git a/changelog-entries/683.md b/changelog-entries/683.md new file mode 100644 index 000000000..6c084d0c1 --- /dev/null +++ b/changelog-entries/683.md @@ -0,0 +1 @@ +- Fixed running without coupling and VTK output in macro-dumux [#683](https://github.com/precice/tutorials/pull/683) diff --git a/two-scale-heat-conduction/macro-dumux/appl/main.cc b/two-scale-heat-conduction/macro-dumux/appl/main.cc index 5932dd22e..5bbd1dd9c 100644 --- a/two-scale-heat-conduction/macro-dumux/appl/main.cc +++ b/two-scale-heat-conduction/macro-dumux/appl/main.cc @@ -108,8 +108,7 @@ int main(int argc, char **argv) // verify that dimensions match const int preciceDim = couplingParticipant.getMeshDimensions(meshName); const int dim = int(leafGridView.dimension); - std::cout << "coupling Dims = " << dim << " , leafgrid dims = " << dim - << std::endl; + std::cout << "Coupling dims = " << preciceDim << " , leafgrid dims = " << dim << std::endl; if (preciceDim != dim) DUNE_THROW(Dune::InvalidStateException, "Dimensions do not match"); } @@ -135,11 +134,15 @@ int main(int argc, char **argv) } } - std::cout << "Number of Coupled Cells:" << coupledElementIdxs.size() - << std::endl; + std::cout << "Number of Coupled Cells:" << coupledElementIdxs.size() << std::endl; + + int numberOfElements; + if (runWithCoupling) { + numberOfElements = coords.size() / couplingParticipant.getMeshDimensions(meshName); + } else { + numberOfElements = coupledElementIdxs.size(); + } - auto numberOfElements = - coords.size() / couplingParticipant.getMeshDimensions(meshName); if (runWithCoupling) { couplingParticipant.setMesh(meshName, coords); @@ -228,14 +231,15 @@ int main(int argc, char **argv) } // time loop parameters - const auto tEnd = getParam("TimeLoop.TEnd"); - double preciceDt = couplingParticipant.getMaxTimeStepSize(); + const auto tEnd = getParam("TimeLoop.TEnd"); + double preciceDt; double solverDt; double dt; if (runWithCoupling) { - solverDt = getParam("TimeLoop.InitialDt"); - dt = std::min(preciceDt, solverDt); + preciceDt = couplingParticipant.getMaxTimeStepSize(); + solverDt = getParam("TimeLoop.InitialDt"); + dt = std::min(preciceDt, solverDt); } else { dt = getParam("TimeLoop.InitialDt"); } @@ -266,7 +270,7 @@ int main(int argc, char **argv) NewtonSolver nonLinearSolver(assembler, linearSolver); // time loop - int n = 0; // counts timesteps for the output interval + int n_out = 0; // counts timesteps for the output interval std::cout << "Time Loop starts" << std::endl; timeLoop->start(); do { @@ -305,7 +309,7 @@ int main(int argc, char **argv) // set new dt as suggested by the Newton solver or by preCICE timeLoop->setTimeStepSize(dt); - std::cout << "Solver starts with target dt: " << dt << std::endl; + std::cout << "nonLinearSolver starts with target dt: " << dt << std::endl; // linearize & solve nonLinearSolver.solve(x, *timeLoop); @@ -319,15 +323,6 @@ int main(int argc, char **argv) timeLoop->reportTimeStep(); xOld = x; - // Vtk output - // TODO: output interval does not work seamlessly when subcycling - n += 1; - if (n == vtkOutputInterval) { - problem->updateVtkOutput(x); - vtkWriter.write(timeLoop->time()); - n = 0; - } - if (runWithCoupling) { int solIdx = 0; for (const auto &element : elements(leafGridView, Dune::Partitions::interior)) { @@ -363,6 +358,25 @@ int main(int argc, char **argv) } } + if (runWithCoupling) { + // if coupling, write VTK output only when time window is complete + if (couplingParticipant.isTimeWindowComplete()) { + n_out += 1; + if (n_out == vtkOutputInterval) { + problem->updateVtkOutput(x); + vtkWriter.write(timeLoop->time()); + n_out = 0; + } + } + } else { + n_out += 1; + if (n_out == vtkOutputInterval) { + problem->updateVtkOutput(x); + vtkWriter.write(timeLoop->time()); + n_out = 0; + } + } + std::cout << "Time: " << timeLoop->time() << std::endl; } while (!timeLoop->finished());