UniformBlockMesh_Mesh1D Subroutine

public subroutine UniformBlockMesh_Mesh1D(this, nElem, x)

Arguments

TypeIntentOptionalAttributesName
class(Mesh1D), intent(out) :: this
integer, intent(in) :: nElem
real(kind=prec), intent(in) :: x(1:2)

Contents


Source Code

  subroutine UniformBlockMesh_Mesh1D(this,nElem,x)
    implicit none
    class(Mesh1D),intent(out) :: this
    integer,intent(in) :: nElem
    real(prec),intent(in) :: x(1:2)
    ! Local
    integer :: iel,ngeo
    integer :: nid,nNodes
    integer :: i
    real(prec) :: xU(1:nElem+1)
    type(Lagrange),target :: linearInterp
    type(Lagrange),target :: nGeoInterp
    type(Scalar1D) :: xLinear
    type(Scalar1D) :: xGeo

    ngeo = 1

    nNodes = nElem*(nGeo+1)
    call this%Init(nElem,nNodes,2)
    this%quadrature = GAUSS_LOBATTO

    ! Set the hopr_nodeCoords
    xU = UniformPoints(x(1),x(2),1,nElem+1)

    call linearInterp%Init(1,GAUSS_LOBATTO, &
                           nGeo,GAUSS_LOBATTO)

    call nGeoInterp%Init(nGeo,GAUSS_LOBATTO, &
                         nGeo,GAUSS_LOBATTO)

    ! Create a linear interpolant to interpolate to nGeo grid
    call xLinear%Init(linearInterp,1,nElem)
    call xGeo%Init(nGeoInterp,1,nElem)

    do iel = 1,nElem
      xLinear%interior(1:2,iel,1) = xU(iel:iel+1)
    enddo

    call xLinear%GridInterp(xGeo%interior)

    ! Set the element information
    nid = 1
    do iel = 1,nElem
      this%elemInfo(1,iel) = selfLineLinear ! Element Type
      this%elemInfo(2,iel) = 1 ! Element Zone
      this%elemInfo(3,iel) = nid ! Node Index Start
      do i = 1,nGeo+1
        this%nodeCoords(nid) = xGeo%interior(i,iel,1)
        nid = nid+1
      enddo
      this%elemInfo(4,iel) = nid-1 ! Node Index End
    enddo

    call xLinear%Free()
    call xGeo%Free()
    call linearInterp%Free()
    call nGeoInterp%Free()

  endsubroutine UniformBlockMesh_Mesh1D