Skip to main content

InterpolationPoint_Triangle

This routine returns the interpolation points on triangle.

Interface

INTERFACE
MODULE FUNCTION InterpolationPoint_Triangle(order, ipType, &
& layout, xij) RESULT(ans)
INTEGER(I4B), INTENT(IN) :: order
!! order
INTEGER(I4B), INTENT(IN) :: ipType
!! interpolation point type
REAL(DFP), OPTIONAL, INTENT(IN) :: xij(:, :)
!! Coord of domain in xij format
CHARACTER(*), INTENT(IN) :: layout
!! local node numbering layout, always VEFC
REAL(DFP), ALLOCATABLE :: ans(:, :)
!! xij coordinates
END FUNCTION InterpolationPoint_Triangle
END INTERFACE
xij
  • xij contains nodal coordinates of triangle in xij format.
  • SIZE(xij,1) = nsd, and SIZE(xij,2)=3
  • If xij is absent then unit triangle is assumed
ipType
  • ipType is interpolation point type, it can take following values

  • Equidistance, uniformly/evenly distributed points

  • BlythPozChebyshev

  • BlythPozLegendre

  • GaussLegendreLobatto, which is same as IsaacLegendre

  • GaussChebyshevLobatto, which is same as IsaacChebyshev

  • GaussJacobiLobatto

  • GaussUltrasphericalLobatto

  • IsaacChebyshev

  • IsaacLegendre

  • ChenBabuska TODO

  • Hesthaven TODO

  • Feket TODO

layout
  • layout specifies the arrangement of points. The nodes are always returned in “VEFC” format (vertex, edge, face, cell). 1:3 are vertex points, then edge, and then internal nodes. The internal nodes also follow the same convention. Please read Gmsh manual on this topic.

Example of Equidistance points

Click here to see example
program main
use easifembase
implicit none
integer( i4b ) :: i1, i2, order, ipType
real( dfp ), allocatable :: x(:,:), xij(:, :)
CHARACTER( 20 ) :: layout
order=5
xij = RefTriangleCoord("UNIT")
ipType=Equidistance
layout="VEFC"
x =InterpolationPoint_Triangle( order=order, &
& ipType=ipType, layout=layout, xij=xij)
call display( Mdencode(TRANSPOSE(x)), "xij (order="//tostring(order)//")=" )
call blanklines(nol=2)
end program main
See results

xij (order=5)=

00
10
01
0.20
0.40
0.60
0.80
0.80.2
0.60.4
0.40.6
0.20.8
00.8
00.6
00.4
00.2
0.20.2
0.60.2
0.20.6
0.40.2
0.40.4
0.20.4

BlythPozChebyshev

Click here to see example
program main
use easifembase
implicit none
integer( i4b ) :: i1, i2, order, ipType
real( dfp ), allocatable :: x(:,:), xij(:, :)
CHARACTER( 20 ) :: layout
order=5
xij = RefTriangleCoord("UNIT")
ipType=BlythPozChebyshev
layout="VEFC"
x =InterpolationPoint_Triangle( order=order, &
& ipType=ipType, layout=layout, xij=xij)
call display( Mdencode(TRANSPOSE(x)), "xij (order="//tostring(order)//")=" )
call blanklines(nol=2)
end program main
See results

xij (order=5)=

00
10
01
9.54915E-02-2.02431E-14
0.34549-3.27886E-14
0.65451-3.28071E-14
0.90451-2.02431E-14
0.904519.54915E-02
0.654510.34549
0.345490.65451
9.54915E-020.90451
-2.02431E-140.90451
-3.28071E-140.65451
-3.27886E-140.34549
-2.02431E-149.54915E-02
0.146990.14699
0.706010.14699
0.146990.70601
0.416670.16667
0.416670.41667
0.166670.41667

BlythPozLegendre

Click here to see example
program main
use easifembase
implicit none
integer( i4b ) :: i1, i2, order, ipType
real( dfp ), allocatable :: x(:,:), xij(:, :)
CHARACTER( 20 ) :: layout
order=5
xij = RefTriangleCoord("UNIT")
ipType=BlythPozLegendre
layout="VEFC"
x =InterpolationPoint_Triangle( order=order, &
& ipType=ipType, layout=layout, xij=xij)
call display( Mdencode(TRANSPOSE(x)), "xij (order="//tostring(order)//")=" )
call blanklines(nol=2)
end program main
See results

xij (order=5)=

-1.85037E-16-1.85037E-16
1-1.85037E-16
-1.85037E-161
0.11747-1.4803E-16
0.35738-3.70074E-17
0.64262-3.70074E-17
0.88253-1.29526E-16
0.882530.11747
0.642620.35738
0.357380.64262
0.117470.88253
-1.29526E-160.88253
-3.70074E-170.64262
-3.70074E-170.35738
-1.4803E-160.11747
0.158290.15829
0.683430.15829
0.158290.68343
0.41330.17339
0.41330.4133
0.173390.4133

IsaacChebyshev

Click here to see example
program main
use easifembase
implicit none
integer( i4b ) :: i1, i2, order, ipType
real( dfp ), allocatable :: x(:,:), xij(:, :)
CHARACTER( 20 ) :: layout
order=5
xij = RefTriangleCoord("UNIT")
ipType=IsaacLegendre
layout="VEFC"
x =InterpolationPoint_Triangle( order=order, &
& ipType=ipType, layout=layout, xij=xij)
call display( Mdencode(TRANSPOSE(x)), "xij (order="//tostring(order)//")=" )
call blanklines(nol=2)
end program main
See results

xij (order=5)=

-8.32667E-17-8.32667E-17
1-8.32667E-17
-8.32667E-171
0.11747-4.89901E-17
0.35738-1.58335E-17
0.64262-1.58335E-17
0.88253-4.89901E-17
0.882530.11747
0.642620.35738
0.357380.64262
0.117470.88253
-4.89901E-170.88253
-1.58335E-170.64262
-1.58335E-170.35738
-4.89901E-170.11747
0.155990.15599
0.688020.15599
0.155990.68802
0.418070.16387
0.418070.41807
0.163870.41807

IssacLegendre

Click here to see example
program main
use easifembase
implicit none
integer( i4b ) :: i1, i2, order, ipType
real( dfp ), allocatable :: x(:,:), xij(:, :)
CHARACTER( 20 ) :: layout
order=5
xij = RefTriangleCoord("UNIT")
ipType=IsaacLegendre
layout="VEFC"
x =InterpolationPoint_Triangle( order=order, &
& ipType=ipType, layout=layout, xij=xij)
call display( Mdencode(TRANSPOSE(x)), "xij (order="//tostring(order)//")=" )
call blanklines(nol=2)
end program main
See results

xij (order=5)=

-8.32667E-17-8.32667E-17
1-8.32667E-17
-8.32667E-171
0.11747-4.89901E-17
0.35738-1.58335E-17
0.64262-1.58335E-17
0.88253-4.89901E-17
0.882530.11747
0.642620.35738
0.357380.64262
0.117470.88253
-4.89901E-170.88253
-1.58335E-170.64262
-1.58335E-170.35738
-4.89901E-170.11747
0.155990.15599
0.688020.15599
0.155990.68802
0.418070.16387
0.418070.41807
0.163870.41807