get_numerical_hessian_api Subroutine

subroutine get_numerical_hessian_api(verror, vmol, vdisp, & 
                                   & vparam, c_hessian) &
      & bind(C, name=namespace//"get_numerical_hessian")
   !DEC$ ATTRIBUTES DLLEXPORT :: get_numerical_hessian_api
   type(c_ptr), value :: verror
   type(vp_error), pointer :: error
   type(c_ptr), value :: vmol
   type(vp_structure), pointer :: mol
   type(c_ptr), value :: vdisp
   type(vp_model), pointer :: disp
   type(c_ptr), value :: vparam
   type(vp_param), pointer :: param
   real(c_double), intent(out) :: c_hessian(*)
   real(wp), allocatable :: hessian(:, :, :, :)
   integer :: nat_sq


   if (debug) print'("[Info]",1x, a)', "get_numerical_hessian"

   if (.not.c_associated(verror)) return
   call c_f_pointer(verror, error)

   if (.not.c_associated(vmol)) then
      call fatal_error(error%ptr, "Molecular structure data is missing")
      return
   end if
   call c_f_pointer(vmol, mol)
   nat_sq = mol%ptr%nat*mol%ptr%nat

   if (.not.c_associated(vdisp)) then
      call fatal_error(error%ptr, "Dispersion model is missing")
      return
   end if
   call c_f_pointer(vdisp, disp)

   if (.not.c_associated(vparam)) then
      call fatal_error(error%ptr, "Damping parameters are missing")
      return
   end if
   call c_f_pointer(vparam, param)

   if (.not.allocated(param%ptr)) then
      call fatal_error(error%ptr, "Damping parameters are not initialized")
      return
   end if

   ! Evaluate hessian numerically 
   hessian = reshape(c_hessian(:9*nat_sq), &
                    &(/3, mol%ptr%nat, 3, mol%ptr%nat/))
   call get_dispersion_hessian(mol%ptr, disp%ptr, param%ptr, &
    & realspace_cutoff(), hessian)
   c_hessian(:9*nat_sq) = reshape(hessian, (/9*nat_sq/))

end subroutine get_numerical_hessian_api