Skip to content

Commit

Permalink
Split expressions into static and dynamic
Browse files Browse the repository at this point in the history 
Closes #1269
  • Loading branch information
@dweindl
dweindl committed Feb 23, 2024
1 parent 4334557 commit f2afbdf
Show file tree
Hide file tree
Showing 48 changed files with 722 additions and 515 deletions.
6 changes: 5 additions & 1 deletion include/amici/abstract_model.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -820,11 +820,15 @@ class AbstractModel {
* @param h Heaviside vector
* @param tcl total abundances for conservation laws
* @param spl spline value vector
* @param include_static Whether to (re-)evaluate only dynamic expressions
* (false) or also static expressions (true).
* Dynamic expressions are those that depend directly or indirectly on time,
* static expressions are those that don't.
*/
virtual void
fw(realtype* w, realtype const t, realtype const* x, realtype const* p,
realtype const* k, realtype const* h, realtype const* tcl,
realtype const* spl);
realtype const* spl, bool include_static = true);

/**
* @brief Model-specific sparse implementation of dwdp
Expand Down
14 changes: 13 additions & 1 deletion include/amici/model.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -239,6 +239,14 @@ class Model : public AbstractModel, public ModelDimensions {
bool computeSensitivities, std::vector<int>& roots_found
);

/**
* @brief Re-initialize model properties after changing simulation context.
* @param x Reference to state variables
* @param t Timepoint
*/
void reinitialize(AmiVector& x, realtype t);


/**
* @brief Initialize model properties.
* @param xB Adjoint state variables
Expand Down Expand Up @@ -1828,8 +1836,12 @@ class Model : public AbstractModel, public ModelDimensions {
* @brief Compute recurring terms in xdot.
* @param t Timepoint
* @param x Array with the states
* @param include_static Whether to (re-)evaluate only dynamic expressions
* (false) or also static expressions (true).
* Dynamic expressions are those that depend directly or indirectly on time,
* static expressions are those that don't.
*/
void fw(realtype t, realtype const* x);
void fw(realtype t, realtype const* x, bool include_static = true);

/**
* @brief Compute parameter derivative for recurring terms in xdot.
Expand Down
2 changes: 1 addition & 1 deletion matlab/@amifun/getArgs.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -108,7 +108,7 @@
case 'dJrzdsigma'
this.argstr = '(double *dJrzdsigma, const int iz, const realtype *p, const realtype *k, const double *rz, const double *sigmaz)';
case 'w'
this.argstr = '(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *spl)';
this.argstr = '(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *spl, bool include_static)';
case 'dwdp'
this.argstr = '(realtype *dwdp, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *tcl, const realtype *stcl, const realtype *spl, const realtype *sspl)';
case 'dwdx'
Expand Down
3 changes: 2 additions & 1 deletion matlab/@amimodel/generateC.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -209,7 +209,7 @@ function generateC(this)
end
fprintf(fid,'};\n\n');
fprintf(fid,['} // namespace model_' this.modelname '\n\n']);
fprintf(fid,'} // namespace amici \n\n');
fprintf(fid,'} // namespace amici\n\n');
fprintf(fid,['#endif /* _amici_' this.modelname '_h */\n']);
fclose(fid);

Expand Down Expand Up @@ -253,6 +253,7 @@ function generateC(this)

argstr = strrep(argstr,'realtype','');
argstr = strrep(argstr,'int','');
argstr = strrep(argstr,'bool','');
argstr = strrep(argstr,'const','');
argstr = strrep(argstr,'double','');
argstr = strrep(argstr,'SUNMatrixContent_Sparse','');
Expand Down
37 changes: 24 additions & 13 deletions models/model_calvetti/CMakeLists.txt
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,11 @@
# Build AMICI model
cmake_minimum_required(VERSION 3.15)
cmake_policy(VERSION 3.15...3.27)

# cmake >=3.27
if(POLICY CMP0144)
cmake_policy(SET CMP0144 NEW)
endif(POLICY CMP0144)

project(model_calvetti)

Expand All @@ -14,7 +20,7 @@ if(CMAKE_CXX_COMPILER_ID STREQUAL "GNU")
endif()
foreach(flag ${MY_CXX_FLAGS})
unset(CUR_FLAG_SUPPORTED CACHE)
check_cxx_compiler_flag(-Werror ${flag} CUR_FLAG_SUPPORTED)
check_cxx_compiler_flag(${flag} CUR_FLAG_SUPPORTED)
if(${CUR_FLAG_SUPPORTED})
string(APPEND CMAKE_CXX_FLAGS " ${flag}")
endif()
Expand All @@ -33,6 +39,23 @@ find_package(Amici REQUIRED HINTS
${CMAKE_CURRENT_LIST_DIR}/../../build)
message(STATUS "Found AMICI ${Amici_DIR}")

# Debug build?
if("$ENV{ENABLE_AMICI_DEBUGGING}" OR "$ENV{ENABLE_GCOV_COVERAGE}")
add_compile_options(-UNDEBUG)
if(MSVC)
add_compile_options(-DEBUG)
else()
add_compile_options(-O0 -g)
endif()
set(CMAKE_BUILD_TYPE "Debug")
endif()

# coverage options
if($ENV{ENABLE_GCOV_COVERAGE})
string(APPEND CMAKE_CXX_FLAGS_DEBUG " --coverage")
string(APPEND CMAKE_EXE_LINKER_FLAGS_DEBUG " --coverage")
endif()

set(MODEL_DIR ${CMAKE_CURRENT_LIST_DIR})

set(SRC_LIST_LIB ${MODEL_DIR}/JSparse.cpp
Expand Down Expand Up @@ -73,18 +96,6 @@ if(NOT "${AMICI_PYTHON_BUILD_EXT_ONLY}")
target_link_libraries(simulate_${PROJECT_NAME} ${PROJECT_NAME})
endif()

# Debug build?
if("$ENV{ENABLE_AMICI_DEBUGGING}" OR "$ENV{ENABLE_GCOV_COVERAGE}")
add_compile_options(-UNDEBUG -O0 -g)
set(CMAKE_BUILD_TYPE "Debug")
endif()

# coverage options
if($ENV{ENABLE_GCOV_COVERAGE})
string(APPEND CMAKE_CXX_FLAGS_DEBUG " --coverage")
string(APPEND CMAKE_EXE_LINKER_FLAGS_DEBUG " --coverage")
endif()

# SWIG
option(ENABLE_SWIG "Build swig/python library?" ON)
if(ENABLE_SWIG)
Expand Down
59 changes: 28 additions & 31 deletions models/model_calvetti/main.cpp
View file
Original file line number Diff line number Diff line change
@@ -1,14 +1,13 @@
#include <iostream>

#include <amici/amici.h> /* AMICI base functions */
#include "wrapfunctions.h" /* model-provided functions */
#include "wrapfunctions.h" /* model-provided functions */
#include <amici/amici.h> /* AMICI base functions */

template < class T >
std::ostream& operator << (std::ostream& os, const std::vector<T>& v)
{
template <class T>
std::ostream& operator<<(std::ostream& os, std::vector<T> const& v) {

Check warning on line 7 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L7 

Added line #L7 was not covered by tests
os << "[";
for (typename std::vector<T>::const_iterator ii = v.begin(); ii != v.end(); ++ii)
{
for (typename std::vector<T>::const_iterator ii = v.begin(); ii != v.end();
++ii) {

Check warning on line 10 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L10 

Added line #L10 was not covered by tests
os << " " << *ii;
}
os << "]";
Expand All @@ -21,9 +20,9 @@ std::ostream& operator << (std::ostream& os, const std::vector<T>& v)
*/

int main() {
std::cout<<"********************************"<<std::endl;
std::cout<<"** Running forward simulation **"<<std::endl;
std::cout<<"********************************"<<std::endl<<std::endl;
std::cout << "********************************" << std::endl;
std::cout << "** Running forward simulation **" << std::endl;
std::cout << "********************************" << std::endl << std::endl;

Check warning on line 25 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L23-L25 

Added lines #L23 - L25 were not covered by tests

// Create a model instance
auto model = amici::generic_model::getModel();
Expand All @@ -44,21 +43,20 @@ int main() {

// Print observable time course
auto observable_ids = model->getObservableIds();
std::cout<<"Simulated observables for timepoints "<<rdata->ts<<"\n\n";
for(int i_observable = 0; i_observable < rdata->ny; ++i_observable) {
std::cout<<observable_ids[i_observable]<<":\n\t";
for(int i_time = 0; i_time < rdata->nt; ++i_time) {
std::cout << "Simulated observables for timepoints " << rdata->ts << "\n\n";

Check warning on line 46 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L46 

Added line #L46 was not covered by tests
for (int i_observable = 0; i_observable < rdata->ny; ++i_observable) {
std::cout << observable_ids[i_observable] << ":\n\t";

Check warning on line 48 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L48 

Added line #L48 was not covered by tests
for (int i_time = 0; i_time < rdata->nt; ++i_time) {
// rdata->y is a flat 2D array in row-major ordering
std::cout<<rdata->y[i_time * rdata->ny + i_observable]<<" ";
std::cout << rdata->y[i_time * rdata->ny + i_observable] << " ";

Check warning on line 51 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L51 

Added line #L51 was not covered by tests
}
std::cout<<std::endl<<std::endl;
std::cout << std::endl << std::endl;

Check warning on line 53 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L53 

Added line #L53 was not covered by tests
}


std::cout<<std::endl;
std::cout<<"**********************************"<<std::endl;
std::cout<<"** Forward sensitivity analysis **"<<std::endl;
std::cout<<"**********************************"<<std::endl<<std::endl;
std::cout << std::endl;
std::cout << "**********************************" << std::endl;
std::cout << "** Forward sensitivity analysis **" << std::endl;
std::cout << "**********************************" << std::endl << std::endl;

Check warning on line 59 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L56-L59 

Added lines #L56 - L59 were not covered by tests

// Enable first-order sensitivity analysis
solver->setSensitivityOrder(amici::SensitivityOrder::first);
Expand All @@ -78,18 +76,17 @@ int main() {
auto state_ids = model->getStateIds();
auto parameter_ids = model->getParameterIds();

std::cout<<"State sensitivities for timepoint "
<<rdata->ts[i_time]
<<std::endl;//nt x nplist x nx
for(int i_state= 0; i_state < rdata->nx; ++i_state) {
std::cout<<"\td("<<state_ids[i_state]<<")/d("
<<parameter_ids[model->plist(i_nplist)]<<") = ";
std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
<< std::endl; // nt x nplist x nx

Check warning on line 80 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L79-L80 

Added lines #L79 - L80 were not covered by tests
for (int i_state = 0; i_state < rdata->nx; ++i_state) {
std::cout << "\td(" << state_ids[i_state] << ")/d("
<< parameter_ids[model->plist(i_nplist)] << ") = ";

Check warning on line 83 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L82-L83 

Added lines #L82 - L83 were not covered by tests

// rdata->sx is a flat 3D array in row-major ordering
std::cout<<rdata->sx[i_time * rdata->nplist * rdata->nx
+ i_nplist * rdata->nx
+ i_state];
std::cout<<std::endl;
std::cout << rdata->sx
[i_time * rdata->nplist * rdata->nx
+ i_nplist * rdata->nx + i_state];
std::cout << std::endl;

Check warning on line 89 in models/model_calvetti/main.cpp

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/main.cpp#L86-L89 

Added lines #L86 - L89 were not covered by tests
}

return 0;
Expand Down
10 changes: 5 additions & 5 deletions models/model_calvetti/model_calvetti.h
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
#ifndef _amici_model_calvetti_h
#define _amici_model_calvetti_h
/* Generated by amiwrap (R2017b) 223c5075608273b17f304556d0e8ccb41233bd21 */
/* Generated by amiwrap (R2017b) 8f7983a82d96eb59b83557581712b4741dcf4468 */
#include <cmath>
#include <memory>
#include "amici/defines.h"
Expand All @@ -23,7 +23,7 @@ extern void dwdx_model_calvetti(realtype *dwdx, const realtype t, const realtype
extern void dydx_model_calvetti(double *dydx, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w, const realtype *dwdx);
extern void root_model_calvetti(realtype *root, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *dx);
extern void sigmay_model_calvetti(double *sigmay, const realtype t, const realtype *p, const realtype *k, const realtype *y);
extern void w_model_calvetti(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *spl);
extern void w_model_calvetti(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *spl, bool include_static);
extern void x0_model_calvetti(realtype *x0, const realtype t, const realtype *p, const realtype *k);
extern void xdot_model_calvetti(realtype *xdot, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *dx, const realtype *w);
extern void y_model_calvetti(double *y, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *w);
Expand Down Expand Up @@ -72,7 +72,7 @@ class Model_model_calvetti : public amici::Model_DAE {

amici::Model* clone() const override { return new Model_model_calvetti(*this); };

std::string getAmiciCommit() const override { return "223c5075608273b17f304556d0e8ccb41233bd21"; };
std::string getAmiciCommit() const override { return "8f7983a82d96eb59b83557581712b4741dcf4468"; };

Check warning on line 75 in models/model_calvetti/model_calvetti.h

View check run for this annotation

Codecov / codecov/patch

models/model_calvetti/model_calvetti.h#L75 

Added line #L75 was not covered by tests

void fJSparse(SUNMatrixContent_Sparse JSparse, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype cj, const realtype *dx, const realtype *w, const realtype *dwdx) override {
JSparse_model_calvetti(JSparse, t, x, p, k, h, cj, dx, w, dwdx);
Expand Down Expand Up @@ -185,8 +185,8 @@ class Model_model_calvetti : public amici::Model_DAE {
void fsz(double *sz, const int ie, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *sx, const int ip) override {
}

void fw(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *spl) override {
w_model_calvetti(w, t, x, p, k, h, tcl, spl);
void fw(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *spl, bool include_static) override {
w_model_calvetti(w, t, x, p, k, h, tcl, spl, include_static);
}

void fx0(realtype *x0, const realtype t, const realtype *p, const realtype *k) override {
Expand Down
2 changes: 1 addition & 1 deletion models/model_calvetti/w.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@ namespace amici {

namespace model_model_calvetti{

void w_model_calvetti(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *spl) {
void w_model_calvetti(realtype *w, const realtype t, const realtype *x, const realtype *p, const realtype *k, const realtype *h, const realtype *tcl, const realtype *spl, bool include_static) {
w[0] = 1.0/k[0];
w[1] = k[2]*k[2];
w[2] = 1.0/(x[1]*x[1]);
Expand Down
37 changes: 24 additions & 13 deletions models/model_dirac/CMakeLists.txt
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,11 @@
# Build AMICI model
cmake_minimum_required(VERSION 3.15)
cmake_policy(VERSION 3.15...3.27)

# cmake >=3.27
if(POLICY CMP0144)
cmake_policy(SET CMP0144 NEW)
endif(POLICY CMP0144)

project(model_dirac)

Expand All @@ -14,7 +20,7 @@ if(CMAKE_CXX_COMPILER_ID STREQUAL "GNU")
endif()
foreach(flag ${MY_CXX_FLAGS})
unset(CUR_FLAG_SUPPORTED CACHE)
check_cxx_compiler_flag(-Werror ${flag} CUR_FLAG_SUPPORTED)
check_cxx_compiler_flag(${flag} CUR_FLAG_SUPPORTED)
if(${CUR_FLAG_SUPPORTED})
string(APPEND CMAKE_CXX_FLAGS " ${flag}")
endif()
Expand All @@ -33,6 +39,23 @@ find_package(Amici REQUIRED HINTS
${CMAKE_CURRENT_LIST_DIR}/../../build)
message(STATUS "Found AMICI ${Amici_DIR}")

# Debug build?
if("$ENV{ENABLE_AMICI_DEBUGGING}" OR "$ENV{ENABLE_GCOV_COVERAGE}")
add_compile_options(-UNDEBUG)
if(MSVC)
add_compile_options(-DEBUG)
else()
add_compile_options(-O0 -g)
endif()
set(CMAKE_BUILD_TYPE "Debug")
endif()

# coverage options
if($ENV{ENABLE_GCOV_COVERAGE})
string(APPEND CMAKE_CXX_FLAGS_DEBUG " --coverage")
string(APPEND CMAKE_EXE_LINKER_FLAGS_DEBUG " --coverage")
endif()

set(MODEL_DIR ${CMAKE_CURRENT_LIST_DIR})

set(SRC_LIST_LIB ${MODEL_DIR}/JSparse.cpp
Expand Down Expand Up @@ -73,18 +96,6 @@ if(NOT "${AMICI_PYTHON_BUILD_EXT_ONLY}")
target_link_libraries(simulate_${PROJECT_NAME} ${PROJECT_NAME})
endif()

# Debug build?
if("$ENV{ENABLE_AMICI_DEBUGGING}" OR "$ENV{ENABLE_GCOV_COVERAGE}")
add_compile_options(-UNDEBUG -O0 -g)
set(CMAKE_BUILD_TYPE "Debug")
endif()

# coverage options
if($ENV{ENABLE_GCOV_COVERAGE})
string(APPEND CMAKE_CXX_FLAGS_DEBUG " --coverage")
string(APPEND CMAKE_EXE_LINKER_FLAGS_DEBUG " --coverage")
endif()

# SWIG
option(ENABLE_SWIG "Build swig/python library?" ON)
if(ENABLE_SWIG)
Expand Down
Loading

0 comments on commit f2afbdf

Please sign in to comment.