SymEngine

.github/workflows/ci.yml

@@ -22,17 +22,19 @@ jobs:
           - name: 'Windows static - C++'
             os: windows-latest
             python_bindings: OFF
-            additional_cmake_options: -DLibXml2_DIR="C:\libxml2\libxml2-2.9.10\CMake" -DZLIB_DIR="C:\zlib\lib\cmake\ZLIB-1.2.12"
+            additional_cmake_options: -DLibXml2_DIR="C:\libxml2\libxml2-2.9.10\CMake" -DSYMENGINE_DIR="C:\symengine\CMake" -DZLIB_DIR="C:\zlib\lib\cmake\ZLIB-1.2.12"
           - name: 'Windows shared - C++/Python'
             os: windows-latest
             python_bindings: ON
-            additional_cmake_options: -DLibXml2_DIR="C:\libxml2\libxml2-2.9.10\CMake" -DZLIB_DIR="C:\zlib\lib\cmake\ZLIB-1.2.12"
+            additional_cmake_options: -DLibXml2_DIR="C:\libxml2\libxml2-2.9.10\CMake" -DSYMENGINE_DIR="C:\symengine\CMake" -DZLIB_DIR="C:\zlib\lib\cmake\ZLIB-1.2.12"
           - name: 'Linux static - C++'
             os: ubuntu-latest
             python_bindings: OFF
+            additional_cmake_options: -DSYMENGINE_DIR="$HOME/symengine/build/cmake"
           - name: 'Linux shared - C++/Python'
             os: ubuntu-latest
             python_bindings: ON
+            additional_cmake_options: -DSYMENGINE_DIR="$HOME/symengine/build/cmake"
           - name: 'macOS static - C++ (Intel)'
             os: macos-15-intel
             python_bindings: OFF
@@ -68,6 +70,23 @@ jobs:
           cd C:\
           curl -L https://github.com/cellml/gha/releases/download/gha/libxml2-Windows.tar.gz -o libxml2.tar.gz -s
           tar -xzf libxml2.tar.gz
+      - name: Install SymEngine (Windows only)
+        if: ${{ runner.os == 'Windows' }}
+        run: |
+          cd C:\
+          curl -L https://github.com/cellml/gha/releases/download/gha/symengine-Windows.tar.gz -o symengine.tar.gz -s
+          tar -xzf symengine.tar.gz
+      - name: Install SymEngine (Linux only)
+        if: ${{ runner.os == 'Linux' }}
+        run: |
+          cd $HOME
+          curl -L https://github.com/cellml/gha/releases/download/gha/symengine-Linux.tar.gz -o symengine.tar.gz -s
+          tar -xzf symengine.tar.gz
+      - name: Install SymEngine (macOS only)
+        if: ${{ runner.os == 'macOS' }}
+        run: |
+          brew update
+          brew install symengine
       - name: Install zlib (Windows only)
         if: ${{ runner.os == 'Windows' }}
         run: |
@@ -101,13 +120,17 @@ jobs:
         run: |
           cd $HOME
           wget https://github.com/cellml/gha/releases/download/gha/libxml2-WASM.tar.gz -O - | tar -xz
+      - name: Install SymEngine
+        run: |
+          cd $HOME
+          wget https://github.com/cellml/gha/releases/download/gha/symengine-WASM.tar.gz -O - | tar -xz
       - name: Install zlib
         run: |
           cd $HOME
           wget https://github.com/cellml/gha/releases/download/gha/zlib-WASM.tar.gz -O - | tar -xz
       - name: Configure libCellML
         run: |
-          emcmake cmake -G Ninja -S . -B build-wasm -DBUILD_TYPE=Release -DLIBXML2_INCLUDE_DIR=$HOME/libxml2/include/libxml2 -DLIBXML2_LIBRARY=$HOME/libxml2/lib/libxml2.a -DZLIB_INCLUDE_DIR=$HOME/zlib/include -DZLIB_LIBRARY=$HOME/zlib/lib/libz.a
+          emcmake cmake -G Ninja -S . -B build-wasm -DBUILD_TYPE=Release -DLIBXML2_INCLUDE_DIR=$HOME/libxml2/include/libxml2 -DLIBXML2_LIBRARY=$HOME/libxml2/lib/libxml2.a -DSYMENGINE_INCLUDE_DIR=$HOME/symengine/include -DSYMENGINE_LIBRARY=$HOME/symengine/lib/libsymengine.a -DZLIB_INCLUDE_DIR=$HOME/zlib/include -DZLIB_LIBRARY=$HOME/zlib/lib/libz.a
       - name: Build libCellML
         run: cmake --build build-wasm
       - name: Unit testing

cmake/environmentchecks.cmake

@@ -183,6 +183,13 @@ if(NOT DEFINED _ZLIB_FIND_REPORTED)
   message(STATUS "Found ZLIB: ${ZLIB_LIBRARIES} (found version \"${ZLIB_VERSION_STRING}\").")
 endif()
 
+# Set the minimum policy version to work around SymEngine's outdated CMake requirements.
+set(CMAKE_POLICY_VERSION_MINIMUM 3.10)
+
+# Find SymEngine.
+find_package(SymEngine REQUIRED CONFIG)
+message(STATUS "Found SymEngine: ${SYMENGINE_LIBRARIES} (found version \"${SymEngine_VERSION}\").")
+
 if(BUILDCACHE_EXE OR CLCACHE_EXE OR CCACHE_EXE)
   set(COMPILER_CACHE_AVAILABLE TRUE CACHE INTERNAL "Executable required to cache compilations.")
 endif()

src/3rdparty/symengine/symenginebegin.h

@@ -0,0 +1,21 @@
+/*
+Copyright libOpenCOR contributors.
+
+Licensed under the Apache License, Version 2.0 (the "License");
+you may not use this file except in compliance with the License.
+You may obtain a copy of the License at
+
+    http://www.apache.org/licenses/LICENSE-2.0
+
+Unless required by applicable law or agreed to in writing, software
+distributed under the License is distributed on an "AS IS" BASIS,
+WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+See the License for the specific language governing permissions and
+limitations under the License.
+*/
+
+#ifdef __clang__
+#    pragma clang diagnostic push
+#    pragma clang diagnostic ignored "-Wdeprecated-literal-operator"
+#    pragma clang diagnostic ignored "-Wunused-parameter"
+#endif

src/3rdparty/symengine/symengineend.h

@@ -0,0 +1,19 @@
+/*
+Copyright libOpenCOR contributors.
+
+Licensed under the Apache License, Version 2.0 (the "License");
+you may not use this file except in compliance with the License.
+You may obtain a copy of the License at
+
+    http://www.apache.org/licenses/LICENSE-2.0
+
+Unless required by applicable law or agreed to in writing, software
+distributed under the License is distributed on an "AS IS" BASIS,
+WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+See the License for the specific language governing permissions and
+limitations under the License.
+*/
+
+#ifdef __clang__
+#    pragma clang diagnostic pop
+#endif

src/CMakeLists.txt

@@ -203,6 +203,13 @@ endif()
 
 apply_libxml2_settings(cellml)
 
+target_link_libraries(cellml PRIVATE ${SYMENGINE_LIBRARIES})
+target_include_directories(cellml
+  PRIVATE
+  $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/3rdparty/symengine>
+  ${SYMENGINE_INCLUDE_DIRS}
+)
+
 # Use target compile features to propagate features to consuming projects.
 target_compile_features(cellml PUBLIC cxx_std_17)
 

src/analyser.cpp

@@ -23,6 +23,10 @@ limitations under the License.
 #include <cmath>
 #include <iterator>
 
+#include "symenginebegin.h"
+#include <symengine/solve.h>
+#include "symengineend.h"
+
 #include "libcellml/analyserequation.h"
 #include "libcellml/analyserexternalvariable.h"
 #include "libcellml/generatorprofile.h"
@@ -198,6 +202,133 @@ bool AnalyserInternalEquation::variableOnLhsOrRhs(const AnalyserInternalVariable
            || variableOnRhs(variable);
 }
 
+SymEngineEquationResult AnalyserInternalEquation::symEngineEquation(const AnalyserEquationAstPtr &ast, const SymEngineSymbolMap &symbolMap)
+{
+    if (ast == nullptr) {
+        return {true, SymEngine::null};
+    }
+
+    AnalyserEquationAstPtr leftAst = ast->leftChild();
+    AnalyserEquationAstPtr rightAst = ast->rightChild();
+
+    // Recursively call getConvertedAst on left and right children.
+    auto [leftSuccess, left] = symEngineEquation(leftAst, symbolMap);
+    auto [rightSuccess, right] = symEngineEquation(rightAst, symbolMap);
+
+    if (!leftSuccess || !rightSuccess) {
+        return {false, SymEngine::null};
+    }
+
+    // Analyse mAst current type and value.
+    switch (ast->type()) {
+    case AnalyserEquationAst::Type::EQUALITY:
+        return {true, Eq(left, right)};
+    case AnalyserEquationAst::Type::PLUS:
+        return {true, add(left, right)};
+    case AnalyserEquationAst::Type::CI:
+        // Seems like the voi doesn't exist in mAllVariables, so we don't have an easy means of access.
+        if (symbolMap.find(ast->variable()->name()) == symbolMap.end()) {
+            return {false, SymEngine::null};
+        }
+        return {true, symbolMap.at(ast->variable()->name())};
+    default:
+        // Rearrangement is not possible with this type.
+        return {false, SymEngine::null};
+    }
+}
+
+AnalyserEquationAstPtr AnalyserInternalEquation::parseSymEngineExpression(const SymEngine::RCP<const SymEngine::Basic> &seExpression,
+                                                                          const AnalyserEquationAstPtr &parentAst,
+                                                                          const SymEngineVariableMap &variableMap)
+{
+    auto children = seExpression->get_args();
+
+    AnalyserEquationAstPtr ast = AnalyserEquationAst::create();
+
+    ast->setParent(parentAst);
+
+    switch (seExpression->get_type_code()) {
+    case SymEngine::SYMENGINE_EQUALITY: {
+        ast->setType(AnalyserEquationAst::Type::EQUALITY);
+        break;
+    }
+    case SymEngine::SYMENGINE_ADD: {
+        ast->setType(AnalyserEquationAst::Type::PLUS);
+        break;
+    }
+    case SymEngine::SYMENGINE_MUL: {
+        ast->setType(AnalyserEquationAst::Type::TIMES);
+        break;
+    }
+    case SymEngine::SYMENGINE_SYMBOL: {
+        SymEngine::RCP<const SymEngine::Symbol> symbolExpr = SymEngine::rcp_dynamic_cast<const SymEngine::Symbol>(seExpression);
+        ast->setType(AnalyserEquationAst::Type::CI);
+        ast->setVariable(variableMap.at(symbolExpr)->mVariable);
+        break;
+    }
+    case SymEngine::SYMENGINE_INTEGER: {
+        ast->setType(AnalyserEquationAst::Type::CN);
+        ast->setValue(seExpression->__str__());
+        break;
+    }
+    default:
+        break;
+    }
+
+    // TODO Update to account for symengine expressions with 3 or more children.
+    if (children.size() > 0) {
+        ast->setLeftChild(parseSymEngineExpression(children[0], ast, variableMap));
+        if (children.size() > 1) {
+            ast->setRightChild(parseSymEngineExpression(children[1], ast, variableMap));
+        }
+    }
+
+    return ast;
+}
+
+AnalyserEquationAstPtr AnalyserInternalEquation::rearrangeFor(const AnalyserInternalVariablePtr &variable)
+{
+    SymEngineSymbolMap symbolMap;
+    SymEngineVariableMap variableMap;
+
+    for (const auto &variable : mAllVariables) {
+        SymEngine::RCP<const SymEngine::Symbol> symbol = SymEngine::symbol(variable->mVariable->name());
+        symbolMap[variable->mVariable->name()] = symbol;
+        variableMap[symbol] = variable;
+    }
+
+    auto [success, seEquation] = symEngineEquation(mAst, symbolMap);
+    if (!success) {
+        return nullptr;
+    }
+
+    SymEngine::RCP<const SymEngine::Set> solutionSet = solve(seEquation, symbolMap[variable->mVariable->name()]);
+    SymEngine::vec_basic solutions = solutionSet->get_args();
+
+    // Our system needs to be able to isolate a single solution.
+    if (solutions.size() != 1) {
+        return nullptr;
+    }
+    SymEngine::RCP<const SymEngine::Basic> answer = solutions.front();
+
+    // Rebuild the AST from the rearranged expression.
+    AnalyserEquationAstPtr ast = AnalyserEquationAst::create();
+    AnalyserEquationAstPtr isolatedVariableAst = AnalyserEquationAst::create();
+    AnalyserEquationAstPtr rearrangedEquationAst = parseSymEngineExpression(answer, nullptr, variableMap);
+
+    ast->setType(AnalyserEquationAst::Type::EQUALITY);
+    ast->setLeftChild(isolatedVariableAst);
+    ast->setRightChild(rearrangedEquationAst);
+
+    isolatedVariableAst->setType(AnalyserEquationAst::Type::CI);
+    isolatedVariableAst->setVariable(variable->mVariable);
+    isolatedVariableAst->setParent(ast);
+
+    rearrangedEquationAst->setParent(ast);
+
+    return ast;
+}
+
 bool AnalyserInternalEquation::check(const AnalyserModelPtr &analyserModel, bool checkNlaSystems)
 {
     // Nothing to check if the equation has a known type.
@@ -278,6 +409,15 @@ bool AnalyserInternalEquation::check(const AnalyserModelPtr &analyserModel, bool
                                    mVariables.front() :
                                    nullptr;
 
+    // If we have one variable left, but it's not isolated, try to rearrange it.
+    if ((unknownVariableLeft != nullptr) && !variableOnLhsOrRhs(unknownVariableLeft)) {
+        auto newAst = rearrangeFor(unknownVariableLeft);
+        if (newAst != nullptr) {
+            // TODO Update variables and/or equation type when necessary.
+            mAst = newAst;
+        }
+    }
+
     if (((unknownVariableLeft != nullptr)
          && (checkNlaSystems || variableOnLhsOrRhs(unknownVariableLeft)))
         || !initialisedVariables.empty()) {

src/analyser_p.h

@@ -22,6 +22,12 @@ limitations under the License.
 #include "logger_p.h"
 #include "utilities.h"
 
+namespace SymEngine {
+template<class T> class RCP;
+class Basic;
+class Symbol;
+} // namespace SymEngine
+
 namespace libcellml {
 
 struct AnalyserInternalEquation;
@@ -39,6 +45,10 @@ using AnalyserEquationPtrs = std::vector<AnalyserEquationPtr>;
 using AnalyserVariablePtrs = std::vector<AnalyserVariablePtr>;
 using AnalyserExternalVariablePtrs = std::vector<AnalyserExternalVariablePtr>;
 
+using SymEngineVariableMap = std::map<SymEngine::RCP<const SymEngine::Symbol>, AnalyserInternalVariablePtr, SymEngine::RCPBasicKeyLess>;
+using SymEngineSymbolMap = std::map<std::string, SymEngine::RCP<const SymEngine::Symbol>>;
+using SymEngineEquationResult = std::tuple<bool, SymEngine::RCP<const SymEngine::Basic>>;
+
 struct AnalyserInternalVariable
 {
     enum struct Type
@@ -128,6 +138,10 @@ struct AnalyserInternalEquation
     bool variableOnRhs(const AnalyserInternalVariablePtr &variable);
     bool variableOnLhsOrRhs(const AnalyserInternalVariablePtr &variable);
 
+    SymEngineEquationResult symEngineEquation(const AnalyserEquationAstPtr &ast, const SymEngineSymbolMap &symbolMap);
+    AnalyserEquationAstPtr parseSymEngineExpression(const SymEngine::RCP<const SymEngine::Basic> &seExpression, const AnalyserEquationAstPtr &parentAst, const SymEngineVariableMap &variableMap);
+    AnalyserEquationAstPtr rearrangeFor(const AnalyserInternalVariablePtr &variable);
+
     bool check(const AnalyserModelPtr &analyserModel, bool checkNlaSystems);
 };
 

tests/analyser/analyser.cpp

@@ -1060,3 +1060,17 @@ TEST(Analyser, unsuitablyConstrainedNlaSystem)
 
     EXPECT_EQ(libcellml::AnalyserModel::Type::UNSUITABLY_CONSTRAINED, analyser->analyserModel()->type());
 }
+
+TEST(Analyser, rearrangeAlgebraicEquation)
+{
+    auto parser = libcellml::Parser::create();
+    auto model = parser->parseModel(fileContents("analyser/unarranged_algebraic_equation.cellml"));
+
+    EXPECT_EQ(size_t(0), parser->issueCount());
+
+    auto analyser = libcellml::Analyser::create();
+
+    analyser->analyseModel(model);
+
+    EXPECT_EQ(libcellml::AnalyserModel::Type::ALGEBRAIC, analyser->analyserModel()->type());
+}

tests/generator/generator.cpp

@@ -1602,30 +1602,6 @@ TEST(Generator, analyserModelScopeTest)
     EXPECT_EQ_FILE_CONTENTS("generator/hodgkin_huxley_squid_axon_model_1952/model.c", generator->implementationCode(analyserModel));
 }
 
-TEST(Generator, daeModel)
-{
-    auto parser = libcellml::Parser::create(false);
-    auto model = parser->parseModel(fileContents("generator/dae_cellml_1_1_model/model.cellml"));
-
-    EXPECT_EQ(size_t(0), parser->errorCount());
-
-    auto analyser = libcellml::Analyser::create();
-
-    analyser->analyseModel(model);
-
-    EXPECT_EQ(size_t(0), analyser->errorCount());
-
-    auto analyserModel = analyser->analyserModel();
-    auto generator = libcellml::Generator::create();
-
-    EXPECT_EQ_FILE_CONTENTS("generator/dae_cellml_1_1_model/model.h", generator->interfaceCode(analyserModel));
-    EXPECT_EQ_FILE_CONTENTS("generator/dae_cellml_1_1_model/model.c", generator->implementationCode(analyserModel));
-
-    auto profile = libcellml::GeneratorProfile::create(libcellml::GeneratorProfile::Profile::PYTHON);
-
-    EXPECT_EQ_FILE_CONTENTS("generator/dae_cellml_1_1_model/model.py", generator->implementationCode(analyserModel, profile));
-}
-
 TEST(Generator, variableInitialisedUsingAnotherVariable)
 {
     // Note: this should be in sync with the corresponding Analyser test.