Init_Lagrange Subroutine

  subroutine Init_Lagrange(this,N,controlNodeType,M,targetNodeType)
    !! Initialize an instance of the Lagrange class
    !! On output, all of the attributes for the Lagrange class are allocated and values are initialized according to the number of
    !! control points, number of target points, and the types for the control and target nodes.
    !! If a GPU is available, device pointers for the Lagrange attributes are allocated and initialized.
    implicit none
    class(Lagrange),intent(out) :: this
    !! Lagrange class instance
    integer,intent(in)          :: N
    !! The number of control points for interpolant
    integer,intent(in)          :: M
    !! The number of target points for the interpolant
    integer,intent(in)          :: controlNodeType
    !! The integer code specifying the type of control points. Parameters are defined in SELF_Constants.f90. One of GAUSS(=1),
    !! GAUSS_LOBATTO(=2), or UNIFORM(=3)
    integer,intent(in)          :: targetNodeType
    !! The integer code specifying the type of target points. Parameters are defined in SELF_Constants.f90. One of GAUSS(=1),
    !! GAUSS_LOBATTO(=2), or UNIFORM(=3)
    ! -------!
    ! Local
    real(prec) :: q(0:M)

    if(.not. GPUAvailable()) then
      print*,__FILE__,':',__LINE__,' : Error : Attempt to use GPU extension, but GPU is not available.'
      stop 1
    endif

    this%N = N
    this%M = M
    this%controlNodeType = controlNodeType
    this%targetNodeType = targetNodeType
    allocate(this%controlPoints(1:N+1), &
             this%targetPoints(1:M+1), &
             this%bWeights(1:N+1), &
             this%qWeights(1:N+1), &
             this%iMatrix(1:N+1,1:M+1), &
             this%dMatrix(1:N+1,1:N+1), &
             this%dgMatrix(1:N+1,1:N+1), &
             this%bMatrix(1:N+1,1:2))

    if(controlNodeType == GAUSS .or. controlNodeType == GAUSS_LOBATTO) then

      call LegendreQuadrature(N, &
                              this%controlPoints, &
                              this%qWeights, &
                              controlNodeType)

    elseif(controlNodeType == CHEBYSHEV_GAUSS .or. controlNodeType == CHEBYSHEV_GAUSS_LOBATTO) then

      call ChebyshevQuadrature(N, &
                               this%controlPoints, &
                               this%qWeights, &
                               controlNodeType)

    elseif(controlNodeType == UNIFORM) then

      this%controlPoints = UniformPoints(-1.0_prec,1.0_prec,0,N)
      this%qWeights = 2.0_prec/real(N,prec)

    endif

    ! Target Points
    if(targetNodeType == GAUSS .or. targetNodeType == GAUSS_LOBATTO) then

      call LegendreQuadrature(M, &
                              this%targetPoints, &
                              q, &
                              targetNodeType)

    elseif(targetNodeType == UNIFORM) then

      this%targetPoints = UniformPoints(-1.0_prec,1.0_prec,0,M)

    endif

    call this%CalculateBarycentricWeights()
    call this%CalculateInterpolationMatrix()
    call this%CalculateDerivativeMatrix()
    this%bMatrix(1:N+1,1) = this%CalculateLagrangePolynomials(-1.0_prec)
    this%bMatrix(1:N+1,2) = this%CalculateLagrangePolynomials(1.0_prec)

    print*,"Lagrange malloc"
    call gpuCheck(hipMalloc(this%iMatrix_gpu,sizeof(this%iMatrix)))
    call gpuCheck(hipMalloc(this%dMatrix_gpu,sizeof(this%dMatrix)))
    call gpuCheck(hipMalloc(this%dgMatrix_gpu,sizeof(this%dgMatrix)))
    call gpuCheck(hipMalloc(this%bMatrix_gpu,sizeof(this%bMatrix)))
    call gpuCheck(hipMalloc(this%qWeights_gpu,sizeof(this%qWeights)))

    print*,"Lagrange memcpy"

    print*,c_loc(this%iMatrix)
    call gpuCheck(hipMemcpy(this%iMatrix_gpu, &
                            c_loc(this%iMatrix), &
                            sizeof(this%iMatrix), &
                            hipMemcpyHostToDevice))
    print*,"Lagrange memcpy"

    call gpuCheck(hipMemcpy(this%dMatrix_gpu, &
                            c_loc(this%dMatrix), &
                            sizeof(this%dMatrix), &
                            hipMemcpyHostToDevice))
    print*,"Lagrange memcpy"

    call gpuCheck(hipMemcpy(this%dgMatrix_gpu, &
                            c_loc(this%dgMatrix), &
                            sizeof(this%dgMatrix), &
                            hipMemcpyHostToDevice))
    print*,"Lagrange memcpy"

    call gpuCheck(hipMemcpy(this%bMatrix_gpu, &
                            c_loc(this%bMatrix), &
                            sizeof(this%bMatrix), &
                            hipMemcpyHostToDevice))
    print*,"Lagrange memcpy"

    call gpuCheck(hipMemcpy(this%qWeights_gpu, &
                            c_loc(this%qWeights), &
                            sizeof(this%qWeights), &
                            hipMemcpyHostToDevice))

  endsubroutine Init_Lagrange