diff --git a/tdescreensaver/kdesavers/rkodesolver.cpp b/tdescreensaver/kdesavers/rkodesolver.cpp index 2127a110..4c581416 100644 --- a/tdescreensaver/kdesavers/rkodesolver.cpp +++ b/tdescreensaver/kdesavers/rkodesolver.cpp @@ -10,6 +10,7 @@ // //============================================================================ +#include #include #include "rkodesolver.h" diff --git a/tdescreensaver/kdesavers/rotation.cpp b/tdescreensaver/kdesavers/rotation.cpp index 9d9a1442..a00f7bbc 100644 --- a/tdescreensaver/kdesavers/rotation.cpp +++ b/tdescreensaver/kdesavers/rotation.cpp @@ -123,9 +123,9 @@ std::valarray EulerOdeSolver::f( ypr[2] = -(B-A)/C * omega_body[0] * omega_body[1]; // r' // e1', e2', e3' - ypr[std::slice(3,3,1)] = vec3::crossprod(omega, e1); - ypr[std::slice(6,3,1)] = vec3::crossprod(omega, e2); - ypr[std::slice(9,3,1)] = vec3::crossprod(omega, e3); + ypr[std::slice(3,3,1)] = std::valarray(vec3::crossprod(omega, e1)); + ypr[std::slice(6,3,1)] = std::valarray(vec3::crossprod(omega, e2)); + ypr[std::slice(9,3,1)] = std::valarray(vec3::crossprod(omega, e3)); return ypr; } @@ -529,15 +529,17 @@ void KRotationSaver::initData() e3_body.rotate(-m_initEulerTheta*e1_body); // omega_body = L_body * J_body^(-1) vec3 omega_body = e3_body * m_Lz; - omega_body /= J; + std::valarray &omega_body_ = omega_body; + std::valarray &J_ = J; + omega_body_ /= J_; // assemble initial y for solver std::valarray y(12); - y[std::slice(0,3,1)] = omega_body; + y[std::slice(0,3,1)] = std::valarray(omega_body); // 3 basis vectors of body system in fixed coordinates - y[std::slice(3,3,1)] = e1t; - y[std::slice(6,3,1)] = e2t; - y[std::slice(9,3,1)] = e3t; + y[std::slice(3,3,1)] = std::valarray(e1t); + y[std::slice(6,3,1)] = std::valarray(e2t); + y[std::slice(9,3,1)] = std::valarray(e3t); // initial rotation vector omega diff --git a/tdescreensaver/kdesavers/vec3.cpp b/tdescreensaver/kdesavers/vec3.cpp index c66730e5..1715937e 100644 --- a/tdescreensaver/kdesavers/vec3.cpp +++ b/tdescreensaver/kdesavers/vec3.cpp @@ -59,7 +59,9 @@ vec3& vec3::rotate(const vec3& r) // part of vector which is parallel to r vec3 par(r*(*this)/(r*r) * r); // part of vector which is perpendicular to r - vec3 perp(*this - par); + vec3 perp(*this); + std::valarray &perp_ = perp; + perp -= std::valarray(par); // rotation direction, size of perp vec3 rotdir(norm(perp) * normalized(crossprod(r,perp))); *this = par + cos(phi)*perp + sin(phi)*rotdir; diff --git a/tdescreensaver/kdesavers/vec3.h b/tdescreensaver/kdesavers/vec3.h index 62aeb8ed..879c55f7 100644 --- a/tdescreensaver/kdesavers/vec3.h +++ b/tdescreensaver/kdesavers/vec3.h @@ -148,7 +148,9 @@ inline vec3::vec3(const std::slice_array& a) template inline vec3 operator+(vec3 a, const vec3& b) { - a += b; /* valarray::operator+=(const valarray&) */ + std::valarray &a_ = a; + const std::valarray &b_ = b; + a_ += b_; /* valarray::operator+=(const valarray&) */ return a; } @@ -158,7 +160,9 @@ inline vec3 operator+(vec3 a, const vec3& b) template inline T operator*(vec3 a, const vec3& b) { - a *= b; /* valarray::operator*=(const T&) */ + std::valarray &a_ = a; + const std::valarray &b_ = b; + a_ *= b_; /* valarray::operator*=(const valarray&) */ return a.sum(); }