muscle_tissue_martins.h

#pragma once

#include <funcy/examples/volumetric_penalty_functions.h>
#include <funcy/funcy.h>

namespace funcy
{
    namespace detail
    {
        template < linalg::Matrix Mat, int n = linalg::dim< Mat >() >
        auto generate_incompressible_muscle_tissue_martins( double c, double b, double d, double e,
                                                            const Mat& A, const Mat& F )
        {
            using namespace linalg;
            auto S = strain_tensor( F );
            auto si1 = mi1< Mat, n >( S() ) - n;
            auto si6 = mi6< Mat, Mat, n >( S(), A ) - 1;

            auto f = c * ( exp( b * si1 ) - 1 ) + d * ( exp( e * squared( si6 ) ) - 1 );
            return f( S );
        }
    } // namespace detail

    template < linalg::Matrix Mat, int offset = linalg::dim< Mat >() >
    auto incompressible_muscle_tissue_martins( double c, double b, double d, double e, const Mat& A,
                                               const Mat& F )
    {
        return finalize( detail::generate_incompressible_muscle_tissue_martins< Mat, offset >(
            c, b, d, e, A, F ) );
    }

    template < linalg::Matrix Mat, int offset = linalg::dim< Mat >() >
    auto incompressible_muscle_tissue_martins( const Mat& A, const Mat& F )
    {
        return incompressible_muscle_tissue_martins< Mat, offset >( 0.387, 23.46, 0.584, 12.43, A,
                                                                    F );
    }

    template < class Inflation, class Compression, linalg::Matrix Mat,
               int offset = linalg::dim< Mat >() >
    auto compressible_muscle_tissue_martins( double c, double b, double A, double a, double d0,
                                             double d1, const Mat& M, const Mat& F )
    {
        return finalize( detail::generate_incompressible_muscle_tissue_martins< Mat, offset >(
                             c, b, A, a, M, F ) +
                         volumetric_penalty< Inflation, Compression >( d0, d1, F ) );
    }

    template < class Inflation, class Compression, linalg::Matrix Mat,
               int offset = linalg::dim< Mat >() >
    auto compressible_muscle_tissue_martins( double d0, double d1, const Mat& M, const Mat& F )
    {
        return compressible_muscle_tissue_martins< Inflation, Compression, Mat, offset >(
            0.387, 23.46, 0.584, 12.43, d0, d1, M, F );
    }
} // namespace funcy