lawmurray / LibBi

Bayesian state-space modelling on high-performance hardware, including multicore, GPUs and distributed clusters.
http://www.libbi.org
Other
96 stars 22 forks source link

Dot product within ode block fails #16

Closed BlackEdder closed 6 years ago

BlackEdder commented 7 years ago

My model involves taking the dot product in to ODE equations, but that currently fails to compile. The relevant lines is:

dS[age]/dt = - (polymod*(I1+I2)) * S[age] - v[age] * S[age]

where polymod is a matrix and I1 and I2 are vectors. Moving this outside of the ode block works. The same problem has been reported on the mailinglist. @sbfnk told me that it did use to work in previous versions.

The relevant part of make.log

depbase=`echo src/model/ModelSIR.o | sed 's|[^/]*$|.deps/&|;s|\.o$||'`;\
g++ -DPACKAGE_NAME=\"LibBi\" -DPACKAGE_TARNAME=\"libbi\" -DPACKAGE_VERSION=\"1.3.0\" -DPACKAGE_STRING=\"LibBi\ 1.3.0\" -DPACKAGE_BUGREPORT=\"bug-report@libbi.org\" -DPACKAGE_URL=\"http://www.libbi.org\" -DHAVE_OMP_H=1 -DTHRUST_DEVICE_SYSTEM=THRUST_DEVICE_SYSTEM_CPP -DTHRUST_HOST_SYSTEM=THRUST_HOST_SYSTEM_OMP -DHAVE_LIBM=1 -DHAVE_LIBGFORTRAN=1 -DHAVE_LIBQRUPDATE=1 -DHAVE_LIBGSL=1 -DHAVE_LIBNETCDF=1 -DHAVE_NETCDF_H=1 -DHAVE_CBLAS_H=1 -DHAVE_GSL_GSL_CBLAS_H=1 -DHAVE_BOOST_MPL_IF_HPP=1 -DHAVE_BOOST_RANDOM_BERNOULLI_DISTRIBUTION_HPP=1 -DHAVE_BOOST_RANDOM_GAMMA_DISTRIBUTION_HPP=1 -DHAVE_BOOST_RANDOM_MERSENNE_TWISTER_HPP=1 -DHAVE_BOOST_RANDOM_NORMAL_DISTRIBUTION_HPP=1 -DHAVE_BOOST_RANDOM_POISSON_DISTRIBUTION_HPP=1 -DHAVE_BOOST_RANDOM_UNIFORM_INT_HPP=1 -DHAVE_BOOST_RANDOM_UNIFORM_REAL_HPP=1 -DHAVE_BOOST_RANDOM_VARIATE_GENERATOR_HPP=1 -DHAVE_BOOST_TYPEOF_TYPEOF_HPP=1 -DHAVE_THRUST_ADJACENT_DIFFERENCE_H=1 -DHAVE_THRUST_BINARY_SEARCH_H=1 -DHAVE_THRUST_COPY_H=1 -DHAVE_THRUST_DEVICE_PTR_H=1 -DHAVE_THRUST_DISTANCE_H=1 -DHAVE_THRUST_EXTREMA_H=1 -DHAVE_THRUST_FILL_H=1 -DHAVE_THRUST_FOR_EACH_H=1 -DHAVE_THRUST_FUNCTIONAL_H=1 -DHAVE_THRUST_GATHER_H=1 -DHAVE_THRUST_INNER_PRODUCT_H=1 -DHAVE_THRUST_ITERATOR_COUNTING_ITERATOR_H=1 -DHAVE_THRUST_ITERATOR_DETAIL_NORMAL_ITERATOR_H=1 -DHAVE_THRUST_ITERATOR_DISCARD_ITERATOR_H=1 -DHAVE_THRUST_ITERATOR_PERMUTATION_ITERATOR_H=1 -DHAVE_THRUST_ITERATOR_TRANSFORM_ITERATOR_H=1 -DHAVE_THRUST_ITERATOR_ZIP_ITERATOR_H=1 -DHAVE_THRUST_LOGICAL_H=1 -DHAVE_THRUST_REDUCE_H=1 -DHAVE_THRUST_SCAN_H=1 -DHAVE_THRUST_SEQUENCE_H=1 -DHAVE_THRUST_SORT_H=1 -DHAVE_THRUST_TRANSFORM_H=1 -DHAVE_THRUST_TRANSFORM_REDUCE_H=1 -DHAVE_THRUST_TRANSFORM_SCAN_H=1 -DHAVE_THRUST_TUPLE_H=1 -DHAVE_GSL_GSL_MULTIMIN_H=1 -DENABLE_DIAGNOSTICS=no -DHAVE_STDBOOL_H=1 -I.  -Isrc  -I/usr/local/cuda/include      -DENABLE_OPENMP    -fopenmp -O3 -g3 -funroll-loops   -MT src/model/ModelSIR.o -MD -MP -MF $depbase.Tpo -c -o src/model/ModelSIR.o src/model/ModelSIR.cpp &&\
mv -f $depbase.Tpo $depbase.Po
In file included from src/model/block/Block6.hpp:18:0,
                 from src/model/block/Block5.hpp:22,
                 from src/model/ModelSIR.hpp:1556,
                 from src/model/ModelSIR.cpp:9:
src/model/block/../action/Action35.hpp: In static member function ‘static void Action35::dfdt(T1, const bi::State<ModelSIR, L>&, int, const CX&, const PX&, T2&)’:
src/model/block/../action/Action35.hpp:73:27: error: ‘un0’ was not declared in this scope
   const VarCoord38 cox0(  un0);
                           ^
In file included from src/model/block/Block10.hpp:18:0,
                 from src/model/block/Block9.hpp:24,
                 from src/model/ModelSIR.hpp:1559,
                 from src/model/ModelSIR.cpp:9:
src/model/block/../action/Action58.hpp: In static member function ‘static void Action58::dfdt(T1, const bi::State<ModelSIR, L>&, int, const CX&, const PX&, T2&)’:
src/model/block/../action/Action58.hpp:73:27: error: ‘un0’ was not declared in this scope
   const VarCoord44 cox0(  un0);
                           ^
In file included from src/model/block/Block6.hpp:18:0,
                 from src/model/block/Block5.hpp:22,
                 from src/model/ModelSIR.hpp:1556,
                 from src/sample_cpu.cpp:21:
src/model/block/../action/Action35.hpp: In static member function ‘static void Action35::dfdt(T1, const bi::State<ModelSIR, L>&, int, const CX&, const PX&, T2&)’:
src/model/block/../action/Action35.hpp:73:27: error: ‘un0’ was not declared in this scope
   const VarCoord38 cox0(  un0);
                           ^
In file included from src/model/block/Block10.hpp:18:0,
                 from src/model/block/Block9.hpp:24,
                 from src/model/ModelSIR.hpp:1559,
                 from src/sample_cpu.cpp:21:
src/model/block/../action/Action58.hpp: In static member function ‘static void Action58::dfdt(T1, const bi::State<ModelSIR, L>&, int, const CX&, const PX&, T2&)’:
src/model/block/../action/Action58.hpp:73:27: error: ‘un0’ was not declared in this scope
   const VarCoord44 cox0(  un0);
                           ^
Makefile:1011: recipe for target 'src/model/ModelSIR.o' failed
make: *** [src/model/ModelSIR.o] Error 1
make: *** Waiting for unfinished jobs....
Makefile:1011: recipe for target 'src/sample_cpu.o' failed
make: *** [src/sample_cpu.o] Error 1
/tmp/Rtmp2p8rQ8/.SIR/build_assert_openmp/make.log (END)
lawmurray commented 7 years ago

Thanks @BlackEdder. Probably related to your other report as you say. Likewise, can you try polymodI1 +polymodI2 as a workaround? Without having looked deeper yet, I suspect the issue is related to the subexpression appearing in the matrix-vector product.

BlackEdder commented 7 years ago

The issue also occurred with a simpler model, where I just used polymod*I

sbfnk commented 6 years ago

Closing this for inactivity. Matrix/vector operations of this sort are currently not supported in ode blocks.