From 8a7dc587a50805584bd6d7e41c46fabf5bdc9c2c Mon Sep 17 00:00:00 2001
From: Daniel Weindl <daniel.weindl@uni-bonn.de>
Date: Mon, 21 Oct 2024 16:25:32 +0200
Subject: [PATCH] Run simulations with json I/O

Adds an executable to evaluate the objective function and its gradient based on parameters provided in a JSON file.
The output will be written to another JSON file. This makes it easier to compare specific simulation results to some other tool.
---
 templates/CMakeLists.template.txt |  27 +++++++
 templates/CMakeLists.txt          |   2 +-
 templates/main_json.cpp           | 122 ++++++++++++++++++++++++++++++
 3 files changed, 150 insertions(+), 1 deletion(-)
 create mode 100644 templates/main_json.cpp

diff --git a/templates/CMakeLists.template.txt b/templates/CMakeLists.template.txt
index d842a809..58cf9e19 100644
--- a/templates/CMakeLists.template.txt
+++ b/templates/CMakeLists.template.txt
@@ -99,6 +99,33 @@ target_link_libraries(${PROJECT_NAME}
 )
 target_include_directories(${PROJECT_NAME} PUBLIC "model")
 
+#################################
+# simulateJSON executable
+#################################
+
+option(BUILD_JSON_SIMULATOR "Build JSON simulator" ON)
+
+if(BUILD_JSON_SIMULATOR)
+    project(simulateJSON_${MODEL_NAME})
+
+    set(SRC_LIST main_json.cpp)
+
+    add_executable(${PROJECT_NAME} ${SRC_LIST})
+    add_dependencies(${PROJECT_NAME} ${MODEL_NAME})
+
+    include(FetchContent)
+
+    FetchContent_Declare(json URL https://github.com/nlohmann/json/releases/download/v3.11.3/json.tar.xz)
+    FetchContent_MakeAvailable(json)
+
+    target_link_libraries(
+        ${PROJECT_NAME}
+        PRIVATE
+        model
+        Upstream::parpe
+        nlohmann_json::nlohmann_json
+    )
+endif()
 #################################
 
 project(${MODEL_NAME}) # for IDE, will use last project as label
diff --git a/templates/CMakeLists.txt b/templates/CMakeLists.txt
index 0eca0608..70069c4d 100644
--- a/templates/CMakeLists.txt
+++ b/templates/CMakeLists.txt
@@ -1,7 +1,7 @@
 # To show in IDE
 add_custom_target(
   templates SOURCES CMakeLists.template.txt main.cpp main_debug.cpp
-                    main_nominal.cpp main_simulate.cpp)
+                    main_nominal.cpp main_simulate.cpp main_json.cpp)
 
 set_target_properties(
   templates
diff --git a/templates/main_json.cpp b/templates/main_json.cpp
new file mode 100644
index 00000000..aeca58b4
--- /dev/null
+++ b/templates/main_json.cpp
@@ -0,0 +1,122 @@
+/*
+ * Executable for evaluating the objective function for parameters specified
+ * in a JSON file. The results will be written to a JSON file.
+ *
+ * The input is a list of objects with the following fields:
+ * - data_file: the name of the file containing the conditions
+ * - parameters: object with parameterId => parameterValue
+ * - gradient: boolean indicating whether the gradient should be computed
+ *
+ * The output is a list of objects with the following fields:
+ * - data_file: the name of the file containing the conditions
+ * - parameters: object with parameterId => parameterValue
+ * - fval: the objective function value
+ * - gradient: the gradient (object with parameterId => gradient entry),
+ *     or `false` if not computed
+ */
+
+#include <parpecommon/parpeConfig.h>
+
+#include <parpeamici/multiConditionDataProvider.h>
+#include <parpeamici/multiConditionProblem.h>
+#include <parpeamici/standaloneSimulator.h>
+#include <parpecommon/misc.h>
+#include <parpecommon/parpeConfig.h>
+
+#include <cstdio> // remove
+#include <iostream>
+#include <fstream>
+#include <string_view>
+#include <nlohmann/json.hpp>
+#ifdef PARPE_ENABLE_MPI
+#include <mpi.h>
+#endif
+
+using json = nlohmann::json;
+
+// fields in the input/output json
+constexpr std::string_view DATA_FILE = "data_file";
+constexpr std::string_view PARAMETERS = "parameters";
+constexpr std::string_view FVAL = "fval";
+constexpr std::string_view GRADIENT = "gradient";
+constexpr std::string_view STATUS = "status";
+constexpr std::string_view ERROR = "error";
+
+
+namespace amici::generic_model {
+std::unique_ptr<amici::Model> getModel();
+}
+
+
+void printUsage() {
+    std::cerr<<"Error: wrong number of arguments.\n";
+    std::cerr<<"Usage: ... infile outfile\n";
+}
+
+int main(int argc, char **argv) {
+    int status = EXIT_SUCCESS;
+
+    if(argc != 3) {
+        printUsage();
+        return EXIT_FAILURE;
+    }
+
+    std::string inFileArgument = argv[1];
+    std::string outFileArgument = argv[2];
+
+    std::ifstream istream(inFileArgument);
+    json jobs = json::parse(istream);
+
+    json results = json::array();
+    std::unique_ptr<parpe::MultiConditionDataProvider> dataProvider;
+    std::string h5file_name;
+
+    for (auto& job : jobs) {
+        std::cout<<job<<std::endl;
+        auto cur_result = job;
+        try {
+            if (job[DATA_FILE] != h5file_name) {
+                h5file_name = job[DATA_FILE];
+                dataProvider = std::make_unique<parpe::MultiConditionDataProviderHDF5>(
+                    amici::generic_model::getModel(), h5file_name);
+            }
+            parpe::MultiConditionProblem problem(
+                dataProvider.get(), nullptr,
+                std::make_unique<parpe::Logger>(),
+                nullptr);
+
+                   // parameter values
+            auto objective = dynamic_cast<parpe::SummedGradientFunctionGradientFunctionAdapter<int>*>(problem.cost_fun_.get());
+
+            auto parameter_ids = objective->getParameterIds();
+            std::vector<double> parameters;
+            parameters.reserve(parameter_ids.size());
+
+            for (auto& parameter_id: parameter_ids) {
+                parameters.push_back(job[PARAMETERS][parameter_id]);
+            }
+
+            double fval;
+            std::vector<double> gradient(parameter_ids.size());
+            auto fv = objective->evaluate(parameters, fval, job[GRADIENT]?gsl::make_span(gradient):gsl::span<double>(nullptr, 0));
+            cur_result[FVAL] = fval;
+            if (job[GRADIENT]) {
+                cur_result[GRADIENT] = json::object();
+                for (size_t i = 0; i < parameter_ids.size(); ++i) {
+                    cur_result[GRADIENT][parameter_ids[i]] = gradient[i];
+                }
+            }
+            cur_result[STATUS] = fv;
+        } catch(std::exception& e) {
+            std::cerr<<e.what()<<std::endl;
+            cur_result[ERROR] = e.what();
+
+        }
+        results.push_back(cur_result);
+        // not efficient, but let's save it after every simulation for now
+        std::ofstream ostream(outFileArgument);
+        ostream << std::setw(4) << results << std::endl;
+    }
+
+    return status;
+}