! Pre-processing
! ......

module nfe_colvar_mod

implicit none

private
!=================== BASIC DEFINITION ===========================
integer, public, parameter :: COLVAR_ANGLE           = 1
......
! Define the index of CUSTOM_CV
integer, public, parameter :: COLVAR_CUSTOM_CV       = XX

......
......
......
!================================================================
! include the CUSTOM_CV to colvar_value interface
function colvar_value(cv, x) result(value)
......
   select case(cv%type)
   ......
      case(COLVAR_CUSTOM_CV)
         value = v_CUSTOM_CV(cv, x)
   ......
   end select
end function colvar_value

!================================================================
! include the CUSTOM_CV to colvar_force interface
subroutine colvar_force(cv, x, fcv, f)
......
   select case(cv%type)
   ......
      case(COLVAR_CUSTOM_CV)
         call f_CUSTOM_CV(cv, x, fcv, f)
   ......
   end select
end subroutine colvar_force

......
!================================================================
! include the CUSTOM_CV to colvar_print interface
subroutine colvar_print(cv, lun)
......
   select case(cv%type)
   ......
      case(COLVAR_CUSTOM_CV)
         write (unit = lun, fmt = '(a)') 'CUSTOM_CV'''
         call p_CUSTOM_CV(cv, lun)
   ......
   end select
end subroutine colvar_print

!================================================================
! include the CUSTOM_CV to colvar_bootstrap interface
subroutine colvar_bootstrap(cv, cvno, amass)
......
   select case(cv%type)
   ......
      case(COLVAR_CUSTOM_CV)
         call b_CUSTOM_CV(cv, cvno)
   ......
   end select
end subroutine colvar_bootstrap

......
......
......

!================================================================
! calculate the value of CUSTOM_CV
function v_CUSTOM_CV(cv, x) result(value)

   use nfe_lib_mod
   use parallel_dat_mod
!  use other modules needed

   implicit none

   double precision :: value

   type(colvar_t), intent(in) :: cv

   double precision, intent(in) :: x(*)

   integer :: a1, a2

   nfe_assert(cv%type == COLVAR_CUSTOM_CV)

   nfe_assert(associated(cv%i))

#ifdef MPI
   if (mytaskid.eq.0) then
#endif /* MPI */
      value = (formula to get the value of CV)
#ifdef MPI
   else
      value = ZERO
   end if
#endif /* MPI */

end function v_CUSTOM_CV

!=============================================================================
! calculate the force
subroutine f_CUSTOM_CV(cv, x, fcv, f)

   use nfe_lib_mod
   use parallel_dat_mod

   implicit none

   type(colvar_t), intent(in) :: cv

   double precision, intent(in) :: x(*), fcv

   double precision, intent(inout) :: f(*)

   double precision :: d1(3), d2(3)

   integer :: a1, a2
#ifdef MPI
   integer :: error
#endif /* MPI */

   nfe_assert(cv%type == COLVAR_CUSTOM_CV)

   NFE_MASTER_ONLY_BEGIN
   
   d1(1:3) =  ( get the derivatives with
   d2(1:3) =    respect to the coordinates )

   f(a1:a1 + 2) = f(a1:a1 + 2) + fcv*d1
   f(a2:a2 + 2) = f(a2:a2 + 2) + fcv*d2

   NFE_MASTER_ONLY_END
#ifdef MPI
   call mpi_bcast(f(1:3*pmemd_natoms()), 3*pmemd_natoms(), MPI_DOUBLE_PRECISION, 0, pmemd_comm, error)
   nfe_assert(error.eq.0)
#endif /* MPI */

end subroutine f_CUSTOM_CV

!=============================================================================
! bootstrap subroutine to do the pre-check
subroutine b_CUSTOM_CV(cv, cvno)

   use nfe_lib_mod
   use parallel_dat_mod

   implicit none

   type(colvar_t), intent(inout) :: cv
   integer,        intent(in)    :: cvno
   integer :: error

   nfe_assert(cv%type == COLVAR_CUSTOM_CV)
   ! check the number of elements in cv_i
   call check_i(cv%i, cvno, 'CUSTOM_CV', ?)

   ! check the number of elements in cv_r
   if ((size(cv%r).ne.?)) then
      NFE_MASTER_ONLY_BEGIN
         write (unit = ERR_UNIT, fmt = ...) &
           ......
      NFE_MASTER_ONLY_END
      call terminate()
   end if

end subroutine b_CUSTOM_CV

!=============================================================================
! print info of CUSTOM_CV
subroutine p_CUSTOM_CV(cv, lun)

   use nfe_lib_mod

   implicit none

   type(colvar_t), intent(in) :: cv
   integer, intent(in) :: lun

   nfe_assert(is_master())
   nfe_assert(cv%type == COLVAR_CUSTOM_CV)
   nfe_assert(associated(cv%i))

   write (unit = lun, fmt = ...) &
      ......
   call print_i(cv%i, lun)

end subroutine p_CUSTOM_CV

.......
.......
.......
end module nfe_colvar_mod
