diff --git a/src/OdeIntegrators/tests/OdeIntegratores_base_tests.cpp b/src/OdeIntegrators/tests/OdeIntegratores_base_tests.cpp
index af12bcc69cd9c553ab3bffe447e9b39efc8e9c61..ce534d2dab8214e5b3ebb9699860e9a6b63b8fb5 100644
--- a/src/OdeIntegrators/tests/OdeIntegratores_base_tests.cpp
+++ b/src/OdeIntegrators/tests/OdeIntegratores_base_tests.cpp
@@ -86,7 +86,11 @@ void Exponential_struct_ofp2( const state_type_3d_ofp &x , state_type_3d_ofp &dx
void Exponential( const state_type &x , state_type &dxdt , const double t )
{
+ //timer tt;
+ //tt.start();
dxdt = x;
+ //tt.stop();
+ //std::cout<<tt.getwct()<<std::endl;
}
void sigmoid( const state_type &x , state_type &dxdt , const double t )
{
@@ -97,7 +101,7 @@ BOOST_AUTO_TEST_SUITE(odeInt_BASE_tests)
BOOST_AUTO_TEST_CASE(odeint_base_test1)
{
- size_t edgeSemiSize = 40;
+ size_t edgeSemiSize = 512;
const size_t sz[2] = {edgeSemiSize,edgeSemiSize };
Box<2, double> box({ 0, 0 }, { 1.0, 1.0 });
size_t bc[2] = { NON_PERIODIC, NON_PERIODIC };
@@ -109,7 +113,9 @@ BOOST_AUTO_TEST_CASE(odeint_base_test1)
BOOST_TEST_MESSAGE("Init vector_dist...");
vector_dist<2, double, aggregate<double, double,double>> Particles(0, box, bc, ghost);
-
+ double t0=-5,tf=5,t;
+ t=t0;
+ const double dt=0.01;
auto it = Particles.getGridIterator(sz);
while (it.isNext())
{
@@ -121,8 +127,8 @@ BOOST_AUTO_TEST_CASE(odeint_base_test1)
mem_id k1 = key.get(1);
double yp0 = k1 * spacing[1];
Particles.getLastPos()[1] = yp0;
- Particles.getLastProp<0>() = xp0*yp0*exp(0);
- Particles.getLastProp<1>() = xp0*yp0*exp(0.4);
+ Particles.getLastProp<0>() = xp0*yp0*exp(t0);
+ Particles.getLastProp<1>() = xp0*yp0*exp(tf);
++it;
}
@@ -136,14 +142,15 @@ BOOST_AUTO_TEST_CASE(odeint_base_test1)
x0=Init;
// The rhs of x' = f(x)
- double t=0,tf=0.4;
- const double dt=0.1;
-
- //This doesnt work Why?
- //size_t steps=boost::numeric::odeint::integrate(Exponential,x0,0.0,tf,dt);
- size_t steps=boost::numeric::odeint::integrate_const( boost::numeric::odeint::runge_kutta4< state_type >(),Exponential,x0,0.0,tf,dt);
+ //This doesnt work Why?
+ //size_t steps=boost::numeric::odeint::integrate(Exponential,x0,t,tf,dt);
+ timer tt;
+ tt.start();
+ size_t steps=boost::numeric::odeint::integrate_const( boost::numeric::odeint::runge_kutta4< state_type >(),Exponential,x0,t,tf,dt);
+ tt.stop();
+ std::cout<<"Time taken by CPU:"<<tt.getwct()<<std::endl;
OdeSol=x0;
auto it2 = Particles.getDomainIterator();
double worst = 0.0;
@@ -157,17 +164,18 @@ BOOST_AUTO_TEST_CASE(odeint_base_test1)
//std::cout<<worst<<std::endl;
BOOST_REQUIRE(worst < 1e-6);
-
+ t=t0;
x0=Init;
+ //Particles.write("first");
boost::numeric::odeint::runge_kutta4< state_type > rk4;
- while (t<tf)
+ for(int i=0;i<floor((tf-t0)/dt+0.5);i++)
{
rk4.do_step(Exponential,x0,t,dt);
- OdeSol=x0;
t+=dt;
}
-
+ //std::cout<<"Final TIme:"<<t<<std::endl;
OdeSol=x0;
+ //Particles.write("second");
auto it3 = Particles.getDomainIterator();
double worst2 = 0.0;
while (it3.isNext()) {
diff --git a/src/OdeIntegrators/tests/Odeintegrators_test_gpu.cu b/src/OdeIntegrators/tests/Odeintegrators_test_gpu.cu
index a1b2771dc39896409326a3eed476dcaf0bb76de3..afd06528fcf4498c0889472f6396a2f6e471ef11 100644
--- a/src/OdeIntegrators/tests/Odeintegrators_test_gpu.cu
+++ b/src/OdeIntegrators/tests/Odeintegrators_test_gpu.cu
@@ -24,7 +24,11 @@ typedef state_type_1d_ofp_gpu state_type;
void ExponentialGPU( const state_type &x , state_type &dxdt , const double t )
{
+ //timer tt;
+ //tt.startGPU();
dxdt.data.get<0>() = x.data.get<0>();
+ //tt.stopGPU();
+ //std::cout<<"GPU Time:"<<tt.getwctGPU()<<std::endl;
//x.data.get<0>().getVector().deviceToHost<0>();
//dxdt.data.get<0>().getVector().deviceToHost<0>();
}
diff --git a/src/Operators/Vector/vector_dist_operators.hpp b/src/Operators/Vector/vector_dist_operators.hpp
index 3d330a0f93ddc24027c6a7abae4c07f5b8517db5..1577ee21c6b8e1b0a0aa87495db427839f37adef 100644
--- a/src/Operators/Vector/vector_dist_operators.hpp
+++ b/src/Operators/Vector/vector_dist_operators.hpp
@@ -835,14 +835,40 @@ struct vector_dist_expression_comp_proxy_sel
{ vector_dist_op_compute_op<0,false,vector_dist_expression_comp_sel<comp_dev,cond_>::type::value>
::compute_expr(v,v_exp);}
};
+
+template<unsigned int prp, typename vector>
+class vector_dist_expression;
+
+template<typename v_exp>
+struct transform_if_temporal
+{
+ template<typename T>
+ static auto transform(T & v) -> decltype(v)
+ {
+ return v;
+ }
+};
+
+template<typename T>
+struct transform_if_temporal<vector_dist_expression<0,openfpm::vector_gpu<aggregate<T>>>>
+{
+ template<typename T_>
+ static auto transform(T_ & v) -> decltype(v.getVector().toKernel())
+ {
+ return v.getVector().toKernel();
+ }
+};
+
+
template<>
struct vector_dist_expression_comp_proxy_sel<false>
{
template<bool cond, typename v_type, typename exp_type>
static void compute(v_type &v, exp_type &v_exp)
{ auto v_ker=v.toKernel();
+ auto v_exp_transformed = transform_if_temporal<typename std::remove_const<exp_type>::type>::transform(v_exp);
vector_dist_op_compute_op<0,false,vector_dist_expression_comp_sel<comp_dev,cond>::type::value>
- ::compute_expr(v_ker,v_exp);}
+ ::compute_expr(v_ker,v_exp_transformed);}
};
template<typename vector, bool is_ker = has_vector_kernel<vector>::type::value>
@@ -1512,12 +1538,24 @@ public:
{
return base::operator=(v_exp);
}
+
+ vector & operator=(const vector_dist_expression<0,openfpm::vector_gpu<aggregate<T>>> & v_exp)
+ {
+ return base::operator=(v_exp);
+ }
+
+ vector & equalGPU(const vector_dist_expression<0,openfpm::vector_gpu<aggregate<T>>> & v_exp)
+ {
+ return base::operator=(v_exp);
+ }
+
template<typename exp1, typename exp2, unsigned int op>
vector & operator=(const vector_dist_expression_op<exp1,exp2,op> & v_exp)
{
return base::operator=(v_exp);
}
+
};
template<typename T> using texp_v = vector_dist_expression<0,openfpm::vector<aggregate<T>>>;