An example from Anisotropic Elasticity for Inversion-Safety and Element Rehabilitation Eq. 47
The original equation:
I❤️LA implementation:
from linearalgebra: tr
`J₃` = 1₃,₃
`κ_angle(Dₘ)` = 3(√2 v)^(2/3)(7/4‖`Dₘ`‖_F^2-1/4tr(`J₃``Dₘ`ᵀ`Dₘ`))⁻¹
where
`Dₘ` ∈ ℝ^(3×3)
v ∈ ℝ
I❤️LA compiled to C++/Eigen:
/*
from linearalgebra: tr
`J₃` = 1₃,₃
`κ_angle(Dₘ)` = 3(√2 v)^(2/3)(7/4‖`Dₘ`‖_F^2-1/4tr(`J₃``Dₘ`ᵀ`Dₘ`))⁻¹
where
`Dₘ` ∈ ℝ^(3×3)
v ∈ ℝ
*/
#include <Eigen/Core>
#include <Eigen/Dense>
#include <Eigen/Sparse>
#include <iostream>
#include <set>
struct anisotropic_elasticity_47ResultType {
Eigen::Matrix<double, 3, 3> J₃;
double κ_angle_left_parenthesis_Dₘ_right_parenthesis;
anisotropic_elasticity_47ResultType(const Eigen::Matrix<double, 3, 3> & J₃,
const double & κ_angle_left_parenthesis_Dₘ_right_parenthesis)
: J₃(J₃),
κ_angle_left_parenthesis_Dₘ_right_parenthesis(κ_angle_left_parenthesis_Dₘ_right_parenthesis)
{}
};
anisotropic_elasticity_47ResultType anisotropic_elasticity_47(
const Eigen::Matrix<double, 3, 3> & Dₘ,
const double & v)
{
Eigen::Matrix<double, 3, 3> J₃ = Eigen::MatrixXd::Ones(3, 3);
double κ_angle_left_parenthesis_Dₘ_right_parenthesis = 3 * pow((sqrt(2) * v), (2 / double(3))) * 1 / ((7 / double(4) * pow((Dₘ).norm(), 2) - 1 / double(4) * (J₃ * Dₘ.transpose() * Dₘ).trace()));
return anisotropic_elasticity_47ResultType(J₃, κ_angle_left_parenthesis_Dₘ_right_parenthesis);
}
void generateRandomData(Eigen::Matrix<double, 3, 3> & Dₘ,
double & v)
{
v = rand() % 10;
Dₘ = Eigen::MatrixXd::Random(3, 3);
}
int main(int argc, char *argv[])
{
srand((int)time(NULL));
Eigen::Matrix<double, 3, 3> Dₘ;
double v;
generateRandomData(Dₘ, v);
anisotropic_elasticity_47ResultType func_value = anisotropic_elasticity_47(Dₘ, v);
std::cout<<"return value:\n"<<func_value.κ_angle_left_parenthesis_Dₘ_right_parenthesis<<std::endl;
return 0;
}
I❤️LA compiled to Python/NumPy/SciPy:
"""
from linearalgebra: tr
`J₃` = 1₃,₃
`κ_angle(Dₘ)` = 3(√2 v)^(2/3)(7/4‖`Dₘ`‖_F^2-1/4tr(`J₃``Dₘ`ᵀ`Dₘ`))⁻¹
where
`Dₘ` ∈ ℝ^(3×3)
v ∈ ℝ
"""
import numpy as np
import scipy
import scipy.linalg
from scipy import sparse
from scipy.integrate import quad
from scipy.optimize import minimize
class anisotropic_elasticity_47ResultType:
def __init__( self, J3, κ_angle_left_parenthesis_Dₘ_right_parenthesis):
self.J3 = J3
self.κ_angle_left_parenthesis_Dₘ_right_parenthesis = κ_angle_left_parenthesis_Dₘ_right_parenthesis
def anisotropic_elasticity_47(Dₘ, v):
Dₘ = np.asarray(Dₘ, dtype=np.float64)
assert Dₘ.shape == (3, 3)
assert np.ndim(v) == 0
J3 = np.ones((3, 3))
κ_angle_left_parenthesis_Dₘ_right_parenthesis = 3 * np.power((np.sqrt(2) * v), (2 / 3)) * 1 / ((7 / 4 * np.power(np.linalg.norm(Dₘ, 'fro'), 2) - 1 / 4 * np.trace(J3 @ Dₘ.T @ Dₘ)))
return anisotropic_elasticity_47ResultType(J3, κ_angle_left_parenthesis_Dₘ_right_parenthesis)
def generateRandomData():
v = np.random.randn()
Dₘ = np.random.randn(3, 3)
return Dₘ, v
if __name__ == '__main__':
Dₘ, v = generateRandomData()
print("Dₘ:", Dₘ)
print("v:", v)
func_value = anisotropic_elasticity_47(Dₘ, v)
print("return value: ", func_value.κ_angle_left_parenthesis_Dₘ_right_parenthesis)
I❤️LA compiled to MATLAB:
function output = anisotropic_elasticity_47(D_m, v)
% output = anisotropic_elasticity_47(`Dₘ`, v)
%
% from linearalgebra: tr
%
% `J₃` = 1₃,₃
% `κ_angle(Dₘ)` = 3(√2 v)^(2/3)(7/4‖`Dₘ`‖_F^2-1/4tr(`J₃``Dₘ`ᵀ`Dₘ`))⁻¹
%
% where
%
% `Dₘ` ∈ ℝ^(3×3)
% v ∈ ℝ
if nargin==0
warning('generating random input data');
[D_m, v] = generateRandomData();
end
function [D_m, v] = generateRandomData()
v = randn();
D_m = randn(3, 3);
end
assert( isequal(size(D_m), [3, 3]) );
assert(numel(v) == 1);
J3 = ones(3, 3);
kappa_angle_D_m = 3 * (sqrt(2) * v).^(2 / 3) * 1 / ((7 / 4 * norm(D_m, 'fro').^2 - 1 / 4 * trace(J3 * D_m' * D_m)));
output.J3 = J3;
output.kappa_angle_D_m = kappa_angle_D_m;
end
I❤️LA compiled to LaTeX:
\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*}
\text{from linearalgebra import tr}\\
\textit{J₃} & = \mathbb{ 1 }_{ 3,3 } \\
\textit{κ\_angle(Dₘ)} & = 3{\left( \sqrt{2}\mathit{v} \right)}^{\frac{2}{3}}\left( \frac{7}{4}\left\|\textit{Dₘ}\right\|_F^{2} - \frac{1}{4}tr\left( \textit{J₃}{\textit{Dₘ}}^T\textit{Dₘ} \right) \right)^{-1} \\
\intertext{where}
\textit{Dₘ} & \in \mathbb{R}^{ 3 \times 3 } \\
\mathit{v} & \in \mathbb{R} \\
\\
\end{align*}
\end{minipage}
}
\end{center}
\end{document}
I❤️LA LaTeX output: