exit_element Subroutine

public subroutine exit_element(TreeID, line, tree)

This subroutine checks at which face, edge or corner the line leaves the element and calculates the next element.

found intersection, so exit loop

Arguments

Type IntentOptional Attributes Name
integer(kind=long_k), intent(in) :: TreeID
type(tem_line) :: line
type(treelmesh_type), intent(in) :: tree

Calls

proc~~exit_element~~CallsGraph proc~exit_element exit_element proc~tem_baryofid tem_BaryOfId proc~exit_element->proc~tem_baryofid proc~tem_coordofid tem_CoordOfId proc~exit_element->proc~tem_coordofid proc~tem_intersec_line_line tem_intersec_line_line proc~exit_element->proc~tem_intersec_line_line proc~tem_intersec_line_plane tem_intersec_line_plane proc~exit_element->proc~tem_intersec_line_plane proc~tem_intersec_ray_point tem_intersec_ray_point proc~exit_element->proc~tem_intersec_ray_point proc~tem_baryofid->proc~tem_coordofid proc~tem_elemsizelevel tem_ElemSizeLevel proc~tem_baryofid->proc~tem_elemsizelevel proc~tem_levelof tem_LevelOf proc~tem_coordofid->proc~tem_levelof

Source Code

  subroutine exit_element( TreeID, line, tree )
    ! ---------------------------------------------------------------------------
    integer(kind=long_k),intent(in) :: TreeID
    type(tem_line) :: line
    type(treelmesh_type), intent(in) :: tree
    ! ---------------------------------------------------------------------------
    type(tem_plane) :: face
    type(tem_line) :: edge
    real(kind=rk) :: corner(3)
    type(tem_intersec) :: intersection
    logical :: intersects
    type(tem_intersec_elem) :: elem
    integer       :: coord(4)        ! spatial index triple for a given ID
    integer :: iDir
    ! ---------------------------------------------------------------------------


    ! @todo: calculate the size of the element
    ! Kann ich von hier auf diese funktion zugreifen? wie ist das bei einem
    ! gridrefinement?
    ! auch in zeile 108 in tem_geometry.f90
    coord = tem_CoordOfId(TreeID)
    elem%length = tree%global%BoundingCubeLength / real(2**coord(4), kind=rk)
    ! calculate the center of the element
    elem%center = tem_BaryOfId(tree, TreeID)

    ! loop over the faces of the element.
    ! 1 - positive x-direction, 2 - negative x-direction
    ! 3 - positive y-direction, 4 - negative y-direction
    ! 5 - positive z-direction, 6 - negative z-direction
    dirLoop: do iDir = 1, 26

      ! FACE
      if( iDir .eq. q__W )then!.or. q__E .or. q__N .or. q__S .or. q__T .or. q__B ) then
        ! check whether or not it is possible that this face is the exit face
        if( line%direction(1) .lt. 0.0_rk ) then
          ! we need the normal of the face
          face%normal = qOffset( iDir, : )

          !< calculate the center of the face
          face%coord = elem%center + 0.5_rk * elem%length * face%normal

          call tem_intersec_line_plane( face, line, intersects, intersection )
          !> found intersection, so exit loop
          if( intersects ) exit dirLoop
        end if

      ! EDGE
      elseif( iDir .eq. q_BS ) then
        ! check whether or not it is possible that this edge is the exit edge
        if ( line%direction(2) .lt. 0.0_rk .and.                               &
          &  line%direction(3) .lt. 0.0_rk ) then
          ! we need the direction and one point of the edge
          edge%direction = qOffset( iDir, : )
          edge%coordStart = elem%center + 0.5_rk * elem%length                 &
            &             * qOffset( iDir, : )

          call tem_intersec_line_line( edge, line, intersects, intersection )

          !< found intersection, so exit loop
          if( intersects ) exit dirLoop
        end if

      elseif( iDir .eq. qBSW ) then
        ! check whether or not it is possible that this corner is the
        ! exit corner
        if( line%direction(1) .lt. 0.0_rk .and.                                &
          & line%direction(2) .lt. 0.0_rk .and.                                &
          & line%direction(3) .lt. 0.0_rk ) then
          !< calculate the coordinates of the corner
          corner = elem%center + 0.5_rk * elem%length * qOffset( iDir, : )
          !< check whether or not the corner intersects with the ray
          call tem_intersec_ray_point( corner, line, intersects, intersection )

          if( intersects ) exit dirLoop

        end if
      end if
    enddo dirLoop

  end subroutine exit_element