!---------------------------------------------------------------------
! Contains the object and routines required to work with 3D vectors 
!---------------------------------------------------------------------
Module vector3D
  !-----------------------------------------------------------
  ! Specify which modules to use here. There none in this case
  !-----------------------------------------------------------
  !  Use module_xyz, Only: functions_in_xyz, variables_in_xyz

  Implicit None
  Save

  !------------------------------------------------------------------
  ! Make everything private by default, then specify which variables 
  ! and functions are allowed for Public Use by other modules 
  !------------------------------------------------------------------
  Private  
  Public :: VecType, Operator(+), Operator(*), Assignment(=), &
      vector3D_init, vector3D_display, vector3D_getnormsq, &
      vector3D_mult_scalar

  !---------------------------------------------------------------
  ! Data Type Definition of the object associated with his module
  !---------------------------------------------------------------
  Type VecType
    Real*8    :: x, y, z
  End Type VecType

  !-----------------------------------
  ! Variable Declarations
  !-----------------------------------
  ! None required in this module

  !-----------------------------------------------------------
  ! Interface Specification for the Subroutines and Functions
  !-----------------------------------------------------------
  ! Static Polymorphism
  Interface vector3D_init
    Module Procedure vector3D_init1
    Module Procedure vector3D_init2
  End Interface

  ! Operator Overloading 
  Interface Operator(+)
    Module Procedure vector3D_add
  End Interface

  Interface Operator(*)
    Module Procedure vector3D_mult_scalar
  End Interface

  ! Operator Overloading 
  Interface Assignment(=)
    Module Procedure vector3D_vec_eq_scalar
  End Interface


  !------------------------------------------------------
  ! Definition of Associated Subroutines and Functions
  !------------------------------------------------------ 
Contains
  !------------------------------------------------------------------
  ! Initialize the vec "object" to the values stored in the array arr
  !------------------------------------------------------------------
  Subroutine vector3D_init1(vec, arr)
    Type(VecType), Intent(out) :: vec
    Real*8, Dimension(:), Intent(in) :: arr

    ! Assign the value from the array
    vec%x = arr(1)
    vec%y = arr(2)
    vec%z = arr(3)
  End Subroutine vector3D_init1

  !---------------------------------------------------------------
  ! Initialize all the components of the vec "object" to the value
  ! in "num"
  !---------------------------------------------------------------
  Subroutine vector3D_init2(vec, num)
    Type(VecType), Intent(out) :: vec
    Real*8, Intent(in)         :: num

    ! Assign "num" to elements of "vec"
    vec%x = num
    vec%y = num
    vec%z = num
  End Subroutine vector3D_init2

  !------------------------------------------------------
  !Performs " vec1 + vec2 ".  This is operator overloaded
  !------------------------------------------------------
  Type(VecType) Function vector3D_add(vec1, vec2)
    Type(VecType), Intent(in) :: vec1, vec2

    ! Do the addition
    vector3D_add%x = vec1%x + vec2%x
    vector3D_add%y = vec1%y + vec2%y
    vector3D_add%z = vec1%z + vec2%z
  End Function vector3D_add


  !------------------------------------------------------
  ! Multiply vector "vec" by scalar "num" 
  !------------------------------------------------------
  Type(VecType) Function vector3D_mult_scalar(vec, num)
    Type(VecType), Intent(in) :: vec
    Real, Intent(in) :: num
    vector3D_mult_scalar%x=num*vec%x
    vector3D_mult_scalar%y=num*vec%y
    vector3D_mult_scalar%z=num*vec%z
  End Function vector3D_mult_scalar

  !--------------------------------------------------------
  ! returns the square of the magnitude of the vector "vec" 
  !--------------------------------------------------------
  Real Function vector3D_getnormsq(vec )
    Type(VecType), Intent(in) :: vec
    vector3D_getnormsq=(vec%x**2)+(vec%y**2)+vec%z**2
  End Function vector3D_getnormsq


  !-------------------------------------------------
  ! Initialize all of a vector's components to "num"
  !-------------------------------------------------
  Subroutine vector3D_vec_eq_scalar(vec, num)
    Type(VecType), Intent(out) :: vec
    Real*8, Intent(in) :: num
    vec%x=num
    vec%y=num
    vec%z=num
  End Subroutine vector3D_vec_eq_scalar
  !------------------------------------------------------------------
  ! Write contents of "vec" structure to the logical unit no "unitno"
  !------------------------------------------------------------------
  Subroutine vector3D_display(vec, optunitno)
    Type(VecType), Intent(in) :: vec
    Integer, Optional, Intent(in) :: optunitno
    Integer :: unitno
    If (Present(optunitno)) Then
      unitno = optunitno
    Else
      unitno = 6 ! Write to standard output
    End If
    Write(unitno,*) "Components of vector : ", vec%x, vec%y, vec%z
  End Subroutine vector3D_display
End Module vector3D