//
// meam_implementation.hpp
//
// LGPL Version 2.1 HEADER START
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
//
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
//
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
// MA 02110-1301  USA
//
// LGPL Version 2.1 HEADER END
//

//
// Copyright (c) 2020--2023, Regents of the University of Minnesota.
// All rights reserved.
//
// Contributors:
//    Yaser Afshar
//

#ifndef MEAM_IMPLEMENTATION_HPP
#define MEAM_IMPLEMENTATION_HPP

#include <cstdio>
#include <memory>
#include <string>
#include <vector>

#include "KIM_ModelDriverHeaders.hpp"
#include "helper.hpp"
#include "meam_c.hpp"
#include "meam_spline.hpp"
#include "meam_sw_spline.hpp"

/*! \class MEAMImplementation
 * \brief MEAM model driver Implementation class
 */
class MEAMImplementation {
 public:
  /*!
   * \brief Construct a new MEAMImplementation object
   *
   * \param model_driver_create A %KIM API Model object
   * \param requested_length_unit User reuested length unit
   * \param requested_energy_unit User reuested energy unit
   * \param requested_charge_unit User reuested charge unit
   * \param requested_temperature_unit User reuested temperature unit
   * \param requested_time_unit User reuested time unit
   * \param ier 0|false if everything goes well and 1|true if it fails
   */
  MEAMImplementation(KIM::ModelDriverCreate *const model_driver_create,
                     KIM::LengthUnit const requested_length_unit,
                     KIM::EnergyUnit const requested_energy_unit,
                     KIM::ChargeUnit const requested_charge_unit,
                     KIM::TemperatureUnit const requested_temperature_unit,
                     KIM::TimeUnit const requested_time_unit, int *const ier);

  /*!
   * \brief Destroy the MEAMImplementation object
   *
   */
  ~MEAMImplementation() = default;

 private:
  /*!
   * \brief Open the MEAM parameter files
   *
   * \param model_driver_create A %KIM API Model object
   * \param number_parameter_files Number of parameter files to open
   * \param parameter_file_pointers FILE pointer to the opened files
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int OpenParameterFiles(KIM::ModelDriverCreate *const model_driver_create,
                         int const number_parameter_files,
                         std::FILE **parameter_file_pointers);

  /*!
   * \brief Parse the input files and process the read parameters
   *
   * \param model_driver_create A %KIM API Model object
   * \param number_parameter_files Number of parameter files to open
   * \param parameter_file_pointers FILE pointer to the opened files
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int ProcessParameterFiles(KIM::ModelDriverCreate *const model_driver_create,
                            int const number_parameter_files,
                            std::FILE *const *parameter_file_pointers);

  /*!
   * \brief Close all the opened files
   *
   * \param number_parameter_files Number of parameter files to close
   * \param parameter_file_pointers FILE pointer to the opened files
   */
  void CloseParameterFiles(int const number_parameter_files,
                           std::FILE *const *parameter_file_pointers);

  /*!
   * \brief Convert units of the parameters from the input files
   *
   * \param model_driver_create A %KIM API Model object
   * \param requested_length_unit User reuested length unit
   * \param requested_energy_unit User reuested energy unit
   * \param requested_charge_unit User reuested charge unit
   * \param requested_temperature_unit User reuested temperature unit
   * \param requested_time_unit User reuested time unit
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int ConvertUnits(KIM::ModelDriverCreate *const model_driver_create,
                   KIM::LengthUnit const &requested_length_unit,
                   KIM::EnergyUnit const &requested_energy_unit,
                   KIM::ChargeUnit const &requested_charge_unit,
                   KIM::TemperatureUnit const &requested_temperature_unit,
                   KIM::TimeUnit const &requested_time_unit);

  /*!
   * \brief Set the mutable values
   *
   * \param model_compute_arguments A %KIM API ComputeArguments object
   * \param is_compute_energy ComputeEnergy flag
   * \param is_compute_forces ComputeForces flag
   * \param is_compute_particle_energy ComputeParticleEnergy flag
   * \param is_compute_virial ComputeVirial flag
   * \param is_compute_particle_virial ComputeParticleVirial flag
   * \param particle_species_codes Particle species code
   * \param particle_contributing Particle contirubuting flag list
   * \param coordinates Particles' coordinates
   * \param energy System energy
   * \param forces Particles' forces
   * \param particle_energy Particles' energy
   * \param virial System virial
   * \param particle_virial Particles' virial
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int SetComputeMutableValues(
      KIM::ModelComputeArguments const *const model_compute_arguments,
      bool &is_compute_energy, bool &is_compute_forces,
      bool &is_compute_particle_energy, bool &is_compute_virial,
      bool &is_compute_particle_virial, int const *&particle_species_codes,
      int const *&particle_contributing, VectorOfSizeDIM const *&coordinates,
      double *&energy, VectorOfSizeDIM *&forces, double *&particle_energy,
      VectorOfSizeSix *&virial, VectorOfSizeSix *&particle_virial);

  /*!
   * \brief Get the compute index
   *
   * \param is_compute_energy ComputeEnergy flag
   * \param is_compute_forces ComputeForces flag
   * \param is_compute_particle_energy ComputeParticleEnergy flag
   * \param is_compute_virial ComputeVirial flag
   * \param is_compute_particle_virial ComputeParticleVirial flag
   *
   * \return int The comput index
   */
  int GetComputeIndex(bool const is_compute_energy,
                      bool const is_compute_forces,
                      bool const is_compute_particle_energy,
                      bool const is_compute_virial,
                      bool const is_compute_particle_virial) const;

  /*!
   * \brief Get the Compute Index object
   *
   * \param is_compute_energy ComputeEnergy flag
   * \param is_compute_forces ComputeForces flag
   * \param is_compute_particle_energy ComputeParticleEnergy flag
   * \param is_compute_virial ComputeVirial flag
   * \param is_compute_particle_virial ComputeParticleVirial flag
   * \param are_knots_on_regular_grid are_knots_on_regular_grid flag
   * \return int The comput index
   */
  int GetComputeIndex(bool const is_compute_energy,
                      bool const is_compute_forces,
                      bool const is_compute_particle_energy,
                      bool const is_compute_virial,
                      bool const is_compute_particle_virial,
                      bool const are_knots_on_regular_grid) const;

  /*!
   * \brief Use (possibly) new values of parameters to compute other quantities
   *
   * \tparam ModelObj A %KIM API Model object
   *
   * \param model_obj It is a %KIM API ModelDriverCreate object during
   * initialization or a ModelRefresh object when the Model's parameters have
   * been altered
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  template <class ModelObj>
  int SetRefreshMutableValues(ModelObj *const model_obj);

  /*!
   * \brief Set the Model's particle Numbering
   *
   * \param model_driver_create A %KIM API Model object
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int RegisterKIMModelSettings(
      KIM::ModelDriverCreate *const model_driver_create) const;

  /*!
   * \brief Register %KIM API parameters
   *
   * \param model_driver_create A %KIM API Model object
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int RegisterKIMParameters(KIM::ModelDriverCreate *const model_driver_create);

  /*!
   * \brief Register %KIM API Functions
   *
   * \param model_driver_create A %KIM API Model object
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int RegisterKIMFunctions(
      KIM::ModelDriverCreate *const model_driver_create) const;

  /*!
   * \brief Register %KIM API arguments settings
   *
   * \param model_compute_arguments_create
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int RegisterKIMComputeArgumentsSettings(
      KIM::ModelComputeArgumentsCreate *const model_compute_arguments_create)
      const;

  /*!
   * \brief Main compute routine for `meam/c` potential
   *
   * \tparam IsComputeEnergy ComputeEnergy flag
   * \tparam IsComputeForces ComputeForces flag
   * \tparam IsComputeParticleEnergy ComputeParticleEnergy flag
   * \tparam IsComputeVirial ComputeVirial flag
   * \tparam IsComputeParticleVirial ComputeParticleVirial flag
   *
   * \param model_compute A %KIM API Model object
   * \param model_compute_arguments A %KIM API ComputeArguments object
   * \param particle_species_codes Species code
   * \param particle_contributing Contribution array
   * \param coordinates Particles coordinates
   * \param energy Global energy value
   * \param forces Particles forces
   * \param particle_energy Particles energy
   * \param virial Global virial term
   * \param particle_virial Particles virial
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  template <bool IsComputeEnergy, bool IsComputeForces,
            bool IsComputeParticleEnergy, bool IsComputeVirial,
            bool IsComputeParticleVirial>
  int MeamCCompute(
      KIM::ModelCompute const *const model_compute,
      KIM::ModelComputeArguments const *const model_compute_arguments,
      int const *const particle_species_codes,
      int const *const particle_contributing,
      const VectorOfSizeDIM *const coordinates, double *const energy,
      VectorOfSizeDIM *const forces, double *const particle_energy,
      VectorOfSizeSix virial, VectorOfSizeSix *const particle_virial) const;

  /*!
   * \brief Main compute routine for `meam/spline` potential
   *
   * \tparam IsComputeEnergy ComputeEnergy flag
   * \tparam IsComputeForces ComputeForces flag
   * \tparam IsComputeParticleEnergy ComputeParticleEnergy flag
   * \tparam IsComputeVirial ComputeVirial flag
   * \tparam IsComputeParticleVirial ComputeParticleVirial flag
   *
   * \param model_compute A %KIM API Model object
   * \param model_compute_arguments A %KIM API ComputeArguments object
   * \param particle_species_codes Species code
   * \param particle_contributing Contribution array
   * \param coordinates Particles coordinates
   * \param energy Global energy value
   * \param forces Particles forces
   * \param particle_energy Particles energy
   * \param virial Global virial term
   * \param particle_virial Particles virial
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  template <bool IsComputeEnergy, bool IsComputeForces,
            bool IsComputeParticleEnergy, bool IsComputeVirial,
            bool IsComputeParticleVirial, bool AreKnotsOnRegularGid>
  int MeamSplineCompute(
      KIM::ModelCompute const *const model_compute,
      KIM::ModelComputeArguments const *const model_compute_arguments,
      int const *const particle_species_codes,
      int const *const particle_contributing,
      const VectorOfSizeDIM *const coordinates, double *const energy,
      VectorOfSizeDIM *const forces, double *const particle_energy,
      VectorOfSizeSix virial, VectorOfSizeSix *const particle_virial) const;

  /*!
   * \brief Main compute routine for `meam/sw/spline` potential
   *
   * \tparam IsComputeEnergy ComputeEnergy flag
   * \tparam IsComputeForces ComputeForces flag
   * \tparam IsComputeParticleEnergy ComputeParticleEnergy flag
   * \tparam IsComputeVirial ComputeVirial flag
   * \tparam IsComputeParticleVirial ComputeParticleVirial flag
   *
   * \param model_compute A %KIM API Model object
   * \param model_compute_arguments A %KIM API ComputeArguments object
   * \param particle_species_codes Species code
   * \param particle_contributing Contribution array
   * \param coordinates Particles coordinates
   * \param energy Global energy value
   * \param forces Particles forces
   * \param particle_energy Particles energy
   * \param virial Global virial term
   * \param particle_virial Particles virial
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  template <bool IsComputeEnergy, bool IsComputeForces,
            bool IsComputeParticleEnergy, bool IsComputeVirial,
            bool IsComputeParticleVirial, bool AreKnotsOnRegularGid>
  int MeamSWSplineCompute(
      KIM::ModelCompute const *const model_compute,
      KIM::ModelComputeArguments const *const model_compute_arguments,
      int const *const particle_species_codes,
      int const *const particle_contributing,
      const VectorOfSizeDIM *const coordinates, double *const energy,
      VectorOfSizeDIM *const forces, double *const particle_energy,
      VectorOfSizeSix virial, VectorOfSizeSix *const particle_virial) const;

  /*!
   * \brief Return the total number of (effective half list) neighbors
   *
   * \param model_compute_arguments A %KIM API ComputeArguments object
   * \param particle_contributing Contribution array
   * \return int Total number of neighbors
   */
  std::size_t TotalNumberOfNeighbors(
      KIM::ModelComputeArguments const *const model_compute_arguments,
      int const *const particle_contributing) const;

  /*!
   * \brief Return the maximum number of neighbors a single atom has
   *
   * \param model_compute_arguments A %KIM API ComputeArguments object
   * \param particle_contributing Contribution array
   * \return int Maximum number of neighbors
   */
  int MaxNumberOfNeighbors(
      KIM::ModelComputeArguments const *const model_compute_arguments,
      int const *const particle_contributing) const;

 public:
  // no explicit Destroy() needed here

  /*!
   * \brief A refresh routine
   *
   * \param model_refresh A refresh routine
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int Refresh(KIM::ModelRefresh *const model_refresh);

  /*!
   * \brief A routine to write the parameterized model files
   *
   * \param model_write_parameterized_model A %KIM API Model object
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int WriteParameterizedModel(KIM::ModelWriteParameterizedModel const
                                  *const model_write_parameterized_model) const;

  /*!
   * \brief Compute routine
   *
   * \param model_compute A %KIM API Model object
   * \param model_compute_arguments A %KIM API ComputeArguments object
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int Compute(KIM::ModelCompute const *const model_compute,
              KIM::ModelComputeArguments const *const model_compute_arguments);

  /*!
   * \brief Create a compute arguments routine
   *
   * \param model_compute_arguments_create
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int ComputeArgumentsCreate(KIM::ModelComputeArgumentsCreate
                                 *const model_compute_arguments_create) const;

  /*!
   * \brief Destroy a compute arguments routine
   *
   * \param model_compute_arguments_destroy
   *
   * \return int 0|false if everything goes well and 1|true if it fails
   */
  int ComputeArgumentsDestroy(KIM::ModelComputeArgumentsDestroy
                                  *const model_compute_arguments_destroy) const;

 private:
  /*!
   * \brief Flag indicating we do not need to request neighbors from
   *        non-contributing particles for this model driver
   */
  int const model_will_not_request_neighbors_of_non_contributing_particles_{1};

  /*!
   * \brief Number of particles
   *
   * \note
   * This is a Mutable value that can change with each
   * call to \b Refresh() and \b Compute().
   */
  int cached_number_of_particles_{0};

  /*! # of unique elements */
  int number_of_elements_{0};

  /*! meam/c style flag */
  int is_meam_c_{0};

  /*! meam/spline style flag */
  int is_meam_spline_{0};

  /*! meam/sw/spline style flag */
  int is_meam_sw_spline_{0};

  /*! max cutoff for all elements */
  double max_cutoff_{0.0};

  /*! max cutoff squared for all elements */
  double max_cutoff_squared_{0.0};

  /*!
   * \brief Cutoff value in %KIM API object
   *
   * \note
   * This is a Mutable value that can change with each
   * call to \b Refresh() and \b Compute().
   */
  double influence_distance_{0.0};

  /*! names of unique elements */
  std::vector<std::string> element_name_;

  /*! MEAMC object */
  std::unique_ptr<MEAMC> meam_c_;

  /*! MEAMSpline object */
  std::unique_ptr<MEAMSpline> meam_spline_;

  /*! MEAMSWSpline object */
  std::unique_ptr<MEAMSWSpline> meam_sw_spline_;
};

#endif  // MEAM_IMPLEMENTATION_HPP