ScalarDataReconstruction Subroutine

public subroutine ScalarDataReconstruction(IBM, elements, STLNum, integralType, iter, autosave, dt)

Uses

    • OrientedBoundingBox
    • MappedGeometryClass
    • mpi
    • MPI_Process_Info
    • IBMClass
    • KDClass
    • MPI_IBMUtilities
    • TessellationTypes

Arguments

Type IntentOptional Attributes Name
type(IBM_type), intent(inout) :: IBM
type(element), intent(inout) :: elements(:)
integer, intent(in) :: STLNum
integer, intent(in) :: integralType
integer, intent(in) :: iter
logical :: autosave
real(kind=RP), intent(in) :: dt

Source Code

   subroutine ScalarDataReconstruction( IBM, elements, STLNum, integralType, iter, autosave, dt ) 
      use TessellationTypes
      use MappedGeometryClass
      use IBMClass
      use OrientedBoundingBox
      use KDClass
      use MPI_Process_Info
      use MPI_IBMUtilities
#ifdef _HAS_MPI_
      use mpi
#endif
!
!        -----------------------------------------------------------------------------------------
!           This function computes Scalar integrals, that is, those
!           in the form:
!                 val = \int \vec{v}ยท\vec{n}dS
!           The data at the boundary point (BP) is computed through a Inverse Distance Weight
!           procedure. 
!        -----------------------------------------------------------------------------------------
      implicit none
      !-arguments--------------------------------------------------------
      type(IBM_type), intent(inout) :: IBM
      type(element),  intent(inout) :: elements(:)
      integer,        intent(in)    :: integralType, STLNum, iter
      real(kind=RP),  intent(in)    :: dt
      !-local-variables-------------------------------------------------
      real(kind=rp), allocatable  :: Qsurf(:,:), U_xsurf(:,:), U_ysurf(:,:), U_zsurf(:,:)
      integer                     :: i, j
      logical                     :: found, autosave
       
      if( .not. IBM% Integral(STLNum)% compute ) return
      
      allocate( Qsurf(NCONS,IBM% NumOfInterPoints) )
      call IBM% BandPoint_state( elements, STLNum, .true. )
 
      if( IBM% stlSurfaceIntegrals(STLNum)% move ) then
         if( IBM% stlSurfaceIntegrals(STLNum)% motionType .eq. ROTATION ) then
            call IBM% stlSurfaceIntegrals(STLNum)% getRotationaMatrix( dt )
            call OBB(STLNum)% STL_rotate( IBM% stlSurfaceIntegrals(STLNum), .true. )
         elseif( IBM% stlSurfaceIntegrals(STLNum)% motionType .eq. LINEAR ) then
            call IBM% stlSurfaceIntegrals(STLNum)% getDisplacement( dt )
            call OBB(STLNum)% STL_translate( IBM% stlSurfaceIntegrals(STLNum), .true. )
         end if
      end if

      if( .not. MPI_Process% isRoot ) return 
!$omp parallel
!$omp do schedule(runtime) private(j,found)
      do i = 1, IBM% stlSurfaceIntegrals(STLNum)% NumOfObjs

         do j = 1, NumOfVertices + 4 
            call GetSurfaceState( IBM, IBM% stlSurfaceIntegrals(STLNum)% ObjectsList(i), IBM% stlSurfaceIntegrals(STLNum)% ObjectsList(i)% vertices(j), STLNum ) 
            Qsurf = IBM% BandRegion(STLNum)% Q(:,IBM% stlSurfaceIntegrals(STLNum)% ObjectsList(i)% vertices(j)% nearestPoints)   
         end do

         do j = 1, NumOfVertices + 4
            IBM% stlSurfaceIntegrals(STLNum)% ObjectsList(i)% vertices(j)% ScalarValue = IntegratedScalarValue( Q                 = Qsurf,                                                         &
                                                                                                                vertex            = IBM% stlSurfaceIntegrals(STLNum)% ObjectsList(i)% vertices(j), &
                                                                                                                normal            = IBM% stlSurfaceIntegrals(STLNum)% ObjectsList(i)% normal,      &
                                                                                                                integralType      = integralType,                                                  &
                                                                                                                InterpolationType = IBM% InterpolationType                                         )   
         end do
      end do
!$omp end do
!$omp end parallel
      if( IBM% stl(STLNum)% move ) then
         IBM% Integral(STLNum)% ListComputed = .false.
      else 
         IBM% Integral(STLNum)% ListComputed = .true.
      end if 
      
      if( autosave ) call GenerateScalarmonitorTECfile( IBM, STLNum, integralType, iter )

   end subroutine ScalarDataReconstruction