//#####################################################################
// Copyright 2004-2007, Geoffrey Irving, Igor Neverov, Eftychios Sifakis, Joseph Teran.
// This file is part of PhysBAM whose distribution is governed by the license contained in the accompanying file PHYSBAM_COPYRIGHT.txt.
//#####################################################################
// Class GENERAL_STRESS_DERIVATIVE
//#####################################################################
#ifndef __GENERAL_STRESS_DERIVATIVE__
#define __GENERAL_STRESS_DERIVATIVE__

#include <Core/Matrices/MATRIX.h>
#include <Core/Matrices/SYMMETRIC_MATRIX.h>
#include <Core/Vectors/VECTOR.h>
namespace PhysBAM{
template<class TV>
class GENERAL_STRESS_DERIVATIVE
{
    typedef typename TV::SCALAR T;
public:
    // Aijkl = A(i,j)(k,l) = Aklij
    MATRIX<MATRIX<T,TV::m>,TV::m> A;

//#####################################################################
    MATRIX<T,TV::m> Contract_03(const MATRIX<T,TV::m>& M) const; // solving Cjk = Mijkl * Mil
    MATRIX<T,TV::m> Contract_13(const MATRIX<T,TV::m>& M) const; // solving Cik = Mijkl * Mjl
    MATRIX<T,TV::m> Contract_23(const MATRIX<T,TV::m>& M) const; // solving Cij = Mijkl * Mkl
    MATRIX<T,TV::m> Differential(const MATRIX<T,TV::m>& M) const // solving Cij = Mijkl * Mkl
    {return Contract_23(M);}
//#####################################################################
};
// Compute Sigma, given SVD of F (U*F*Transpose(V))
template<class T,int d> MATRIX<T,d>
Sigma(const MATRIX<T,d>& P,const MATRIX<T,d>& F);
// Compute dSigma_dF:M, given SVD of F (U*F*Transpose(V))
template<class TV,class T=typename TV::SCALAR,int d=TV::m> MATRIX<T,d>
dSigma_dF_Times(const GENERAL_STRESS_DERIVATIVE<TV>& dP_dF,
    const MATRIX<T,TV::m>& P,const MATRIX<T,TV::m>& F,const MATRIX<T,TV::m>& M);
}
#endif
