fix a silly initialization bug

src/language/algorithms/ParabolicHeat.cpp

@@ -235,8 +235,9 @@ ParabolicHeatScheme<Dimension>::ParabolicHeatScheme(
 
     CellValue<double> Tj =
       InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(T_id, mesh_data.xj());
-    CellValue<double> Temperature{mesh->connectivity()};
-    Temperature.fill(1);
+    CellValue<double> Temperature =
+      InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(T_id, mesh_data.xj());
+    //{mesh->connectivity()};
     FaceValue<double> Tl =
       InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::face>(T_id, mesh_data.xl());
     FaceValue<double> Temperature_face =
@@ -254,7 +255,14 @@ ParabolicHeatScheme<Dimension>::ParabolicHeatScheme(
     }
     double time = 0;
     Assert(dt > 0, "The time-step have to be positive!");
+    const CellValue<const double> primal_Vj = mesh_data.Vj();
     do {
+      {
+        VTKWriter vtk_writer("Temperature_" + std::to_string(Dimension) + std::to_string(time), 0.01);
+
+        vtk_writer.write(mesh, {NamedItemValue{"T", Temperature}, NamedItemValue{"Vj", primal_Vj}}, 0,
+                         true);   // forces last output
+      }
       double deltat = std::min(dt, Tf - time);
       std::cout << "Current time = " << time << " time-step = " << deltat << " final time = " << Tf << "\n";
       ScalarDiamondScheme<Dimension>(i_mesh, bc_descriptor_list, kappa_id, f_id, Temperature, Temperature_face, Tf,
@@ -267,7 +275,6 @@ ParabolicHeatScheme<Dimension>::ParabolicHeatScheme(
       Vector<double> face_error{mesh->numberOfFaces()};
       double error_max                    = 0.;
       size_t cell_max                     = 0;
-      const CellValue<const double> primal_Vj = mesh_data.Vj();
       const FaceValue<const double> mes_l = [&] {
         if constexpr (Dimension == 1) {
           FaceValue<double> compute_mes_l{mesh->connectivity()};