diff --git a/src/pew/camera/depthCamera.hpp b/src/pew/camera/depthCamera.hpp
index e1ffc1534305e298146ba2c3edcddc83323fa2c4..d9e68523831dbaf69542c3c8d2bb81dde1736d52 100644
--- a/src/pew/camera/depthCamera.hpp
+++ b/src/pew/camera/depthCamera.hpp
@@ -31,7 +31,7 @@ namespace pew {
 
             origin = lookFrom;
 
-            double f_focusDistance;
+            vec3 f_focusDistance;
             if (focusDistance.has_value()) {
                 f_focusDistance = focusDistance.value();
             } else {
diff --git a/src/pew/hittable/collection.hpp b/src/pew/hittable/collection.hpp
index 30cc803760b12184c4fedb28f6f0fc72f5188b4b..0fb6a6a942400b6b56c179ca4b7ac3e58996e0b6 100644
--- a/src/pew/hittable/collection.hpp
+++ b/src/pew/hittable/collection.hpp
@@ -22,7 +22,7 @@ namespace pew {
                 if (auto result = hittable->hit(rayCopy)) {
                     if (!closest || closest->t > result->t) {
                         closest = result;
-                        rayCopy.setRange({ rayCopy.getRange().first, result->t });
+                        rayCopy.setRange({ rayCopy.getRange().first, result->t.x() });
                     }
                 }
             }
diff --git a/src/pew/hittable/hittable.hpp b/src/pew/hittable/hittable.hpp
index 241ce3b0036458165e9dfc71ecbcbb43d205bc79..ca277cfe96343edae7429b1fe744f3e7bd320c46 100644
--- a/src/pew/hittable/hittable.hpp
+++ b/src/pew/hittable/hittable.hpp
@@ -47,7 +47,7 @@ namespace pew {
 
     struct Hit {
         Ray                       ray;
-        double                    t;
+        vec3                      t;
         std::shared_ptr<Material> material;
         bool                      frontFace = false;
         vec3                      normal {};
diff --git a/src/pew/hittable/sphere.hpp b/src/pew/hittable/sphere.hpp
index 66a2209663cc8693500cd997505d2bcf7857a9ed..3e47fbfee80d65ada2319584e902a7c196fe0038 100644
--- a/src/pew/hittable/sphere.hpp
+++ b/src/pew/hittable/sphere.hpp
@@ -31,7 +31,7 @@ namespace pew {
 
             auto discriminant = half_b * half_b - a * c;
             if (discriminant < 0) return std::nullopt;
-            auto sqrtd = sqrt(discriminant);
+            auto sqrtd = discriminant.sqrt();
 
             // Find the nearest root that lies in the acceptable range.
             auto root = (-half_b - sqrtd) / a;
diff --git a/src/pew/hittable/triangle.hpp b/src/pew/hittable/triangle.hpp
index 59a51c2372a63d4b381f82414bdad2541b50a5a1..ec79d1e46dd418fb5702bc6932a50da93adbeb4b 100644
--- a/src/pew/hittable/triangle.hpp
+++ b/src/pew/hittable/triangle.hpp
@@ -34,14 +34,14 @@ namespace pew {
 
             if (a > -0.00001 && a < 0.00001) return std::nullopt;
 
-            auto f = 1 / a;
+            auto f = vec3(1) / a;
             auto s = ray.getOrigin() - v0.point;
-            auto u = f * s.dot(h);
+            double u = (f * s.dot(h)).x();
 
             if (u < 0.0 || u > 1.0) return std::nullopt;
 
             auto q = s.cross(e1);
-            auto v = f * ray.getDirection().dot(q);
+            double v = (f * ray.getDirection().dot(q)).x();
 
             if (v < 0.0 || u + v > 1.0) return std::nullopt;
 
diff --git a/src/pew/material/dielectric.hpp b/src/pew/material/dielectric.hpp
index 096e46a8a446fa10b8027bb77a6832607873a2a3..c7fc542859791c900d6e79e4fa5e89ff0264d77d 100644
--- a/src/pew/material/dielectric.hpp
+++ b/src/pew/material/dielectric.hpp
@@ -20,12 +20,12 @@ namespace pew::material {
             const double refractionRatio = hit.frontFace ? (1.0 / ir) : ir;
             const auto   unitDirection   = hit.ray.getDirection().unit();
 
-            const double cosTheta = fmin((-unitDirection).dot(hit.normal), 1.0);
-            const double sinTheta = sqrt(1.0 - cosTheta * cosTheta);
+            const auto cosTheta = std::min((-unitDirection).dot(hit.normal), vec3(1.0));
+            const auto sinTheta = (vec3(1.0) - cosTheta * cosTheta).sqrt();
 
             bool shouldReflect = refractionRatio * sinTheta > 1.0;
 
-            if (shouldReflect || reflectance(cosTheta, refractionRatio) > PCG_T::randomDouble()) {
+            if (shouldReflect || reflectance(cosTheta.x(), refractionRatio) > PCG_T::randomDouble()) {
                 return Ray(hit.getPoint(), unitDirection.reflect(hit.normal));
             } else {
                 return Ray(hit.getPoint(),
diff --git a/src/pew/ray/ray.hpp b/src/pew/ray/ray.hpp
index f54ee3604b43780d8f5880bc27d3ec24237246c4..0d2ff1733ca4a1db334642f7a2a4780361b96c67 100644
--- a/src/pew/ray/ray.hpp
+++ b/src/pew/ray/ray.hpp
@@ -38,7 +38,7 @@ namespace pew {
         }
         
         constexpr const auto
-        operator()(double t) const {
+        operator()(vec3 t) const {
             return origin + t * direction;
         }
     };
diff --git a/src/pew/test_int.cpp b/src/pew/test_int.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..72d5af01d336327ea6199d49ccf93c5ceb709604
--- /dev/null
+++ b/src/pew/test_int.cpp
@@ -0,0 +1,61 @@
+#include <pew/vector/intrinsic.hpp>
+using simd = pew::intrinsic<double>;
+
+int main() {
+{
+        auto reg1 = simd::set(-2, 20, 5);
+        auto reg2 = simd::set(11, 7, 36);
+        auto reg3 = simd::dot(reg1, reg2);
+        printf("reg1: %s\n", simd::dump(reg1).c_str());
+        printf("reg2: %s\n", simd::dump(reg2).c_str());
+        printf("reg3: %s\n", simd::dump(reg3).c_str());
+    }
+    {
+        auto reg1 = simd::set(-2, 54, -89);
+        auto reg2 = simd::abs(reg1);
+        printf("vec1: %s\n", simd::dump(reg2).c_str());
+    }
+    // {
+    //     auto reg1 = simd::set(0);
+    //     auto reg2 = simd::set(1, 0, 0);
+    //     auto reg3 = simd::set(0.000001);
+    //     printf("Reg1 ~~ 0: %i\n", simd::zero(reg1));
+    //     printf("Reg2 ~~ 0: %i\n", simd::zero(reg2));
+    //     printf("Reg3 ~~ 0: %i\n", simd::zero(reg3));
+    // }
+    // {
+    //     auto reg1 = simd::set(2);
+    //     auto reg2 = simd::set(2.000001);
+    //     auto reg3 = simd::set(1);
+    //     printf("Reg1 == Reg1: %i\n", simd::equal(reg1, reg1));
+    //     printf("Reg1 ~~ Reg1: %i\n", simd::equal(reg1, reg1));
+    //     printf("Reg1 == Reg2: %i\n", simd::equal(reg1, reg2));
+    //     printf("Reg1 ~~ Reg2: %i\n", simd::cmp(reg1, reg2));
+    //     printf("Reg1 == Reg3: %i\n", simd::equal(reg1, reg3));
+    //     printf("Reg1 ~~ Reg3: %i\n", simd::cmp(reg1, reg3));
+    // }
+    // {
+    //     auto reg1 = simd::set(200, 0, 0);
+    //     auto reg2 = simd::set(0, 1, 0);
+    //     auto reg3 = simd::cross(reg1, reg2);
+    //     printf("reg3: %s\n", simd::dump(reg3).c_str());
+    // }
+    {
+        auto reg1 = simd::set(3, 4, 20);
+        auto reg2 = simd::len(reg1);
+        // auto reg3 = simd::len3(reg1);
+        printf("reg2: %s\n", simd::dump(reg2).c_str());
+        // printf("reg3: %s\n", simd::dump(reg3).c_str());
+    }
+    // {
+    //     auto reg1 = simd::set(2, 30, 5.6);
+    //     auto reg2 = simd::rmax(reg1);
+    //     printf("reg2: %s\n", simd::dump(reg2).c_str());
+    // }
+    // {
+    //     auto reg1 = simd::set(2, 30, 5.6);
+    //     auto reg2 = simd::set(-2, 13, -3);
+    //     auto reg3 = simd::ang(reg1, reg2);
+    //     printf("reg3: %s\n", simd::dump(reg3).c_str());
+    // }
+}
\ No newline at end of file
diff --git a/src/pew/vector/intrinsic.hpp b/src/pew/vector/intrinsic.hpp
index 43c0dd3b82893d511827031cad61aa1698d5fb62..f274114bfe2fae1f51b932882711072e84e2906b 100644
--- a/src/pew/vector/intrinsic.hpp
+++ b/src/pew/vector/intrinsic.hpp
@@ -5,6 +5,7 @@
 #include <emmintrin.h>
 #include <immintrin.h>
 #include <smmintrin.h>
+#include <sstream>
 #include <tmmintrin.h>
 #include <xmmintrin.h>
 
@@ -20,7 +21,7 @@ namespace pew {
 
         inline static auto
         set(double x, double y, double z) {
-            return _mm256_set_pd(0, z, y, x);
+            return _mm256_set_pd(20, z, y, x);
         }
 
         inline static auto
@@ -29,8 +30,9 @@ namespace pew {
         }
 
         inline static auto
-        store(double* ptr, type vec) {
-            return _mm256_store_pd(ptr, vec);
+        store(double* arr, type vec) {
+            _mm256_store_pd(arr, vec);
+            return arr;
         }
 
         inline static auto
@@ -38,43 +40,84 @@ namespace pew {
             return _mm256_mul_pd(left, right);
         }
 
+        inline static auto
+        add(const type& left, const type& right) {
+            return _mm256_add_pd(left, right);
+        }
+
+        inline static auto
+        sub(const type& left, const type& right) {
+            return _mm256_sub_pd(left, right);
+        }
+
+        inline static auto
+        div(const type& left, const type& right) {
+            return _mm256_div_pd(left, right);
+        }
+
+        inline static auto
+        div(const type& vec, double s) {
+            return div(vec, set(s));
+        }
+
         inline static auto
         mul(const type& vec, double s) {
             return mul(vec, set(s));
         }
 
         inline static auto
-        sign(const type& vec) {
-            return mul(vec, set(-1.0));
+        dot(const type& left, const type& right) -> type {
+            auto reg_xy = mul(left, right);
+
+            // clang-format off
+            return _mm256_permute2f128_pd(
+                _mm256_permute_pd(
+                    add(
+                        _mm256_hadd_pd(reg_xy, reg_xy),
+                        _mm256_permute2f128_pd(reg_xy, reg_xy, 0x81)
+                    ), 0),
+                left, 0);
+            // clang-format on
         }
 
         inline static auto
-        add(const type& left, const type& right) {
-            return _mm256_add_pd(left, right);
+        sqrt(const type& reg) {
+            return _mm256_sqrt_pd(reg);
         }
 
         inline static auto
-        sub(const type& left, const type& right) {
-            return _mm256_sub_pd(left, right);
+        len(const type& reg) {
+            return sqrt(dot(reg, reg));
         }
 
-        inline static double
-        dot(const type& left, const type& right) {
-            alignas(alignment) double arr[4];
-            store(arr, mul(left, right));
-            return arr[0] + arr[1] + arr[2];
+        inline static auto
+        unit(const type& reg) {
+            return div(reg, len(reg));
         }
 
         inline static auto
-        div(const type& vec, double s) {
-            return _mm256_div_pd(vec, set(s));
+        sign(const type& vec) {
+            return mul(vec, set(-1.0));
+        }
+
+        inline static auto
+        cmplt(const type& left, const type& right) -> bool {
+            auto reg_ge = _mm256_cmp_pd(left, right, _CMP_GE_OQ);
+            return _mm256_testz_pd(reg_ge, reg_ge);
         }
 
         inline static auto
         abs(const type& abs) {
+            return _mm256_blendv_pd(abs, sign(abs), _mm256_cmp_pd(abs, set(0), _CMP_LT_OQ));
+        }
+
+        inline static auto
+        dump(const type& reg) {
             alignas(alignment) double arr[4];
-            store(arr, abs);
-            return set(std::abs(arr[0]), std::abs(arr[1]), std::abs(arr[2]));
+            store(arr, reg);
+            std::stringstream stream;
+            stream << "[" << arr[0] << " " << arr[1] << " " << arr[2] << " " << arr[3] << "]";
+            return stream.str();
         }
     };
 
diff --git a/src/pew/vector/vector.hpp b/src/pew/vector/vector.hpp
index 29a0f36b2f5086aabab6fa9d51bde04f532300af..2674dc7fa89f2353e485ad15b06ec2e1c3e42616 100644
--- a/src/pew/vector/vector.hpp
+++ b/src/pew/vector/vector.hpp
@@ -2,6 +2,7 @@
 #define PEW_VECTOR_VECTOR
 
 #include <cmath>
+#include <optional>
 #include <pew/random/random.hpp>
 #include <pew/tools/sqrt.hpp>
 #include <pew/vector/intrinsic.hpp>
@@ -33,7 +34,7 @@ namespace pew {
                 v = intrinsic<T>::set(x, y, z);
             }
         }
-        explicit constexpr Vector(T c) {
+        constexpr Vector(T c) {
             if (std::is_constant_evaluated()) {
                 v = std::array { c, c, c };
             } else {
@@ -256,12 +257,12 @@ namespace pew {
             }
         }
 
-        constexpr T
+        constexpr Vector
         dot(const Vector& other) const {
             if (std::is_constant_evaluated()) {
                 const auto& aRef = std::get<std::array<T, 3>>(v);
                 const auto& bRef = std::get<std::array<T, 3>>(other.v);
-                return aRef[0] * bRef[0] + aRef[1] * bRef[1] + aRef[2] * bRef[2];
+                return Vector { aRef[0] * bRef[0] + aRef[1] * bRef[1] + aRef[2] * bRef[2] };
             } else {
                 return intrinsic<T>::dot(asRegister(), other.asRegister());
             }
@@ -302,23 +303,49 @@ namespace pew {
             }
         }
 
-        constexpr T
+        constexpr Vector
+        operator/(const Vector& s) const {
+            if (std::is_constant_evaluated()) {
+                const auto& vecArr = std::get<std::array<T, 3>>(v);
+                const auto& vecS   = std::get<std::array<T, 3>>(s.v);
+                return { vecArr[0] / vecS[0], vecArr[1] / vecS[1], vecArr[2] / vecS[2] };
+            } else {
+                return intrinsic<T>::div(asRegister(), s.asRegister());
+            }
+        }
+
+        constexpr Vector&
+        operator/=(const Vector& s) {
+            if (std::is_constant_evaluated()) {
+                auto&       vecArr = std::get<std::array<T, 3>>(v);
+                const auto& vecS   = std::get<std::array<T, 3>>(s.v);
+                vecArr[0] /= vecS[0];
+                vecArr[1] /= vecS[1];
+                vecArr[2] /= vecS[2];
+                return *this;
+            } else {
+                v = intrinsic<T>::div(asRegister(), s.asRegister());
+                return *this;
+            }
+        }
+
+        constexpr Vector
         squareLength() const {
             return dot(*this);
         }
 
-        constexpr T
+        constexpr Vector
         length() const {
             if (std::is_constant_evaluated()) {
-                return pew::sqrt(squareLength());
+                return Vector { pew::sqrt(squareLength().x()) };
             } else {
-                return std::sqrt(squareLength());
+                return intrinsic<T>::len(asRegister());
             }
         }
 
         constexpr Vector
         unit() const {
-            return *this / this->length();
+            return intrinsic<T>::unit(asRegister());
         }
 
         constexpr Vector
@@ -327,11 +354,20 @@ namespace pew {
                 auto& vecArr = std::get<std::array<T, 3>>(v);
                 return { std::abs(vecArr[0]), std::abs(vecArr[1]), std::abs(vecArr[2]) };
             } else {
-                auto vecArr = clone();
                 return intrinsic<T>::abs(asRegister());
             }
         }
 
+        constexpr Vector
+        sqrt() const {
+            if (std::is_constant_evaluated()) {
+                auto& vecArr = std::get<std::array<T, 3>>(v);
+                return { std::sqrt(vecArr[0]), std::sqrt(vecArr[1]), std::sqrt(vecArr[2]) };
+            } else {
+                return intrinsic<T>::sqrt(asRegister());
+            }
+        }
+
         constexpr T
         max() const {
             if (std::is_constant_evaluated()) {
@@ -416,12 +452,19 @@ namespace pew {
             }
         }
 
+        inline friend int
+        operator<=>(const Vector& left, const Vector& right) {
+            return intrinsic<T>::cmplt(left.asRegister(), right.asRegister())
+                       ? -1
+                       : (intrinsic<T>::cmplt(right.asRegister(), left.asRegister()) ? 1 : 0);
+        }
+
         inline static Vector
         randomInUnitSphere() {
             if constexpr (std::is_floating_point_v<T>) {
                 while (true) {
                     auto p = random(-1.0, 1.0);
-                    if (p.squareLength() >= 1) continue;
+                    if (p.squareLength() >= Vector { 1 }) continue;
                     return p;
                 }
             } else {
@@ -434,7 +477,7 @@ namespace pew {
             while (true) {
                 auto p
                     = Vector(dyn_random::randomDouble(-1, 1), dyn_random::randomDouble(-1, 1), 0);
-                if (p.squareLength() >= 1) continue;
+                if (p.squareLength() >= Vector { 1 }) continue;
                 return p;
             }
         }
@@ -457,16 +500,16 @@ namespace pew {
         }
 
         constexpr Vector
-        refract(const Vector& normal,
-                double        etaiOverEtat,
-                const double* cosThetaRet = nullptr) const {
-            double cosTheta;
-            if (cosThetaRet)
-                cosTheta = *cosThetaRet;
-            else
-                cosTheta = fmin(dot(-normal), 1.0);
-            const auto rOutPerp     = etaiOverEtat * (*this + cosTheta * normal);
-            const auto rOutParallel = -sqrt(fabs(1.0 - rOutPerp.squareLength())) * normal;
+        refract(const Vector&          normal,
+                double                 etaiOverEtat,
+                std::optional<const Vector*> cosThetaRet = std::nullopt) const {
+            Vector cosTheta
+                = cosThetaRet ? *cosThetaRet.value() : std::min(dot(-normal), Vector { 1 });
+            const auto rOutPerp = etaiOverEtat * (*this + cosTheta * normal);
+            // const auto rOutParallel = -sqrt(fabs(Vector{1} - rOutPerp.squareLength())) * normal;
+            const auto rOutParallel
+                = (rOutPerp.squareLength() - Vector { 1 }).abs().sqrt() * -normal;
+
             return rOutPerp + rOutParallel;
         }
     };