diff options
Diffstat (limited to 'p-spin.hpp')
| -rw-r--r-- | p-spin.hpp | 23 |
1 files changed, 23 insertions, 0 deletions
@@ -121,6 +121,29 @@ Tensor<Scalar, 3> ddzDot(const Vector<Scalar>& z, const Vector<Scalar>& dH) { } template <class Scalar> +Tensor<Scalar, 3> dcdzDot(const Vector<Scalar>& z, const Vector<Scalar>& dH, const Matrix<Scalar>& ddH) { + Tensor<Scalar, 1> zT = Eigen::TensorMap<Tensor<const Scalar, 1>>(z.data(), {z.size()}); + Tensor<Scalar, 1> dHT = Eigen::TensorMap<Tensor<const Scalar, 1>>(dH.data(), {dH.size()}); + Tensor<Scalar, 2> ddHT = Eigen::TensorMap<Tensor<const Scalar, 2>>(ddH.data(), {dH.size(), dH.size()}); + + Matrix<Scalar> I = Matrix<Real>::Identity(z.size(), z.size()); + Tensor<Scalar, 2> IT = Eigen::TensorMap<Tensor<const Scalar, 2>>(I.data(), {z.size(), z.size()}); + + Eigen::array<Eigen::IndexPair<int>, 0> ei = {}; + + Scalar z² = z.squaredNorm(); + + return ddHT.conjugate().contract(zT.conjugate(), ei) / z² + +IT.contract(dHT.conjugate(), ei).shuffle((std::array<int, 3>){0,2,1}) / z² + -zT.contract(dHT.conjugate(), ei).contract(zT.conjugate(), ei) / pow(z², 2) + -ddHT.conjugate().contract(zT, ip00).contract(zT.conjugate(), ei).contract(zT.conjugate(), ei) / pow(z², 2) + -IT.contract(zT.conjugate(), ei) * (dH.dot(z) / pow(z², 2)) + -IT.contract(zT.conjugate(), ei).shuffle((std::array<int, 3>){0,2,1}) * (dH.dot(z) / pow(z², 2)) + +zT.contract(zT.conjugate(), ei).contract(zT.conjugate(), ei) * (2.0 * dH.dot(z) / pow(z², 3)) + ; +} + +template <class Scalar> Tensor<Scalar, 3> ddzDotConjugate(const Vector<Scalar>& z, const Vector<Scalar>& dH, const Matrix<Scalar>& ddH, const Tensor<Scalar, 3>& dddH) { Tensor<Scalar, 1> zT = Eigen::TensorMap<Tensor<const Scalar, 1>>(z.data(), {z.size()}); Tensor<Scalar, 1> dHT = Eigen::TensorMap<Tensor<const Scalar, 1>>(dH.data(), {dH.size()}); |