4 files changed, 79 insertions, 8 deletions
diff --git a/libprakpp/include/prakcommon.hpp b/libprakpp/include/prakcommon.hpp
index 074e290..dcafecf 100644
--- a/libprakpp/include/prakcommon.hpp
+++ b/libprakpp/include/prakcommon.hpp
@@ -1,5 +1,7 @@
 #pragma once
 
+#include <cstddef>
+#include <cstdint>
 #include <cmath>
 #include <functional>
 #include <vector>
@@ -132,3 +134,7 @@ template <typename T> struct pvalue<T> operator-(const T &a, const struct pvalue
 template <typename T> std::ostream &operator<<(std::ostream &os, const struct pvalue<T> &p) { return os << p.val << "±" << p.err; }
 
 } // namespace prak
+
+namespace std {
+template <typename T> struct prak::pvalue<T> abs(const struct prak::pvalue<T> &p) { return {abs(p.val), abs(p.err)}; }
+}
diff --git a/libprakpp/include/prakmath.hpp b/libprakpp/include/prakmath.hpp
index bb09e74..d43c11f 100644
--- a/libprakpp/include/prakmath.hpp
+++ b/libprakpp/include/prakmath.hpp
@@ -321,7 +321,6 @@ void least_squares_linear(
 			  << ", s.size() = " << ss.size() << std::endl;
 		return;
 	}
-	[[assume(xs.size() == ys.size() && ys.size() == ss.size())]];
 	size_t sz = xs.size();
 
 	prak::vector<T> ones(sz);
@@ -364,6 +363,27 @@ void least_squares_linear(
 	b.err = sqrt(d1 / D);
 }
 
+template <std::floating_point T>
+void least_squares_prop(
+	const prak::vector<T> &xs, 
+	const prak::vector<T> &ys, 
+	const prak::vector<T> &ss, 
+	struct pvalue<T> &a) {
+	assert("Different-sized vectors passed to least squares method" && xs.size() == ys.size());
+	size_t s = xs.size();
+	prak::vector<T> buf(s), si(s);
+	vmul(ss.data(), ss.data(), si.data(), s);
+	
+	vmul(xs.data(), ys.data(), buf.data(), s);
+	vdiv(buf.data(), ss.data(), buf.data(), s);
+	T D1 = vsum(buf.data(), s);
+	
+	vmul(xs.data(), xs.data(), buf.data(), s);
+	vdiv(buf.data(), ss.data(), buf.data(), s);
+	T D2 = vsum(buf.data(), s);
+	a.val = D1/D2;
+	a.err = 1/D2;
+}
 
 /// May throw std::bad_optional_access
 template <typename T>
diff --git a/libprakpp/include/prakmatrix.hpp b/libprakpp/include/prakmatrix.hpp
index f1ccddc..2f87401 100644
--- a/libprakpp/include/prakmatrix.hpp
+++ b/libprakpp/include/prakmatrix.hpp
@@ -2,6 +2,7 @@
 
 #include "prakcommon.hpp"
 
+#include <iomanip>
 #include <cstring>
 #include <ostream>
 #include <optional>
@@ -75,7 +76,7 @@ struct matrix : public std::vector<T> {
 	}
 	/// multiply matrix by a scalar
 	template <typename Arg>
-	std::enable_if_t<std::is_arithmetic_v<Arg>, struct matrix<T>>::type
+	typename std::enable_if_t<std::is_arithmetic_v<Arg>, struct matrix<T>>::type
 	operator*(const Arg other) const {
 		struct matrix<T> ret(rows, cols);
 		for (size_t i = 0; i < this->size(); ++i)
@@ -161,12 +162,27 @@ struct matrix : public std::vector<T> {
 		}
 		return ret;
 	}
+	
+	void add_row(size_t index, std::vector<T> row_data) {
+		++rows;
+		this->resize(rows * cols);
+		memmove(this->data() + cols * (index+1), this->data() + cols * index, (rows - index - 1) * cols * sizeof (T));
+		for (size_t i = 0; i < cols; ++i)
+			this->operator[](index * cols + i) = row_data[i];
+	}
+	void add_column(size_t index, std::vector<T> column_data) {
+		this->reserve(rows * (cols + 1));
+		for (size_t i = 0; i < rows; ++i)
+			this->insert(this->begin() + (cols + 1) * i + index, column_data[i]);
+		++cols;
+	}
 
 	std::ostream &print(std::ostream &out = std::cout, size_t width = 10) const {
+		std::cout.setf(std::ios::left);
 		for (size_t row = 0; row < rows; ++row) {
 			for (size_t col = 0; col < cols; ++col)
-				out << "|" << std::setw(width) << SUBSCR_OPRTR(row, col); 
-			std::cout << "|\n";
+				out /*<< "|"*/ << std::setw(width) << SUBSCR_OPRTR(row, col); 
+			std::cout << "\n";
 		}
 		return out;
 	}
diff --git a/libprakpp/include/praktable.hpp b/libprakpp/include/praktable.hpp
index 1870ec6..aabcb8e 100644
--- a/libprakpp/include/praktable.hpp
+++ b/libprakpp/include/praktable.hpp
@@ -105,10 +105,10 @@ public:
 		iterator operator--(int) { iterator ret = *this; --(*this); return ret; }
 		iterator &operator+(int x) { data_index += columns * x; return *this; }
 		iterator &operator-(int x) { data_index -= columns * x; return *this; }
-		std::strong_ordering operator<=>(const iterator& other) {
+		std::strong_ordering operator<=>(const iterator& other) const {
 			return parent == other.parent && col_index == other.col_index && data_index <=> other.data_index;
 		}
-		bool operator==(const iterator& other) {
+		bool operator==(const iterator& other) const {
 			return parent == other.parent && col_index == other.col_index && data_index == other.data_index;
 		}
 		reference operator*() { return parent->data[data_index]; };
@@ -270,6 +270,18 @@ public:
 
 	// apply a function to several columns and store result in another column
 	// function must accept std::vector or arguments
+	table &apply(function_type function, std::initializer_list<std::string> _args, std::optional<std::string> result) {
+		std::vector<std::string> args = _args;
+		size_t result_index = result.has_value() ? index(*result) : 0;
+		for (size_t i = 0; i < rows; ++i) {
+			std::vector<dtype> v(args.size());
+			for (size_t j = 0; j < args.size(); ++j) 
+				v[j] = SUBSCR_OPRTR(args[j], i);
+			if (result.has_value()) data[columns * i + result_index] = function(v);
+			else (void)function(v);
+		}
+		return *this;
+	}
 	table &apply(function_type function, stringvec args, std::optional<std::string> result) {
 		size_t result_index = result.has_value() ? index(*result) : 0;
 		for (size_t i = 0; i < rows; ++i) {
@@ -281,7 +293,6 @@ public:
 		}
 		return *this;
 	}
-
 	table& apply(function_type function, const std::string& arg, std::optional<std::string> result) {
 		size_t result_index = result.has_value() ? index(*result) : 0;
 		for (size_t i = 0; i < rows; ++i) {
@@ -411,11 +422,12 @@ public:
 
 
 	// Appends a row to the table. if name is set, appends it to `opt_rownames`
-	void add_row(std::vector<dtype> values, std::optional<std::string> name = std::nullopt) {
+	table& add_row(std::vector<dtype> values, std::optional<std::string> name = std::nullopt) {
 		if (values.size() == 0) values = std::vector<dtype>(columns, dtype{});
 		data.resize(columns * (++rows));
 		std::copy_n(values.cbegin(), columns, data.end() - columns);
 		if (name.has_value()) opt_rownames.push_back(*name);
+		return *this;
 	}
 
 
@@ -518,6 +530,23 @@ public:
 		prak::least_squares_linear<dtype>(_x, _y, _s, ret.first, ret.second);
 		return ret;
 	}
+	prak::pvalue<dtype>
+	least_squares_prop(std::string x, std::string y, std::optional<std::string> sigma, std::optional<dtype> sigma_fixed) const {
+		if (sigma.has_value() == sigma_fixed.has_value()) 
+			throw std::invalid_argument("sigma and sigma_fixed can't both have (no) value");
+		prak::vector<dtype> _x(rows);
+		prak::vector<dtype> _y(rows);
+		prak::vector<dtype> _s(rows);
+
+		std::copy(cbegin(x), cend(x), _x.begin());
+		std::copy(cbegin(y), cend(y), _y.begin());
+		if (sigma.has_value()) std::copy(cbegin(*sigma), cend(*sigma), _s.begin());
+		else _s = prak::vector<dtype>(rows, static_cast<dtype>(*sigma_fixed));
+		
+		prak::pvalue<dtype> ret;
+		prak::least_squares_prop(_x, _y, _s, ret);
+		return ret;
+	}
 
 	// calculate an average of the column 
 	dtype col_avg(const std::string &column) const {