diff --git a/include/vectorwrapper/vectorwrapper.hpp b/include/vectorwrapper/vectorwrapper.hpp
index 2feb29d..1654b9e 100644
--- a/include/vectorwrapper/vectorwrapper.hpp
+++ b/include/vectorwrapper/vectorwrapper.hpp
@@ -169,6 +169,8 @@ namespace vwr {
 		private:
 			template <size_type... I, typename... Args>
 			void assign_values (const bt::number_seq<size_type, I...>&, Args... parArgs);
+			template <typename Op, typename V2, size_type... I>
+			void assign_values_op (Op parOp, const bt::number_seq<size_type, I...>& parSeq, const VecBase<V2>& parOther);
 
 			vector_type m_wrapped;
 		};
diff --git a/include/vectorwrapper/vectorwrapper.inl b/include/vectorwrapper/vectorwrapper.inl
index ccca6e1..7fe12e2 100644
--- a/include/vectorwrapper/vectorwrapper.inl
+++ b/include/vectorwrapper/vectorwrapper.inl
@@ -55,6 +55,12 @@ namespace vwr {
 			static_cast<void>(t);
 		}
 
+		template <typename V>
+		template <typename Op, typename V2, size_type... I>
+		void VecBase<V>::assign_values_op (Op parOp, const bt::number_seq<size_type, I...>& parSeq, const VecBase<V2>& parOther) {
+			this->assign_values(parSeq, parOp((*this)[I], parOther[I])...);
+		}
+
 		template <typename V>
 		auto VecBase<V>::operator[] (size_type parIndex) -> scalar_type& {
 			return VecGetter<V>::get_at(m_wrapped, parIndex);
@@ -95,36 +101,28 @@ namespace vwr {
 		template <typename V2>
 		VecBase<V>& VecBase<V>::operator+= (const VecBase<V2>& parOther) {
 			static_assert(static_cast<int>(VectorWrapperInfo<V>::dimensions) == static_cast<int>(VectorWrapperInfo<V2>::dimensions), "Dimensions mismatch");
-			for (size_type z = 0; z < VectorWrapperInfo<V>::dimensions; ++z) {
-				(*this)[z] += parOther[z];
-			}
+			this->assign_values_op(std::plus<scalar_type>(), bt::number_range<size_type, 0, VectorWrapperInfo<V>::dimensions>(), parOther);
 			return *this;
 		}
 		template <typename V>
 		template <typename V2>
 		VecBase<V>& VecBase<V>::operator-= (const VecBase<V2>& parOther) {
 			static_assert(static_cast<int>(VectorWrapperInfo<V>::dimensions) == static_cast<int>(VectorWrapperInfo<V2>::dimensions), "Dimensions mismatch");
-			for (size_type z = 0; z < VectorWrapperInfo<V>::dimensions; ++z) {
-				(*this)[z] -= parOther[z];
-			}
+			this->assign_values_op(std::minus<scalar_type>(), bt::number_range<size_type, 0, VectorWrapperInfo<V>::dimensions>(), parOther);
 			return *this;
 		}
 		template <typename V>
 		template <typename V2>
 		VecBase<V>& VecBase<V>::operator*= (const VecBase<V2>& parOther) {
 			static_assert(static_cast<int>(VectorWrapperInfo<V>::dimensions) == static_cast<int>(VectorWrapperInfo<V2>::dimensions), "Dimensions mismatch");
-			for (size_type z = 0; z < VectorWrapperInfo<V>::dimensions; ++z) {
-				(*this)[z] *= parOther[z];
-			}
+			this->assign_values_op(std::multiplies<scalar_type>(), bt::number_range<size_type, 0, VectorWrapperInfo<V>::dimensions>(), parOther);
 			return *this;
 		}
 		template <typename V>
 		template <typename V2>
 		VecBase<V>& VecBase<V>::operator/= (const VecBase<V2>& parOther) {
 			static_assert(static_cast<int>(VectorWrapperInfo<V>::dimensions) == static_cast<int>(VectorWrapperInfo<V2>::dimensions), "Dimensions mismatch");
-			for (int z = 0; z < VectorWrapperInfo<V>::dimensions; ++z) {
-				(*this)[z] /= parOther[z];
-			}
+			this->assign_values_op(std::divides<scalar_type>(), bt::number_range<size_type, 0, VectorWrapperInfo<V>::dimensions>(), parOther);
 			return *this;
 		}
 
diff --git a/test/unit/test_operators.cpp b/test/unit/test_operators.cpp
index 4c1cad4..3e839ea 100644
--- a/test/unit/test_operators.cpp
+++ b/test/unit/test_operators.cpp
@@ -137,3 +137,32 @@ TEST(vwr, bin_operators_scalar) {
 		EXPECT_EQ(res, 1000 % a);
 	}
 }
+
+TEST(vwr, bin_assign_op) {
+	using namespace vwr;
+
+	{
+		ivec3 a(2, 4, 8);
+		ivec3 res(2 + 20, 4 + 20, 8 + 20);
+		a += ivec3(20);
+		EXPECT_EQ(res, a);
+	}
+	{
+		ivec3 a(2, 4, 8);
+		ivec3 res(2 - 20, 4 - 20, 8 - 20);
+		a -= ivec3(20);
+		EXPECT_EQ(res, a);
+	}
+	{
+		ivec3 a(2, 4, 8);
+		ivec3 res(2 * 20, 4 * 20, 8 * 20);
+		a *= ivec3(20);
+		EXPECT_EQ(res, a);
+	}
+	{
+		ivec3 a(2, 4, 8);
+		ivec3 res(2 / 2, 4 / 2, 8 / 2);
+		a /= ivec3(2);
+		EXPECT_EQ(res, a);
+	}
+}