UltrasphericalTransform

Discrete Ultraspherical transform.

Interface 1

܀ Interface

INTERFACE
  MODULE PURE FUNCTION UltrasphericalTransform1(n, lambda, coeff, x, w, &
    &  quadType) RESULT(ans)
    INTEGER(I4B), INTENT(IN) :: n
    !! order of jacobi polynomial
    REAL(DFP), INTENT(IN) :: lambda
    !! $\lambda > -0.5, \lambda \ne 0.0$
    REAL(DFP), INTENT(IN) :: coeff(0:n)
    !! nodal value (at quad points)
    REAL(DFP), INTENT(IN) :: x(0:n)
    !! quadrature points
    REAL(DFP), INTENT(IN) :: w(0:n)
    !! weights
    INTEGER(I4B), INTENT(IN) :: quadType
    !! Quadrature type, Gauss, GaussLobatto, GaussRadau, GaussRadauLeft
    !! GaussRadauRight
    REAL(DFP) :: ans(0:n)
    !! modal values  or coefficients
  END FUNCTION UltrasphericalTransform1
END INTERFACE

️܀ See example

This example shows the usage of UltrasphericalTransform method.

This routine performs the forward Ultraspherical transform.

In this example $\lambda=-0.5$ (that is, proportional to Legendre polynomial)

program main
  use easifembase
  use easifemclasses
  implicit none
  integer( i4b ) :: n
  real(dfp), allocatable :: uhat(:), u( : ), pt( : ), wt(:)
  real( dfp ), parameter :: lambda=0.5_DFP, tol=1.0E-10
  type(string) :: astr
  integer( i4b ), parameter :: quadType = GaussLobatto
  type(CSVFile_) :: afile
  character( len=* ), parameter :: fname = "./results/test21"
  type(PLPlot_) :: aplot

  n = 10
  call reallocate( pt, n+1, wt, n+1 )
  call UltrasphericalQuadrature( n=n+1, lambda=lambda, &
    & pt=pt, wt=wt, quadType=quadType )
  u = UltrasphericalEval(n=5_I4B, lambda=lambda, &
    & x=pt)
  uhat = UltrasphericalTransform(n=n, lambda=lambda, coeff=u, &
    & x=pt, w=wt, quadType=quadType)
  uhat( 6 ) = uhat( 6 ) - 1.0_dFP
  CALL ok( SOFTEQ( NORM2(uhat), 0.0_DFP, tol), "test"  )

  n = 10
  call reallocate( pt, n+1, wt, n+1 )
  call UltrasphericalQuadrature( n=n+1, lambda=lambda, &
    & pt=pt, wt=wt, quadType=quadType )
  u = SIN(4.0_DFP * pi * pt)
  uhat = UltrasphericalTransform(n=n, lambda=lambda, coeff=u, x=pt, &
      & w=wt, quadType=quadType)

  CALL afile%initiate(filename=fname // ".csv", &
    & status="NEW", action="WRITE")
  CALL afile%open()
  CALL afile%write(val="n="//tostring(n))
  CALL afile%write(val=uhat)
  CALL afile%write()

  CALL aplot%Initiate()
  CALL aplot%Set( &
    & device="svg", &
    & filename=fname//"-%n.svg")
  CALL aplot%figure()
  CALL aplot%subplot(1,1,1)
  pt = arange(0,n,1)
  CALL aplot%setXYLim([-10.0_DFP, 50.0_DFP], [ -3.0_DFP, 3.0_DFP])
  CALL aplot%setTicks()

  CALL aplot%plot2D(x=pt,y=uhat, &
    & pointType=PS_DOT, &
    & pointSize=4.0_DFP, &
    & pointColor="b", &
    & lineWidth=0.0_DFP)
  !CALL aplot%setLabels("","","")
  !CALL aplot%Show()

  n = 20
  call reallocate( pt, n+1, wt, n+1 )
  call UltrasphericalQuadrature( n=n+1, lambda=lambda, &
    & pt=pt, wt=wt, quadType=quadType )
  u = SIN(4.0_DFP * pi * pt)
  uhat = UltrasphericalTransform(n=n, lambda=lambda,&
    & coeff=u, x=pt, w=wt, quadType=quadType)
  CALL afile%write(val="n="//tostring(n))
  CALL afile%write(val=uhat)
  CALL afile%write()

  pt = arange(0,n,1)
  CALL aplot%plot2D(x=pt,y=uhat, &
    & pointType=PS_DOT, &
    & pointSize=4.0_DFP, &
    & pointColor="r", &
    & lineWidth=0.0_DFP)
  !CALL aplot%setLabels("N","uhat","")
  !CALL aplot%Show()

  n = 30
  call reallocate( pt, n+1, wt, n+1 )
  call UltrasphericalQuadrature( n=n+1, lambda=lambda, &
    & pt=pt, wt=wt, quadType=quadType )
  u = SIN(4.0_DFP * pi * pt)
  uhat = UltrasphericalTransform(n=n, lambda=lambda, &
    & coeff=u, x=pt, w=wt, quadType=quadType)
  CALL afile%write(val="n="//tostring(n))
  CALL afile%write(val=uhat)
  CALL afile%write()

  pt = arange(0,n,1)
  CALL aplot%plot2D(x=pt,y=uhat, &
    & pointType=PS_DOT, &
    & pointSize=10.0_DFP, &
    & pointColor="k", &
    & lineWidth=0.0_DFP)
  CALL aplot%setLabels("N","uhat","")
  CALL aplot%Show()

  n = 40
  call reallocate( pt, n+1, wt, n+1 )
  call UltrasphericalQuadrature( n=n+1, lambda=lambda, &
    & pt=pt, wt=wt, quadType=quadType )
  u = SIN(4.0_DFP * pi * pt)
  uhat = UltrasphericalTransform(n=n, lambda=lambda, &
    & coeff=u, x=pt, w=wt, quadType=quadType)
  CALL afile%write(val="n="//tostring(n))
  CALL afile%write(val=uhat)
  CALL afile%write()

  CALL afile%deallocate()

end program main

↢

Interface 2

܀ Interface

INTERFACE
  MODULE PURE FUNCTION UltrasphericalTransform2(n, lambda, coeff, x, w, &
    & quadType) RESULT(ans)
    INTEGER(I4B), INTENT(IN) :: n
    !! order of polynomial
    REAL(DFP), INTENT(IN) :: lambda
    !! $\lambda > -0.5, \lambda \ne 0.0$
    REAL(DFP), INTENT(IN) :: coeff(0:, 1:)
    !! nodal value (at quad points)
    REAL(DFP), INTENT(IN) :: x(0:n)
    !! quadrature points
    REAL(DFP), INTENT(IN) :: w(0:n)
    !! weights
    INTEGER(I4B), INTENT(IN) :: quadType
    !! Quadrature type, Gauss, GaussLobatto, GaussRadau, GaussRadauLeft
    !! GaussRadauRight
    REAL(DFP) :: ans(0:n, 1:SIZE(coeff, 2))
    !! modal values  or coefficients for each column of val
  END FUNCTION UltrasphericalTransform2
END INTERFACE

️܀ See example

• This example shows the usage of UltrasphericalTransform method.
• This routine performs the forward Ultraspherical transform of column data.

In this example $\lambda=0.5$ (that is, proportional to Legendre polynomial)

program main
  use easifembase
  implicit none
  integer( i4b ) :: n
  real(dfp), allocatable :: uhat(:, :), u( :, : ), pt( : ), wt(:)
  real( dfp ), parameter :: lambda=0.5_DFP, tol=1.0E-10
  type(string) :: astr
  integer( i4b ), parameter :: quadType = GaussLobatto
  n = 10
  call reallocate( pt, n+1, wt, n+1 )
  call UltrasphericalQuadrature( n=n+1, lambda=lambda, &
    & pt=pt, wt=wt, quadType=quadType )
  u = SIN(4.0_DFP * pi * pt) .COLCONCAT. SIN(4.0_DFP * pi * pt)
  uhat = UltrasphericalTransform(n=n, lambda=lambda, coeff=u, &
    & x=pt, w=wt, quadType=quadType)
  CALL Display(MdEncode(uhat), "uhat=")
end program main

See results

$\tilde{u}_{n}$	$\tilde{u}_{n}$
-2.94751E-16	-2.94751E-16
-0.24244	-0.24244
-7.45173E-16	-7.45173E-16
-0.46574	-0.46574
1.71874E-15	1.71874E-15
-0.56276	-0.56276
4.33325E-15	4.33325E-15
1.7997	1.7997
4.13122E-15	4.13122E-15
-0.52876	-0.52876
-1.21279E-15	-1.21279E-15

↢

Interface 3

܀ Interface

INTERFACE
  MODULE FUNCTION UltrasphericalTransform3(n, lambda, f, quadType) &
    & RESULT(ans)
    INTEGER(I4B), INTENT(IN) :: n
    !! order of jacobi polynomial
    REAL(DFP), INTENT(IN) :: lambda
    !! $\lambda > -0.5, \lambda \ne 0.0$
    PROCEDURE(iface_1DFunction), POINTER, INTENT(IN) :: f
    !! 1D space function
    INTEGER(I4B), INTENT(IN) :: quadType
    !! Quadrature type, Gauss, GaussLobatto, GaussRadau, GaussRadauLeft
    !! GaussRadauRight
    REAL(DFP) :: ans(0:n)
    !! modal values  or coefficients
  END FUNCTION UltrasphericalTransform3
END INTERFACE

️܀ See example

This example shows the usage of UltrasphericalTransform method.

In this example $\alpha=\beta=0.0$ (that is, proportional to Legendre polynomial)

program main
  use easifembase
  integer( i4b ) :: n, ii
  real(dfp), allocatable :: uhat(:), u( : ), pt( : ), wt(:), fhat(:)
  real( dfp ), parameter :: lambda=0.5_DFP, tol=1.0E-10
  type(string) :: astr
  integer( i4b ), parameter :: quadType = GaussLobatto
  procedure(iface_1Dfunction), pointer :: f => NULL()

note

To make quadrature points we call UltrasphericalQuadrature routine, and specify QuadType.

  n = 10
  call reallocate( pt, n+1, wt, n+1, u, n+1 )
  call UltrasphericalQuadrature( n=n+1, lambda=lambda, &
    & pt=pt, wt=wt, quadType=quadType )

  f => func1
  !!
  do ii = 1, size(u)
    u(ii) = f(pt(ii))
  end do
  !!
  uhat = UltrasphericalTransform(n=n, lambda=lambda, &
    & coeff=u, x=pt, w=wt, quadType=quadType)
  !!
  fhat = UltrasphericalTransform(n=n, lambda=lambda, &
    & f=f, quadType=quadType)
  !!
  call OK(ALL(SOFTEQ(fhat, uhat, tol)), "n=10")
  !!

  contains
  pure function func1(x) result(ans)
    real(dfp), intent(in) :: x
    real(dfp) :: ans
    ans = SIN(4.0_DFP * pi * x)
  end function func1

end program main

↢