这个东西没啥技术含量,所以就直接上代码了:
#ifndef _MAT_HPP_
#define _MAT_HPP_
#include <memory>
#include <climits>
#include <float.h>
#include <map>
#include <iostream>
#include <iomanip>
#include "ht_memory.h"
#ifdef WITH_BOOST
#include <boost/pool/pool.hpp>
#endif
template<typename val_t>
val_t max_and_swap(const val_t& v1, const val_t& v2)
{
return v1 < v2 ? v2 : v1;
}
template<typename type>
struct destoryer
{
static void do_destory(type& v)
{
//printf("%s\r\n", typeid(type).name());
}
};
template<int i_size, typename val_t>
struct mat_m
{
#ifdef WITH_BOOST
static boost::pool<> s_pool;
#else
val_t sz_ele[i_size];
#endif
val_t* p;
mat_m() :p(nullptr)
{
#ifndef WITH_BOOST
//p = (val_t*)malloc(i_size * sizeof(val_t));
p = sz_ele;
#else
p = (val_t*)(s_pool.malloc());
#endif
for (int i = 0; i < i_size; ++i)
{
new(p + i) val_t(0);
//p[i] = 0;
}
}
~mat_m()
{
if (p)
{
#ifndef WITH_BOOST
/*
for (int i = 0; i < i_size; ++i)
destoryer<val_t>::do_destory(p[i]);
free(p);
*/
#else
for (int i = 0; i < i_size; ++i)
destoryer<val_t>::do_destory(p[i]);
//p[i].~val_t();
s_pool.free(p);
#endif
p = 0;
}
}
val_t& get(const int& len_1d, const int& i_1d_idx, const int& i_2d_idx)
{
val_t& ret = p[i_2d_idx + len_1d * i_1d_idx];
return ret;
}
val_t max_abs() const
{
double d = -1*DBL_MAX;
for (int i = 0; i < i_size; ++i)
{
d = d < abs(p[i]) ? abs(p[i]) : d;
}
return d;
}
val_t max() const
{
val_t d = -1 * DBL_MAX;
for (int i = 0; i < i_size; ++i)
{
//d = d < (p[i]) ? (p[i]) : d;
d = max_and_swap(p[i], d);
}
return d;
}
val_t sum() const
{
val_t d_sum = 0.;
for (int i = 0; i < i_size; ++i)
{
d_sum = d_sum + p[i];
}
return d_sum;
}
template<int len_1d, int i_1d_idx, int i_2d_idx>
inline val_t& get_val()
{
static_assert((i_2d_idx + len_1d * i_1d_idx) < i_size, "ERROR:mat_m over flow!!!");
return p[i_2d_idx + len_1d * i_1d_idx];
}
template<int len_1d, int i_1d_idx, int i_2d_idx>
inline val_t get_val() const
{
return p[i_2d_idx + len_1d * i_1d_idx];
}
};
#ifdef WITH_BOOST
template<int i_size, typename val_t>
boost::pool<> mat_m<i_size, val_t>::s_pool(i_size * sizeof(val_t));
#endif
template<int row_num, int col_num, typename val_t = double>
struct mat
{
using t_type = mat<col_num, row_num, val_t>;
using type = val_t;
typedef val_t vt;
static constexpr int r = row_num;
static constexpr int c = col_num;
using mat_m_t = mat_m<row_num * col_num, val_t>;
std::shared_ptr<mat_m_t> pval;
bool b_t;
mat():b_t(false)
{
pval = std::make_shared<mat_m_t>();
}
mat(const mat<row_num, col_num, val_t>& other) :pval(other.pval), b_t(other.b_t)
{
}
~mat()
{
pval.reset();
}
mat(const val_t&& v):b_t(false)
{
pval = std::make_shared<mat_m_t>();
for (int i = 0; i < row_num; ++i)
{
for (int j = 0; j < col_num; ++j)
{
pval->get(col_num, i, j) = v;
}
}
}
mat(const val_t& v) :b_t(false)
{
pval = std::make_shared<mat_m_t>();
for (int i = 0; i < row_num; ++i)
{
for (int j = 0; j < col_num; ++j)
{
pval->get(col_num, i, j) = v;
}
}
}
#if 0
template<typename val_other_t>
mat(const val_other_t& v) :b_t(false)
{
pval = std::make_shared<mat_m_t>();
for (int i = 0; i < row_num; ++i)
{
for (int j = 0; j < col_num; ++j)
{
pval->get(col_num, i, j) = static_cast<val_t>(v);
}
}
}
#endif
mat(const std::initializer_list<val_t>& lst):b_t(false)
{
pval = std::make_shared<mat_m_t>();
auto itr = lst.begin();
for (int i = 0; i < row_num; ++i)
{
for (int j = 0; j < col_num; ++j)
{
if (itr == lst.end())return;
pval->get(col_num, i, j) = *itr;
itr++;
}
}
}
val_t& get(const int& i_row, const int& i_col)
{
if (!b_t)
return pval->get(col_num, i_row, i_col);
else
return pval->get(row_num, i_col, i_row);
}
val_t get(const int& i_row, const int& i_col) const
{
if (!b_t)
return pval->get(col_num, i_row, i_col);
else
return pval->get(row_num, i_col, i_row);
}
template<int i_1d_idx, int i_2d_idx>
inline val_t& get_val()
{
if (!b_t)
return pval->template get_val<col_num, i_1d_idx, i_2d_idx>();
else
return pval->template get_val<row_num, i_2d_idx, i_1d_idx>();
}
template<int i_1d_idx, int i_2d_idx>
inline val_t get_val() const
{
static_assert(i_1d_idx < row_num && i_2d_idx < col_num, "ERROR: mat::get_val overflow!!!!!");
if (!b_t)
return pval->template get_val<col_num, i_1d_idx, i_2d_idx>();
else
return pval->template get_val<row_num, i_2d_idx, i_1d_idx>();
}
mat<col_num, row_num, val_t> t() const
{
mat<col_num, row_num, val_t> ret;
ret.pval = pval;
ret.b_t = !b_t;
return ret;
}
val_t max_abs() const
{
return pval->max_abs();
}
val_t max() const
{
return pval->max();
}
val_t sum() const
{
return pval->sum();
}
void print() const
{
std::cout << "[" << std::endl;
for (int i = 0; i < row_num; ++i)
{
std::cout << std::setw(3) << "[";
for (int j = 0; j < col_num; ++j)
{
std::cout << (j != 0 ? "," : "") << std::setw(10) << get(i, j);
}
std::cout << std::setw(3) << "]" << std::endl;
}
std::cout << "]" << std::endl;
}
template<int other_col_num>
mat<row_num, other_col_num, val_t> dot(const mat<col_num, other_col_num, val_t>& mt) const
{
return dot<row_num, col_num, col_num, other_col_num, val_t>(*this, mt);
}
mat<row_num, col_num, val_t> rot180() const
{
mat<row_num, col_num, val_t> ret;
for (int r = 0; r < row_num; ++r)
{
for (int c = 0; c < col_num; ++c)
{
ret.get(r, c) = get(row_num-1-r, col_num-1-c);
}
}
return ret;
}
template<int row_base, int col_base, int row_num_other, int col_num_other>
void assign(const mat<row_num_other, col_num_other, val_t>& mt_other)
{
/* 这里不麻烦了,直接写成运行时 */
for (int r = 0; r < row_num_other; ++r)
{
for (int c = 0; c < col_num_other; ++c)
{
if (r + row_base < 0 || c + col_base < 0)
{
continue;
}
if (r + row_base >= row_num || c + col_base >= col_num)
{
break;
}
get(r + row_base, c + col_base) = mt_other.get(r, c);
}
}
}
template<int top_pad, int left_pad, int right_pad, int bottom_pad>
mat<row_num + top_pad + bottom_pad, col_num + left_pad + right_pad, val_t>
pad() const
{
using mat_ret_t = mat<row_num + top_pad + bottom_pad, col_num + left_pad + right_pad, val_t>;
mat_ret_t mt_ret;
mt_ret.template assign<top_pad, left_pad>(*this);
return mt_ret;
}
template<int row_span, int col_span>
mat<row_num + row_span*(row_num - 1), col_num + col_span*(col_num-1)>
span() const
{
using mat_ret_t = mat<row_num + row_span * (row_num - 1), col_num + col_span * (col_num - 1)>;
mat_ret_t mt_ret;
for (int r = 0; r < row_num; ++r)
{
for (int c = 0; c < col_num; ++c)
{
mt_ret.get(r*(row_span + 1), c*(col_span + 1)) = get(r, c);
}
}
return mt_ret;
}
template<int row_base, int col_base, int row_len, int col_len>
val_t region_max(int& i_row, int& i_col) const
{
static_assert(row_base < row_num && col_base < col_num, "region_max overflow!!!");
val_t d_max = -1. * DBL_MAX;
for (int r = row_base; r < row_base + row_len && r < row_num; ++r)
{
for (int c = col_base; c < col_base + col_len && c < col_num; ++c)
{
if (d_max < get(r, c))
{
i_row = r, i_col = c;
d_max = get(r, c);
}
}
}
return d_max;
}
mat<row_num*col_num, 1, val_t> one_col() const
{
mat<row_num*col_num, 1, val_t> ret;
ret.pval = pval;
return ret;
}
template<typename t>
static void print_sub_type(t)
{
std::cout << typeid(t).name();
}
template<int r, int c, typename t>
static void print_sub_type(mat<r, c, t>)
{
mat<r, c, t>::print_type();
}
static void print_type()
{
printf("<matrix %d * %d type: ", row_num, col_num);
print_sub_type<val_t>(val_t());
printf(">\r\n");
}
const val_t& operator[](const int& idx) const
{
return pval->p[idx];
}
val_t& operator[](const int& idx)
{
return pval->p[idx];
}
friend std::ostream& operator<<(std::ostream& ofs, const mat<row_num, col_num, val_t>& mt)
{
std::cout << "[" ;
for (int i = 0; i < row_num; ++i)
{
std::cout << std::setw(3) << "[";
for (int j = 0; j < col_num; ++j)
{
std::cout << (j != 0 ? "," : "") << mt.get(i, j);
}
std::cout << std::setw(3) << "]";
}
std::cout << "]";
return ofs;
}
int size() const
{
return row_num * col_num;
}
mat<row_num - 1, col_num - 1, val_t> algebraic_complement_val(const int& r, const int& c) const
{
mat<row_num - 1, col_num - 1, val_t> mt_ret;
int retr = 0, retc = 0;
for (int mr = 0; mr < row_num; mr++)
{
if (mr == r) continue;
for (int mc = 0; mc < col_num; mc++)
{
if (mc == c) continue;
mt_ret.get(retr, retc++) = get(mr, mc);
}
retc = 0;
retr++;
}
return mt_ret;
}
};
template<typename type>
struct mat_size
{
static constexpr int num = 1;
};
template<int target_col, int target_row, typename target_val_t>
struct mat_size<mat<target_col, target_row, target_val_t> >
{
static constexpr int num = target_col * target_row * mat_size<target_val_t>::num;
};
template<typename unite_type>
struct mat_unite_type
{
using type = unite_type;
};
template<int row_num, int col_num, typename val_t>
struct mat_unite_type<mat<row_num, col_num, val_t> >
{
using type = typename mat_unite_type<val_t>::type;
};
template<typename mat_t>
struct unite_mat
{
using type = mat<mat_size<mat_t>::num, 1, typename mat_unite_type<mat_t>::type>;
};
template<typename mat_t>
using unite_mat_t = typename unite_mat<mat_t>::type;
template<typename target_type>
struct stretch_to_unite
{
template<typename type1>
static target_type cal(type1& v)
{
auto k = v.one_col();
return stretch_to_unite<target_type>::cal(k);
}
static target_type cal(target_type& v)
{
return v;
}
};
template<int i1, int i2, typename val_t>
mat<i1, i2, val_t> max_and_swap(const mat<i1, i2, val_t>& mt_max, const mat<i1, i2, val_t>& mt_optional)
{
using ret_type = mat<i1, i2, val_t>;
ret_type ret;
for (int i = 0; i < i1; ++i)
{
for (int j = 0; j < i2; ++j)
{
ret.get(i, j) = max_and_swap(mt_max.get(i,j), mt_optional.get(i, j));
}
}
return ret;
}
template<typename val_t>
val_t max_and_choose(const val_t& v1, const val_t& v2, const val_t& v3, const val_t& v4)
{
return v1 < v2 ? v3 : v4;
}
template<int i1, int i2, typename val_t>
mat<i1, i2, val_t> max_and_choose(const mat<i1, i2, val_t>& mt_judge1, const mat<i1, i2, val_t>& mt_judge2, const mat<i1, i2, val_t>& mt_optional1, const mat<i1, i2, val_t>& mt_optional2)
{
using ret_type = mat<i1, i2, val_t>;
ret_type ret;
for (int i = 0; i < i1; ++i)
{
for (int j = 0; j < i2; ++j)
{
ret.get(i, j) = max_and_choose(mt_judge1.get(i, j), mt_judge2.get(i, j), mt_optional1.get(i, j), mt_optional2.get(i, j));
}
}
return ret;
}
template<int row_num, int col_num, typename val_t>
struct destoryer<mat<row_num, col_num, val_t> >
{
using type = mat<row_num, col_num, val_t>;
static void do_destory(type& v)
{
v.~type();
}
};
#endif
带逆矩阵的基础运算:
#ifndef _BASE_FUNCTION_HPP_
#define _BASE_FUNCTION_HPP_
#include <algorithm>
#include <math.h>
#include "base_logic.hpp"
// * 计算梯度
template<typename func_t>
auto derivative(func_t&& f, const decltype(f(0))& v)
{
constexpr double SMALL_VAL = 1e-11;
return (f(v + SMALL_VAL) - f(v - SMALL_VAL)) / (2.*SMALL_VAL);
}
/* 点乘运算 */
template<int r1, int c1, int r2, int c2, typename imatt1, typename imatt2, typename vt = double>
inline typename imatt1::type n_dot(const imatt1& mt1, const imatt2& mt2)
{
static_assert(c1 == r2, "[matrix dot error]\tleft matrix column number do not match right matrix's row number.");
if constexpr (c1 != 0 || r2 != 0)
{
return mt1.template get_val<r1, c1>() * mt2.template get_val<r2, c2>() + n_dot<r1, c1 - 1, r2 - 1, c2>(mt1, mt2);
}
if constexpr (c1 == 0 && r2 == 0)
{
return mt1.template get_val<r1, c1>() * mt2.template get_val<r2, c2>();
}
}
template<int r, int c>
class v_dot
{
public:
template<typename imatt1, typename imatt2, typename vt = double>
static typename imatt1::type cal(const imatt1& mt1, const imatt2& mt2)
{
return n_dot<r, mt1.c - 1, mt2.r - 1, c>(mt1, mt2);
}
};
template<int row_num1, int col_num1, int row_num2, int col_num2, typename val_t = double>
mat<row_num1, col_num2, val_t> dot(const mat<row_num1, col_num1, val_t>& mt1, const mat<row_num2, col_num2, val_t>& mt2)
{
using omatt = mat<row_num1, col_num2, val_t>;
using imatt1 = mat<row_num1, col_num1, val_t>;
using imatt2 = mat<row_num2, col_num2, val_t>;
omatt mt_ret;
col_loop<col_num2 - 1, v_dot>(mt_ret, mt1, mt2);
return mt_ret;
}
/* 加法运算 */
template<int r, int c>
class v_add
{
public:
template<typename imatt1, typename imatt2>
static typename imatt1::type cal(const imatt1& mt1, const imatt2& mt2)
{
return mt1.template get_val<r, c>() + mt2.template get_val<r, c>();
}
};
template<int r, int c>
class n_add
{
public:
template<typename imatt2, typename vt = double>
static vt cal(const vt& mt1, const imatt2& mt2)
{
return mt1 + mt2.template get_val<r, c>();
}
};
template<int r, int c>
struct c_add
{
template<int row_num, int col_num, typename vt>
static vt cal(const mat<row_num, col_num, vt>& mt, const mat<row_num, 1, vt>& v)
{
return mt.template get_val<r, c>() + v.template get_val<r, 0>();
}
};
template<int r, int c>
struct r_add
{
template<int row_num, int col_num, typename vt>
static vt cal(const mat<row_num, col_num, vt>& mt, const mat<1, col_num, vt>& v)
{
return mt.template get_val<r, c>() + v.template get_val<0, c>();
}
};
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator+(const mat<row_num, col_num, val_t>& mt1, const mat<row_num, col_num, val_t>& mt2)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, v_add>(mt_ret, mt1, mt2);
return mt_ret;
}
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator+(const val_t& v, const mat<row_num, col_num, val_t>& mt)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_add>(mt_ret, v, mt);
return mt_ret;
}
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator+(const mat<row_num, col_num, val_t>& mt, const val_t& v)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_add>(mt_ret, v, mt);
return mt_ret;
}
#if 0
template<int row_num, int col_num, typename val_t_other, typename val_t = double>
mat<row_num, col_num, val_t> operator+(const mat<row_num, col_num, val_t>& mt, const val_t_other& v)
{
return mt + static_cast<val_t>(v);
}
template<int row_num, int col_num, typename val_t_other, typename val_t = double>
mat<row_num, col_num, val_t> operator+(const val_t_other& v, const mat<row_num, col_num, val_t>& mt)
{
return mt + static_cast<val_t>(v);
}
#endif
template<int row_num, int col_num, typename val_t = double>
void spread_add(mat<row_num, col_num, val_t>& mt_ret, const mat<row_num, col_num, val_t>& mt, const mat<row_num, 1, val_t>& v)
{
col_loop<col_num - 1, c_add>(mt_ret, mt, v);
}
template<int row_num, int col_num, typename val_t = double>
void spread_add(mat<row_num, col_num, val_t>& mt_ret, const mat<row_num, col_num, val_t>& mt, const mat<1, col_num, val_t>& v)
{
col_loop<col_num - 1, r_add>(mt_ret, mt, v);
}
template<int row_num, int col_num, typename val_t = double>
void spread_add(mat<row_num, col_num, val_t>& mt_ret, const mat<row_num, col_num, val_t>& mt, const mat<1, 1, val_t>& v)
{
col_loop<col_num - 1, n_add>(mt_ret, v.template get_val<0, 0>(), mt);
}
template<int row_num, int col_num, typename val_t = double>
void spread_add(mat<row_num, col_num, val_t>& mt_ret, const mat<row_num, 1, val_t>& v, const mat<row_num, col_num, val_t>& mt)
{
col_loop<col_num - 1, c_add>(mt_ret, mt, v);
}
template<int row_num, int col_num, typename val_t = double>
void spread_add(mat<row_num, col_num, val_t>& mt_ret, const mat<1, col_num, val_t>& v, const mat<row_num, col_num, val_t>& mt)
{
col_loop<col_num - 1, r_add>(mt_ret, mt, v);
}
template<int row_num, int col_num, typename val_t = double>
void spread_add(mat<row_num, col_num, val_t>& mt_ret, const mat<1, 1, val_t>& v, const mat<row_num, col_num, val_t>& mt)
{
col_loop<col_num - 1, n_add>(mt_ret, v.template get_val<0, 0>(), mt);
}
template<typename val_t = double>
void spread_add(mat<1, 1, val_t>& mt_ret, const mat<1, 1, val_t>& v, const mat<1, 1, val_t>& mt)
{
col_loop<0, n_add>(mt_ret, v.template get_val<0, 0>(), mt.template get_val<0, 0>);
}
/* 减法运算 */
template<int r, int c>
class v_minus
{
public:
template<typename imatt1, typename imatt2, typename vt = double>
static typename imatt1::type cal(const imatt1& mt1, const imatt2& mt2)
{
return mt1.template get_val<r, c>() - mt2.template get_val<r, c>();
}
};
template<int r, int c>
class n_minus
{
public:
template<int mat_r, int mat_c, typename vt = double>
static vt cal(const vt& v, const mat<mat_r, mat_c, vt>& mt2)
{
return v - mt2.template get_val<r, c>();
}
template<int mat_r, int mat_c, typename vt = double>
static vt cal(const mat<mat_r, mat_c, vt>& mt2, const vt& v)
{
return mt2.template get_val<r, c>() - v;
}
};
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator-(const mat<row_num, col_num, val_t>& mt1, const mat<row_num, col_num, val_t>& mt2)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, v_minus>(mt_ret, mt1, mt2);
return mt_ret;
}
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator-(const val_t& v, const mat<row_num, col_num, val_t>& mt)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_minus>(mt_ret, v, mt);
return mt_ret;
}
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator-(const mat<row_num, col_num, val_t>& mt, const val_t& v)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_minus>(mt_ret, mt, v);
return mt_ret;
}
/*
template<int row_num, int col_num, typename val_t_other, typename val_t = double>
mat<row_num, col_num, val_t> operator-(const mat<row_num, col_num, val_t>& mt, const val_t_other& v)
{
return mt - val_t(v);
}
template<int row_num, int col_num, typename val_t_other, typename val_t = double>
mat<row_num, col_num, val_t> operator-(const val_t_other& v, const mat<row_num, col_num, val_t>& mt)
{
return (val_t(v) - mt);
}
*/
/* 乘法运算 */
template<int r, int c>
class n_mul
{
public:
template<typename imatt2>
static typename imatt2::type cal(const typename imatt2::type& mt1, const imatt2& mt2)
{
return mt1 * mt2.template get_val<r, c>();
}
};
template<int r, int c>
class v_mul
{
public:
template<typename imatt1, typename imatt2>
static typename imatt1::type cal(const imatt1& mt1, const imatt2& mt2)
{
return mt1.template get_val<r, c>() * mt2.template get_val<r, c>();
}
};
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator*(const val_t& v, const mat<row_num, col_num, val_t>& mt)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_mul>(mt_ret, v, mt);
return mt_ret;
}
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator*(const mat<row_num, col_num, val_t>& mt, const val_t& v)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_mul>(mt_ret, v, mt);
return mt_ret;
}
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> operator*(const mat<row_num, col_num, val_t>& mt1, const mat<row_num, col_num, val_t>& mt2)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, v_mul>(mt_ret, mt1, mt2);
return mt_ret;
}
/* 除法 */
template<int r, int c>
class n_div
{
public:
template<int row_num, int col_num, typename vt = double>
static vt cal(const mat<row_num, col_num, vt>& mt, const vt& v)
{
return mt.template get_val<r, c>() / v;
}
template<int row_num, int col_num, typename vt = double>
static vt cal(const vt& v, const mat<row_num, col_num, vt>& mt)
{
return v / mt.template get_val<r, c>();
}
};
template<int r, int c>
class v_div
{
public:
template<int row_num, int col_num, typename vt = double>
static vt cal(const mat<row_num, col_num, vt>& mt1, const mat<row_num, col_num, vt>& mt2)
{
return mt1.template get_val<r, c>() / mt2.template get_val<r, c>();
}
};
template<int row_num, int col_num, typename val_t>
mat<row_num, col_num, val_t> operator/(const mat<row_num, col_num, val_t>& mt, const val_t& v)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_div>(mt_ret, mt, v);
return mt_ret;
}
template<int row_num, int col_num, typename val_t>
mat<row_num, col_num, val_t> operator/(const val_t& v, const mat<row_num, col_num, val_t>& mt)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_div>(mt_ret, v, mt);
return mt_ret;
}
template<int row_num, int col_num, typename vt>
mat<row_num, col_num, vt> operator/(const mat<row_num, col_num, vt>& mt1, const mat<row_num, col_num, vt>& mt2)
{
using omatt = mat<row_num, col_num, vt>;
omatt mt_ret;
col_loop<col_num - 1, v_div>(mt_ret, mt1, mt2);
return mt_ret;
}
double sqrtl(double d)
{
return sqrt(d);
}
template<int i1, int i2, typename val_t>
mat<i1, i2, val_t> sqrtl(const mat<i1, i2, val_t>& mt)
{
using type = mat<i1, i2, val_t>;
type ret;
for (int i = 0; i < i1; ++i)
{
for (int j = 0; j < i2; ++j)
{
ret.get(i, j) = sqrtl(mt.get(i,j));
}
}
return ret;
}
template<int r, int c>
class n_sqrt
{
public:
template<typename imatt>
static typename imatt::type cal(const imatt& mt)
{
return sqrtl(mt.template get_val<r, c>());
}
};
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> sqrtm(const mat<row_num, col_num, val_t>& mt)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_sqrt>(mt_ret, mt);
return mt_ret;
}
/* exp运算 */
template<int r, int c>
struct n_exp
{
template<typename imatt>
static typename imatt::type cal(const imatt& mt )
{
return exp(mt.template get_val<r, c>());
}
};
template<int row_num, int col_num, typename val_t = double>
mat<row_num, col_num, val_t> exp(const mat<row_num, col_num, val_t>& mt)
{
using omatt = mat<row_num, col_num, val_t>;
omatt mt_ret;
col_loop<col_num - 1, n_exp>(mt_ret, mt);
return mt_ret;
}
/* 卷积运算 */
template<int row_base, int col_base, int row_delta, int col_delta, typename imat_origin, typename imat_tpl>
inline auto col_loop_mul(const imat_origin& mt_origin, const imat_tpl& mt_tpl)
{
if constexpr (col_delta != 0)
{
return mt_origin.template get_val<row_base + row_delta, col_base + col_delta>() * mt_tpl.template get_val<row_delta, col_delta>()
+ col_loop_mul<row_base, col_base, row_delta, col_delta - 1, imat_origin, imat_tpl>(mt_origin, mt_tpl);
}
if constexpr (col_delta == 0)
{
return mt_origin.template get_val<row_base + row_delta, col_base + col_delta>() * mt_tpl.template get_val<row_delta, col_delta>();
}
}
template<int row_base, int col_base, int row_delta, int col_delta, typename imat_origin, typename imat_tpl>
inline auto row_loop_add(const imat_origin& mt_origin, const imat_tpl& mt_tpl)
{
if constexpr (row_delta != 0)
{
return col_loop_mul<row_base, col_base, row_delta, col_delta>(mt_origin, mt_tpl)
+ col_loop_mul<row_base, col_base, row_delta - 1, col_delta>(mt_origin, mt_tpl);
}
if constexpr (row_delta == 0)
{
return col_loop_mul<row_base, col_base, row_delta, col_delta>(mt_origin, mt_tpl);
}
}
template<int r, int c>
struct v_inner_conv
{
template<typename imat_origin_t, typename imat_tpl_t>
inline static auto cal(const imat_origin_t& mt_origin, const imat_tpl_t& mt_tpl)
{
return row_loop_add<r, c, imat_tpl_t::r - 1, imat_tpl_t::c - 1>(mt_origin, mt_tpl);
}
};
constexpr int get_step_inner_size(int i_origin, int i_tpl, int i_step)
{
return (i_origin - i_tpl) / i_step + 1;
}
class def_config;
class same_pad;
template<typename op_t = def_config>
constexpr int get_pad_size(int i_origin, int i_tpl, int i_step)
{
return (((i_origin - i_tpl) / i_step) + (((i_origin - i_tpl) % i_step) == 0 ? 0 : 1))*i_step - (i_origin - i_tpl);
}
template<>
constexpr int get_pad_size<same_pad>(int i_origin, int i_tpl, int i_step)
{
//static_assert(i_step == 1, "ERROR:get_pad_size same pad step must equal 1.");
return i_tpl - 1; // a + x - b + 1 = a; x = b-1;
}
constexpr int get_ceil_div(int i_origin, int i_tpl)
{
return (i_origin / i_tpl + ((i_origin%i_tpl)==0?0:1));
}
template<int input_row, int intput_col, int tpl_row, int tpl_col, int row_step, int col_step, typename op_t = def_config>
struct pad_size_t
{
static constexpr int top = get_pad_size<op_t>(input_row, tpl_row, row_step) / 2; // 上方的扩充数量
static constexpr int left = get_pad_size<op_t>(intput_col, tpl_col, col_step) / 2; // 左方的扩充数量
static constexpr int right = get_pad_size<op_t>(intput_col, tpl_col, col_step) - left; // 右方的扩充数量
static constexpr int bottom = get_pad_size<op_t>(input_row, tpl_row, row_step) - top; // 下方的扩充数量
};
// * 内部卷积
template<int row_step, int col_step, int row_num, int col_num, int tpl_row, int tpl_col, typename val_t>
inline mat<get_step_inner_size(row_num, tpl_row, row_step), get_step_inner_size(col_num, tpl_col, col_step), val_t>
inner_conv(const mat<row_num, col_num, val_t>& mt_origin, const mat<tpl_row, tpl_col, val_t>& mt_tpl)
{
using ret_type = mat<get_step_inner_size(row_num, tpl_row, row_step), get_step_inner_size(col_num, tpl_col, col_step), val_t>;
ret_type mt_ret;
col_loop<ret_type::c - 1, v_inner_conv>(mt_ret, mt_origin, mt_tpl);
return mt_ret;
}
template<typename mat_t, typename ...mat_ts>
struct st_one_col
{
static constexpr int all_size = (mat_t::r*mat_t::c) + st_one_col<mat_ts...>::all_size;
};
template<typename mat_t>
struct st_one_col<mat_t>
{
static constexpr int all_size = (mat_t::r*mat_t::c);
};
template<typename mat_t, typename ...mat_ts>
void concat_mat(typename mat_t::type* p, const mat_t& mt, const mat_ts... mts)
{
constexpr int cpy_size = mat_t::r*mat_t::c;
//memcpy(p, mt.pval->p, cpy_size*sizeof(mat_t::type));
for (int i = 0; i < cpy_size; ++i)
{
p[i] = mt.pval->p[i];
}
if constexpr(0!=sizeof...(mat_ts))
concat_mat(p + cpy_size, mts...);
}
template<typename mat_t, typename ...mat_ts>
mat<st_one_col<mat_t, mat_ts...>::all_size, 1> stretch_one_col(const mat_t& mt, const mat_ts&...mts)
{
using ret_type = mat<st_one_col<mat_t, mat_ts...>::all_size, 1>;
ret_type ret;
concat_mat(ret.pval->p, mt, mts...);
return ret;
}
template<typename mat_t, typename ...mat_ts>
void split_mat(typename mat_t::type* p, const mat_t& mt, const mat_ts... mts)
{
constexpr int cpy_size = mat_t::r*mat_t::c;
//memcpy(mt.pval->p, p, cpy_size * sizeof(mat_t::type));
for (int i = 0; i < cpy_size; ++i)
{
mt.pval->p[i] = p[i];
}
if constexpr (0 != sizeof...(mat_ts))
split_mat(p + cpy_size, mts...);
}
template<typename mat_t, typename ...mat_ts>
void split_one_mat(const mat_t& mt, const mat_ts&...mts)
{
split_mat(mt.pval->p, mts...);
}
template<int i1, int i2, typename val_t>
mat<i1, i2, val_t> abs(const mat<i1, i2, val_t>& mt)
{
mat<i1, i2, val_t> ret;
for (int i = 0; i < i1; ++i)
{
for (int j = 0; j < i2; ++j)
{
ret.get(i, j) = abs(mt.get(i, j));
}
}
return ret;
}
#include <vector>
/* 输出改成0均值1均方差的 */
template<int row_num, int col_num, typename val_t>
inline std::vector<mat<row_num, col_num, val_t> > normalize(const std::vector<mat<row_num, col_num, val_t> >& vec_input, mat<row_num, col_num, val_t>& mt_mean, mat<row_num, col_num, val_t>& mt_div)
{
if (1 == vec_input.size())return vec_input;
using type = mat<row_num, col_num, val_t>;
std::vector<type> vec_ret;
//type mt_mean, mt_div;
for (int i = 0; i < vec_input.size(); ++i)
{
mt_mean = mt_mean + vec_input[i] / static_cast<val_t>(vec_input.size());
}
for (int i = 0; i < vec_input.size(); ++i)
{
auto delta = vec_input[i] - mt_mean;
mt_div = mt_div + delta * delta / static_cast<val_t>(vec_input.size());
}
auto mt_s = sqrtl(mt_div);
for (int i = 0; i < vec_input.size(); ++i)
{
vec_ret.push_back((vec_input[i] - mt_mean) / mt_s);
}
return vec_ret;
}
#include "ht_memory.h"
template<typename val_t>
void write_file(const val_t& vt, ht_memory& mry)
{
mry << vt;
}
template<int row_num, int col_num, typename val_t>
void write_file(const mat<row_num, col_num, val_t>& mt, ht_memory& mry)
{
for (int r = 0; r < row_num; ++r)
{
for (int c = 0; c < col_num; ++c)
{
write_file(mt.get(r, c), mry);
}
}
}
template<typename val_t>
void read_file(ht_memory& mry, val_t& vt)
{
mry >> vt;
}
template<int row_num, int col_num, typename val_t>
void read_file(ht_memory& mry, mat<row_num, col_num, val_t>& mt)
{
for (int r = 0; r < row_num; ++r)
{
for (int c = 0; c < col_num; ++c)
{
read_file(mry, mt.get(r, c));
}
}
}
template<int N>
double sub_mul_cal(const mat<N, N, double>& mt, const std::vector<int>& vec_cur)
{
double dmul = 1.;
for (size_t siz_add_itr = 0; siz_add_itr < vec_cur.size(); ++siz_add_itr)
{
dmul *= mt.get(siz_add_itr, vec_cur[siz_add_itr]);
//printf("%4.2lf ", mt.get(siz_add_itr, vec_cur[siz_add_itr]));
}
//printf("=%4.2lf", dmul);
return dmul;
}
template<int N>
void det_cal(double& dret, const mat<N, N, double>& mt, const std::vector<int>& vec_idx, const size_t& siz_cur, const double& dflag = 1.)
{
//printf("\n%4.2lf ", dflag);
dret += (sub_mul_cal(mt, vec_idx)*dflag); // * 把自己算进去
for (size_t siz_itr = siz_cur; siz_itr < vec_idx.size() - 1 ; ++siz_itr)
{
std::vector<int> vec_cur(vec_idx);
double dcurflag = dflag;
for (size_t siz_in_itr = siz_itr + 1; siz_in_itr < vec_cur.size(); ++siz_in_itr)
{
std::swap(vec_cur[siz_itr], vec_cur[siz_in_itr]);
dcurflag *= -1.;
//dret += (sub_mul_cal(mt, vec_cur)*dflag);
det_cal(dret, mt, vec_cur, siz_itr + 1, dcurflag);
}
}
}
template<int N>
double det(const mat<N, N, double>& mt)
{
std::vector<int> vec(N, 0);
int idx = 0;
std::generate(vec.begin(), vec.end(), [&idx]()
{
return idx++;
});
double dret = 0.;
det_cal(dret, mt, vec, 0, 1.);
return dret;
}
template<int N, typename val_t>
mat<N, N, val_t> algebraic_complement(const mat<N, N, val_t>& mt)
{
mat<N, N, val_t> mtret;
double drflag = 1.;
for (int i = 0; i < N; ++i)
{
double dcflag = 1.;
for (int j = 0; j < N; ++j)
{
mtret.get(i, j) = (drflag * dcflag * det(mt.algebraic_complement_val(i, j)));
dcflag *= -1.;
}
drflag *= -1.;
}
return mtret.t();
}
template<int N, typename val_t>
mat<N, N, val_t> inverse(const mat<N, N, val_t>& mt)
{
return algebraic_complement(mt)/det(mt);
}
#endif