Polygon Mesh Processing page 74

`E_LSCM` = ∑_T A_T‖M_T v_T - [0 -1
                              1  0] M_T u_T‖²
where
 
v_T ∈ ℝ^3
u_T ∈ ℝ^3
M_T ∈ ℝ^(2×3)
A_T ∈ ℝ

/*
`E_LSCM` = ∑_T A_T‖M_T v_T - [0 -1
                              1  0] M_T u_T‖²
where
 
v_T ∈ ℝ^3
u_T ∈ ℝ^3
M_T ∈ ℝ^(2×3)
A_T ∈ ℝ
*/
#include <Eigen/Core>
#include <Eigen/Dense>
#include <Eigen/Sparse>
#include <iostream>
#include <set>

struct pmp_74ResultType {
    double E_LSCM;
    pmp_74ResultType(const double & E_LSCM)
    : E_LSCM(E_LSCM)
    {}
};

pmp_74ResultType pmp_74(
    const std::vector<Eigen::Matrix<double, 3, 1>> & v,
    const std::vector<Eigen::Matrix<double, 3, 1>> & u,
    const std::vector<Eigen::Matrix<double, 2, 3>> & M,
    const std::vector<double> & A)
{
    const long dim_0 = A.size();
    assert( v.size() == dim_0 );
    assert( u.size() == dim_0 );
    assert( M.size() == dim_0 );

    double sum_0 = 0;
    for(int T=1; T<=M.size(); T++){
        Eigen::Matrix<double, 2, 2> E_LSCM_0;
        E_LSCM_0 << 0, -1,
        1, 0;
        sum_0 += A.at(T-1) * pow((M.at(T-1) * v.at(T-1) - E_LSCM_0 * M.at(T-1) * u.at(T-1)).lpNorm<2>(), 2);
    }
    double E_LSCM = sum_0;

    return pmp_74ResultType(E_LSCM);
}


void generateRandomData(std::vector<Eigen::Matrix<double, 3, 1>> & v,
    std::vector<Eigen::Matrix<double, 3, 1>> & u,
    std::vector<Eigen::Matrix<double, 2, 3>> & M,
    std::vector<double> & A)
{
    const int dim_0 = rand()%10;
    v.resize(dim_0);
    for(int i=0; i<dim_0; i++){
        v[i] = Eigen::VectorXd::Random(3);
    }
    u.resize(dim_0);
    for(int i=0; i<dim_0; i++){
        u[i] = Eigen::VectorXd::Random(3);
    }
    M.resize(dim_0);
    for(int i=0; i<dim_0; i++){
        M[i] = Eigen::MatrixXd::Random(2, 3);
    }
    A.resize(dim_0);
    for(int i=0; i<dim_0; i++){
        A[i] = rand() % 10;
    }
}


int main(int argc, char *argv[])
{
    srand((int)time(NULL));
    std::vector<Eigen::Matrix<double, 3, 1>> v;
    std::vector<Eigen::Matrix<double, 3, 1>> u;
    std::vector<Eigen::Matrix<double, 2, 3>> M;
    std::vector<double> A;
    generateRandomData(v, u, M, A);
    pmp_74ResultType func_value = pmp_74(v, u, M, A);
    std::cout<<"return value:\n"<<func_value.E_LSCM<<std::endl;
    return 0;
}

"""
`E_LSCM` = ∑_T A_T‖M_T v_T - [0 -1
                              1  0] M_T u_T‖²
where
 
v_T ∈ ℝ^3
u_T ∈ ℝ^3
M_T ∈ ℝ^(2×3)
A_T ∈ ℝ
"""
import numpy as np
import scipy
import scipy.linalg
from scipy import sparse
from scipy.integrate import quad
from scipy.optimize import minimize


class pmp_74ResultType:
    def __init__( self, E_LSCM):
        self.E_LSCM = E_LSCM


def pmp_74(v, u, M, A):
    v = np.asarray(v, dtype=np.float64)
    u = np.asarray(u, dtype=np.float64)
    M = np.asarray(M, dtype=np.float64)
    A = np.asarray(A, dtype=np.float64)

    dim_0 = A.shape[0]
    assert v.shape == (dim_0, 3, )
    assert u.shape == (dim_0, 3, )
    assert M.shape == (dim_0, 2, 3)
    assert A.shape == (dim_0,)

    sum_0 = 0
    for T in range(1, len(M)+1):
        E_LSCM_0 = np.zeros((2, 2))
        E_LSCM_0[0] = [0, -1]
        E_LSCM_0[1] = [1, 0]
        sum_0 += A[T-1] * np.power(np.linalg.norm(M[T-1] @ v[T-1] - E_LSCM_0 @ M[T-1] @ u[T-1], 2), 2)
    E_LSCM = sum_0
    return pmp_74ResultType(E_LSCM)


def generateRandomData():
    dim_0 = np.random.randint(10)
    v = np.random.randn(dim_0, 3, )
    u = np.random.randn(dim_0, 3, )
    M = np.random.randn(dim_0, 2, 3)
    A = np.random.randn(dim_0)
    return v, u, M, A


if __name__ == '__main__':
    v, u, M, A = generateRandomData()
    print("v:", v)
    print("u:", u)
    print("M:", M)
    print("A:", A)
    func_value = pmp_74(v, u, M, A)
    print("return value: ", func_value.E_LSCM)

function output = pmp_74(v, u, M, A)
% output = pmp_74(v, u, M, A)
%
%    `E_LSCM` = ∑_T A_T‖M_T v_T - [0 -1
%                                  1  0] M_T u_T‖²
%    where
%     
%    v_T ∈ ℝ^3
%    u_T ∈ ℝ^3
%    M_T ∈ ℝ^(2×3)
%    A_T ∈ ℝ
    if nargin==0
        warning('generating random input data');
        [v, u, M, A] = generateRandomData();
    end
    function [v, u, M, A] = generateRandomData()
        dim_0 = randi(10);
        v = randn(dim_0,3);
        u = randn(dim_0,3);
        M = randn(dim_0,2,3);
        A = randn(dim_0,1);
    end

    A = reshape(A,[],1);

    dim_0 = size(A, 1);
    assert( isequal(size(v), [dim_0, 3]) );
    assert( isequal(size(u), [dim_0, 3]) );
    assert( isequal(size(M), [dim_0, 2, 3]) );
    assert( size(A,1) == dim_0 );

    sum_0 = 0;
    for T = 1:size(M, 1)
        E_LSCM_0 = zeros(2, 2);
        E_LSCM_0(1,:) = [0, -1];
        E_LSCM_0(2,:) = [1, 0];
        sum_0 = sum_0 + A(T) * norm(squeeze(M(T,:,:)) * v(T,:)' - E_LSCM_0 * squeeze(M(T,:,:)) * u(T,:)', 2).^2;
    end
    E_LSCM = sum_0;
    output.E_LSCM = E_LSCM;
end

\documentclass[12pt]{article}
\usepackage{mathdots}
\usepackage[bb=boondox]{mathalfa}
\usepackage{mathtools}
\usepackage{amssymb}
\usepackage{libertine}
\DeclareMathOperator*{\argmax}{arg\,max}
\DeclareMathOperator*{\argmin}{arg\,min}
\usepackage[paperheight=8in,paperwidth=4in,margin=.3in,heightrounded]{geometry}
\let\originalleft\left
\let\originalright\right
\renewcommand{\left}{\mathopen{}\mathclose\bgroup\originalleft}
\renewcommand{\right}{\aftergroup\egroup\originalright}
\begin{document}

\begin{center}
\resizebox{\textwidth}{!} 
{
\begin{minipage}[c]{\textwidth}
\begin{align*}
\textit{E\_LSCM} & = \sum_\mathit{T} \mathit{A}_{ \mathit{T} }\left\|\mathit{M}_{ \mathit{T} }\mathit{v}_{ \mathit{T} } - \begin{bmatrix}
0 & -1\\
1 & 0\\
\end{bmatrix}\mathit{M}_{ \mathit{T} }\mathit{u}_{ \mathit{T} }\right\|_2^{2} \\
\intertext{where} 
\mathit{v}_{\mathit{T}} & \in \mathbb{R}^{ 3} \\
\mathit{u}_{\mathit{T}} & \in \mathbb{R}^{ 3} \\
\mathit{M}_{\mathit{T}} & \in \mathbb{R}^{ 2 \times 3 } \\
\mathit{A}_{\mathit{T}} & \in \mathbb{R} \\
\\
\end{align*}
\end{minipage}
}
\end{center}

\end{document}