// #define FUNCTION_NAME <name-of-function-to-define>
// #define TRANS(x,y) x,y or y,x

#pragma GCC diagnostic ignored "-Wsequence-point"
template<class A,class B,
    class T=typename A::SCALAR,
    int am=(TRANS(A::fixed_n,A::fixed_m)),
    int an=(TRANS(A::fixed_m,A::fixed_n)),
    int bm=B::fixed_m,
    int m=am,int n=std::max(an,bm),
    class R=VEC_TYPE<T,mul_vec(A::type,B::type),m>>
std::enable_if_t<IS_MATRIX<A>::value&&IS_VECTOR<B>::value,R>
__attribute__((flatten))
FUNCTION_NAME(const A& a,const B& b)
{
    if constexpr(an<0||bm<0){assert((TRANS(a.Rows(),a.Columns()))==b.Size());}
    if constexpr(R::type==vec_zero) return R{};
    else if constexpr(n==0) return R{};
    else if constexpr(A::type==mat_id) return b;
    else if constexpr(A::type>=mat_diag) return a.x*b;
    else if constexpr(m>=0 && n>=0 && m<=size_limit && n<=size_limit)
    {
        return []<int... i,int... k>(const auto& a,const auto& b,std::integer_sequence<int,i...>,std::integer_sequence<int,k...>)
            {
                int j{};(void)j;
                return R{(j=i,((a(TRANS(j,k))*b(k))+...))...};
            }(a,b,std::make_integer_sequence<int,m>{},std::make_integer_sequence<int,n>{});
    }
    else
    {
        int mm=(TRANS(a.Columns(),a.Rows()));
        int nn=(TRANS(a.Rows(),a.Columns()));
        R r((INITIAL_SIZE)mm);
        for(int i=0;i<mm;i++)
        {
            T z{};
            for(int j=0;j<nn;j++)
                z+=a(TRANS(i,j))*b(j);
            r(i)=z;
        }
        return r;
    }
}
#undef FUNCTION_NAME
#undef TRANS
