diff --git a/2019/src/day16.cpp b/2019/src/day16.cpp
index 2ead298..22ce998 100644
--- a/2019/src/day16.cpp
+++ b/2019/src/day16.cpp
@@ -3,7 +3,9 @@
 #include <vector>
 #include <cassert>
 #include <iterator>
+#include <numeric>
 #include "days.hpp"
+#include "utils.hpp"
 
 namespace {
     std::array<int, 4> base_pattern{0, 1, 0, -1};
@@ -26,47 +28,38 @@ namespace {
         return result;
     }
 
-    std::vector<int> partial_sum(const std::vector<int> &numbers) {
+    void simulate(std::vector<int> numbers, std::ostream &output) {
+        std::vector<int> new_numbers(numbers.size());
         std::vector<int> partial_sums(numbers.size());
-        partial_sums[0] = numbers[0];
-        for (int j = 1; j < numbers.size(); ++j) {
-            partial_sums[j] = partial_sums[j - 1] + numbers[j];
-        }
 
-        return partial_sums;
-    }
-
-    void next_phase(const std::vector<int>& numbers, std::vector<int> &new_numbers) {
-        new_numbers.resize(numbers.size());
-
-        for (int rank = 0; rank < numbers.size(); ++rank) {
-            auto partial_sums = partial_sum(numbers);
-            int n = 0;
-            for (int pos = rank; pos < numbers.size(); pos += rank + 1) {
-                int run = std::min(rank + 1, (int) numbers.size() - pos);
-                if (int modifier = get_modifier(rank, pos); modifier) {
-                    n += modifier * (partial_sums[pos + run - 1] - partial_sums[pos] + numbers[pos]);
+        for (int i = 0; i < 100; ++i) {
+            for (int rank = 0; rank < numbers.size(); ++rank) {
+                std::partial_sum(numbers.begin() + rank, numbers.end(), partial_sums.begin() + rank);
+                int n = 0;
+                for (int pos = rank; pos < numbers.size(); pos += rank + 1) {
+                    int run = std::min(rank + 1, (int) numbers.size() - pos);
+                    if (int modifier = get_modifier(rank, pos); modifier) {
+                        n += modifier * (partial_sums[pos + run - 1] - partial_sums[pos] + numbers[pos]);
+                    }
                 }
+
+                n = std::abs(n % 10);
+
+                new_numbers[rank] = n;
             }
 
-            n = std::abs(n % 10);
-
-            new_numbers[rank] = n;
+            std::swap(numbers, new_numbers);
         }
+
+        std::copy(numbers.begin(), numbers.begin() + 8, std::ostream_iterator<int>(output));
+        output << std::endl;
     }
 }
 
 void aoc2019::day16_part1(std::istream &input, std::ostream &output) {
     auto numbers = read_input(input);
-    auto new_numbers = numbers;
 
-    for (int i = 0; i < 100; ++i) {
-        next_phase(numbers, new_numbers);
-        std::swap(numbers, new_numbers);
-    }
-
-    std::copy(numbers.begin(), numbers.begin() + 8, std::ostream_iterator<int>(output));
-    output << std::endl;
+    simulate(std::move(numbers), output);
 }
 
 void aoc2019::day16_part2(std::istream &input, std::ostream &output) {
@@ -89,9 +82,10 @@ void aoc2019::day16_part2(std::istream &input, std::ostream &output) {
 
     numbers = std::vector(numbers.begin() + offset, numbers.end());
     std::vector<int> new_numbers(numbers.size());
+    std::vector<int> partial_sums(numbers.size());
 
     for (int i = 0; i < 100; ++i) {
-        std::vector<int> partial_sums = partial_sum(numbers);
+        std::partial_sum(numbers.begin(), numbers.end(), partial_sums.begin());
 
         for (int j = 0; j < numbers.size(); ++j) {
             int n = partial_sums.back() - partial_sums[j] + numbers[j];