diff --git a/examples/plain/copasi.d.ts b/examples/plain/copasi.d.ts new file mode 100644 index 0000000..f3539e4 --- /dev/null +++ b/examples/plain/copasi.d.ts @@ -0,0 +1,99 @@ +import type exp from "constants"; + +export interface SimResult { + num_variables: number; + recorded_steps: number; + titles: string[]; + columns: string[]; +} + +export interface SpeciesInfo { + compartment: string; + concentration: number; + id: string; + name: string; + initial_concentration: number; + initial_particle_number: number; + particle_number: number; + type: string; +} + +export interface CompartmentInfo { + id: string; + name: string; + size: number; + type: string; +} +export interface LocalParamterInfo { + name: string; + value: number; +} + +export interface ReactionInfo { + id: string; + name: string; + reversible: boolean; + scheme: string; + local_paramters: LocalParamterInfo[]; +} + +export interface ModelName +{ + name: string; + notes: string; +} + +export interface GlobalParamterInfo { + name: string; + value: number; + initial_value: number; + id: string; + type: string; +} + +export interface ModelInfo { + species: SpeciesInfo[]; + compartments: CompartmentInfo[]; + reactions: ReactionInfo[]; + global_parameters: GlobalParamterInfo[]; + model : ModelName; + time: number; + status: string; + messages: string; +} + +export default class COPASI { + constructor(module: any); + reset(): void; + readonly version: string; + _vectorToArray(v: any): any[]; + loadExample(path: string) : ModelInfo; + loadModel(modelCode: string): ModelInfo; + simulate() : object; + simulate2D() : number[][]; + simulateEx(startTime : number, endTime : number, numPoints : number) : SimResult; + simulateEx2D(startTime : number, endTime : number, numPoints : number) : number[][]; + simulateYaml(yamlProcessingOptions : string|object) : SimResult; + simulateYaml2D(yamlProcessingOptions : string|object) : SimResult; + readonly floatingSpeciesConcentrations : number[]; + readonly ratesOfChange : number[]; + readonly floatingSpeciesNames : string[]; + readonly floatingSpeciesIds : string[]; + readonly boundarySpeciesConcentrations : number[]; + readonly boundarySpeciesNames : string[]; + readonly boundarySpeciesIds : string[]; + readonly reactionNames : string[]; + readonly reactionIds : string[]; + readonly reactionRates : number[]; + readonly compartmentNames : string[]; + readonly compartmentIds : string[]; + readonly compartmentSizes : number[]; + readonly globalParameterNames : string[]; + readonly globalParameterIds : string[]; + readonly globalParameterValues : number[]; + readonly localParameterNames : string[]; + readonly localParameterValues : number[]; +} + +export {COPASI}; +export default COPASI; \ No newline at end of file diff --git a/src/copasijs.cpp b/src/copasijs.cpp index 6874571..3a3b69b 100755 --- a/src/copasijs.cpp +++ b/src/copasijs.cpp @@ -20,7 +20,7 @@ using namespace nlohmann; // define a struct for information from model, containing ptr to CDataObject, SBML id string, and common name struct CModelElement { - const CDataObject *pObj; + CDataObject *pObj; std::string sbmlId; CRegisteredCommonName cnInitial; CRegisteredCommonName cn; @@ -310,8 +310,75 @@ double getValue(const std::string &nameOrId) return std::numeric_limits::quiet_NaN(); } +bool setModelElement(std::map& map, std::map& idMap, const std::string& name, double value, CModel* model) +{ + auto it = map.find(name); + if (it == map.end()) + { + auto idIt = idMap.find(name); + if (idIt != idMap.end()) + it = map.find(idIt->second); + } + if (it != map.end()) + { + CMetab* pMetab = dynamic_cast(it->second.pObj); + if (pMetab) + { + model->updateInitialValues(pMetab->getInitialConcentrationReference(), false); + + double oldValue = model->getMathContainer().getMathObject(pMetab->getInitialConcentrationReference())->getValue(); + + pMetab->setConcentration(value); + pMetab->setInitialConcentration(value); + model->updateInitialValues(pMetab->getInitialConcentrationReference(), false); + + double newValue = model->getMathContainer().getMathObject(pMetab->getInitialConcentrationReference())->getValue(); + + + return true; + } + + CModelEntity* pEntity = dynamic_cast(it->second.pObj); + if (pEntity) + { + pEntity->setValue(value); + pEntity->setInitialValue(value); + std::set changes = {pEntity->getInitialValueReference(), pEntity->getValueReference()}; + + model->updateInitialValues(changes, true); + return true; + } + + } + return false; +} + void setValue(const std::string& nameOrId, double value) { + if (!pDataModel) + return; + + auto* model = pDataModel->getModel(); + if (!model) + return; + + // look through floating species by name first, falling back + // for id if not found + if (setModelElement(mFloatingSpecies, mFloatingSpeciesIdMap, nameOrId, value, model)) + return; + + // boundary species + if (setModelElement(mBoundarySpecies, mBoundarySpeciesIdMap, nameOrId, value, model)) + return; + + // now compartments + if (setModelElement(mCompartments, mCompartmentsIdMap, nameOrId, value, model)) + return; + + // and global parameters + if (setModelElement(mGlobalParameters, mGlobalParametersIdMap, nameOrId, value, model)) + return; + } @@ -628,7 +695,7 @@ ordered_json buildModelInfo() continue; CModelElement lp = { - obj[0], + const_cast(obj[0]), "", cns[0], cns[0], @@ -930,7 +997,7 @@ std::string simulateJSON(ordered_json& yaml) if (!task.initialize(CCopasiTask::OUTPUT_UI, pDataModel, NULL)) return strdup("Error initializing task"); - if (!task.process(false)) + if (!task.process(true)) return strdup("Error processing task"); if (!task.restore()) diff --git a/test/copasijs_test.cpp b/test/copasijs_test.cpp index 449178b..7207c6e 100644 --- a/test/copasijs_test.cpp +++ b/test/copasijs_test.cpp @@ -100,6 +100,104 @@ TEST_CASE("Load Model", "[copasijs]") } +TEST_CASE("Load SBML Model", "[copasijs][sbml]") +{ + Instance instance; + std::string model = loadFromFile(getTestFile("../example_files/oscli.xml")); + REQUIRE(!model.empty()); + REQUIRE(model != "Error loading model"); + + auto selectionList = getSelectionList(); + REQUIRE_THAT(selectionList, Catch::Matchers::Equals(std::vector{"Time", "S1", "S2"})); + auto selectionValues = getSelectionValues(); + REQUIRE(selectionValues.size() == 3); + REQUIRE_THAT(selectionValues, Catch::Matchers::Approx(std::vector{0, 0, 1})); + + // lets have a look at the time course settings + auto timeCourseSettings = getTimeCourseSettings(); + CAPTURE(timeCourseSettings); + + // now simulate using onestep + double time = 0; + double stepSize = 1; + double numSteps = 11; + std::vector timePoints; + std::vector xPoints; + std::vector yPoints; + while (time < numSteps) { + oneStep(time, stepSize); + timePoints.push_back(time); + auto values = getSelectionValues(); + xPoints.push_back(values[1]); + yPoints.push_back(values[2]); + time += stepSize; + } + REQUIRE(timePoints.size() == 11); + REQUIRE_THAT(timePoints, Catch::Matchers::Approx(std::vector{0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10})); + + reset(); + + // ensure that time is again 0 + selectionValues = getSelectionValues(); + REQUIRE(selectionValues.size() == 3); + REQUIRE_THAT(selectionValues, Catch::Matchers::Approx(std::vector{0, 0, 1})); + + // now run a simulation to look that the json format is in the expected format + auto data = simulateEx(0, 10, 11); + CAPTURE(data); + auto json = nlohmann::json::parse(data); + REQUIRE(json["titles"][0] == "Time"); + REQUIRE(json["columns"][0].size() == json["recorded_steps"].get()); + + reset(); + // now use setValue / getValue to change the value of S1 + double value = getValue("S1"); + CAPTURE(value); + REQUIRE(value == 0); + setValue("S1", 1); + value = getValue("S1"); + CAPTURE(value); + REQUIRE(value == 1); + + // now the global parameter J0_v0 that initially is 8 + value = getValue("J0_v0"); + CAPTURE(value); + REQUIRE(value == 8); + setValue("J0_v0", 10); + value = getValue("J0_v0"); + CAPTURE(value); + REQUIRE(value == 10); + + // and compartment that is initially 1 + value = getValue("compartment"); + CAPTURE(value); + REQUIRE(value == 1); + setValue("compartment", 2); + value = getValue("compartment"); + CAPTURE(value); + REQUIRE(value == 2); + + reset(); + setValue("S1", 1); + + data = simulateEx(0, 10, 11); + CAPTURE(data); + json = nlohmann::json::parse(data); + REQUIRE(json["titles"][0] == "Time"); + REQUIRE(json["titles"][1] == "S1"); + REQUIRE(json["columns"][0].size() == json["recorded_steps"].get()); + REQUIRE(json["columns"][1][0] == 1); + + // reset again and check that the value is back to 0 + reset(); + data = simulateEx(0, 10, 11); + CAPTURE(data); + json = nlohmann::json::parse(data); + REQUIRE(json["titles"][1] == "S1"); + REQUIRE(json["columns"][1][0] == 0); + +} + TEST_CASE("Load COVID Model", "[copasijs][slow]") {