// This file is part of GOAST, a C++ library for variational methods in Geometry Processing
//
// Copyright (C) 2021 Behrend Heeren & Josua Sassen, University of Bonn <goast@ins.uni-bonn.de>
//
// This Source Code Form is subject to the terms of the Mozilla
// Public License v. 2.0. If a copy of the MPL was not distributed
// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
#include <shell-energy/membrane_energy.h>
#include <shell-energy/bending_energy.h>
#include <shell-energy/shell_energy.h>

#include <igl/read_triangle_mesh.h>
#include <igl/edge_flaps.h>

int main( int argc, char *argv[] ) {
  if (argc != 3) {
    std::cerr << "Usage: evaluateEnergy <mesh1> <mesh2>" << std::endl;
    return -1;
  }
  // Data storage of meshes
  Eigen::MatrixXd V1, V2;
  Eigen::MatrixXi F1, F2;

  // Load the two meshes
  igl::read_triangle_mesh( argv[1], V1, F1 );
  igl::read_triangle_mesh( argv[2], V2, F2 );

  // Output number and dimension of vertices
  std::cout << " .. V1 = " << V1.rows() << " x " << V1.cols() << std::endl;
  std::cout << " .. V2 = " << V2.rows() << " x " << V2.cols() << std::endl;

  // Compute edge flaps for bending energy
  Eigen::MatrixXi uE, EF, EI;
  Eigen::VectorXi EMAP;
  igl::edge_flaps( F1, uE, EMAP, EF, EI );

  // Evaluate energies, gradients, and Hessians
  Eigen::MatrixXd grad;
  Eigen::SparseMatrix<double> Hess;

  std::cout << " --- MEMBRANE --- " << std::endl;
  std::cout << " .. membrane_energy = " << shell::membrane_energy( V1, V2, F1 ) << std::endl;

  shell::membrane_deformed_gradient( V1, V2, F1, grad );
  std::cout << " .. membrane_deformed_gradient.norm = " << grad.norm() << std::endl;

  shell::membrane_undeformed_gradient( V1, V2, F1, grad );
  std::cout << " .. membrane_undeformed_gradient.norm = " << grad.norm() << std::endl;

  shell::membrane_deformed_hessian( V1, V2, F1, Hess );
  std::cout << " .. membrane_deformed_hessian.norm = " << Hess.norm() << std::endl;

  shell::membrane_undeformed_hessian( V1, V2, F1, Hess );
  std::cout << " .. membrane_undeformed_hessian.norm = " << Hess.norm() << std::endl;

  shell::membrane_mixed_hessian( V1, V2, F1, Hess );
  std::cout << " .. membrane_mixed_hessian.norm = " << Hess.norm() << std::endl;

  std::cout << std::endl;

  std::cout << " --- BENDING --- " << std::endl;
  std::cout << " .. bending_energy = " << shell::bending_energy( V1, V2, F1, uE, EMAP, EF, EI ) << std::endl;

  shell::bending_deformed_gradient( V1, V2, F1, uE, EMAP, EF, EI, grad );
  std::cout << " .. bending_deformed_gradient.norm = " << grad.norm() << std::endl;

  shell::bending_undeformed_gradient( V1, V2, F1, uE, EMAP, EF, EI, grad );
  std::cout << " .. bending_undeformed_gradient.norm = " << grad.norm() << std::endl;

  shell::bending_deformed_hessian( V1, V2, F1, uE, EMAP, EF, EI, Hess );
  std::cout << " .. bending_deformed_hessian.norm = " << Hess.norm() << std::endl;

  shell::bending_undeformed_hessian( V1, V2, F1, uE, EMAP, EF, EI, Hess );
  std::cout << " .. bending_undeformed_hessian.norm = " << Hess.norm() << std::endl;

  shell::bending_mixed_hessian( V1, V2, F1, uE, EMAP, EF, EI, Hess );
  std::cout << " .. bending_mixed_hessian.norm = " << Hess.norm() << std::endl;

  std::cout << std::endl;

  std::cout << " --- SHELL --- " << std::endl;
  std::cout << " .. shell_energy = " << shell::shell_energy( V1, V2, F1, uE, EMAP, EF, EI, 0.001 ) << std::endl;

  shell::shell_deformed_gradient( V1, V2, F1, uE, EMAP, EF, EI, 0.001, grad );
  std::cout << " .. shell_deformed_gradient.norm = " << grad.norm() << std::endl;

  shell::shell_undeformed_gradient( V1, V2, F1, uE, EMAP, EF, EI, 0.001, grad );
  std::cout << " .. shell_undeformed_gradient.norm = " << grad.norm() << std::endl;

  shell::shell_deformed_hessian( V1, V2, F1, uE, EMAP, EF, EI, 0.001, Hess );
  std::cout << " .. shell_deformed_hessian.norm = " << Hess.norm() << std::endl;

  shell::shell_undeformed_hessian( V1, V2, F1, uE, EMAP, EF, EI, 0.001, Hess );
  std::cout << " .. shell_undeformed_hessian.norm = " << Hess.norm() << std::endl;

  shell::shell_mixed_hessian( V1, V2, F1, uE, EMAP, EF, EI, 0.001, Hess );
  std::cout << " .. shell_mixed_hessian.norm = " << Hess.norm() << std::endl;
}