From 1129dc36f21b682b37f1c363c1f55b6fdd6993ac Mon Sep 17 00:00:00 2001
From: Stefanie Guenther <guenther5@llnl.gov>
Date: Mon, 7 Apr 2025 13:33:24 -0700
Subject: [PATCH] Bugfix for petsc-parallel output.

Expected energy and population need to be computed outside of the
if-statement that check whether the file is open on this processor.
---
 src/output.cpp | 13 +++++++------
 1 file changed, 7 insertions(+), 6 deletions(-)

diff --git a/src/output.cpp b/src/output.cpp
index 746faf48..641bb8f8 100644
--- a/src/output.cpp
+++ b/src/output.cpp
@@ -201,24 +201,25 @@ void Output::writeDataFiles(int timestep, double time, const Vec state, MasterEq
   /* Write output only every <num> time-steps */
   if (timestep % output_frequency == 0) {
 
+
     /* Write expected energy levels to file */
     for (size_t iosc = 0; iosc < expectedfile.size(); iosc++) {
+      double expected = mastereq->getOscillator(iosc)->expectedEnergy(state);
       if (expectedfile[iosc] != NULL) {
-        double expected = mastereq->getOscillator(iosc)->expectedEnergy(state);
         fprintf(expectedfile[iosc], "%.8f %1.14e\n", time, expected);
       }
     }
 
+    double expected_comp = mastereq->expectedEnergy(state);
     if (expectedfile_comp != NULL) {
-      double expected_comp = mastereq->expectedEnergy(state);
       fprintf(expectedfile_comp, "%.8f %1.14e\n", time, expected_comp);
     }
 
     /* Write population to file */
     for (size_t iosc = 0; iosc < populationfile.size(); iosc++) {
+      std::vector<double> pop (mastereq->getOscillator(iosc)->getNLevels(), 0.0);
+      mastereq->getOscillator(iosc)->population(state, pop);
       if (populationfile[iosc] != NULL) {
-        std::vector<double> pop (mastereq->getOscillator(iosc)->getNLevels(), 0.0);
-        mastereq->getOscillator(iosc)->population(state, pop);
         fprintf(populationfile[iosc], "%.8f ", time);
         for (size_t i = 0; i<pop.size(); i++) {
           fprintf(populationfile[iosc], " %1.14e", pop[i]);
@@ -227,9 +228,9 @@ void Output::writeDataFiles(int timestep, double time, const Vec state, MasterEq
       }
     }
 
+    std::vector<double> population_comp; 
+    mastereq->population(state, population_comp);
     if (populationfile_comp != NULL) {
-      std::vector<double> population_comp; 
-      mastereq->population(state, population_comp);
       fprintf(populationfile_comp, "%.8f  ", time);
       for (size_t i=0; i<population_comp.size(); i++){
         fprintf(populationfile_comp, "%1.14e  ", population_comp[i]);