tem_eval_miescatter_magnradial Function

private function tem_eval_miescatter_magnradial(me, coord, time, n) result(res)

Evaluate magnetizing field (radial-component) for Mie-Scattering of electromagnetic wave at dielectric cylinder.

Arguments

Type IntentOptional Attributes Name
type(tem_miescatter_field_type), intent(in) :: me

The function to evaluate
real(kind=rk), intent(in) :: coord(n,3)

barycentric Ids of an elements. 1st index goes over number of elements and 2nd index goes over x,y,z coordinates
real(kind=rk), intent(in) :: time

The time to evaluate the function at.
integer, intent(in) :: n

Number of points to evaluate the function for.

Return Value real(kind=rk), (n)

return value of the function


Calls

proc~~tem_eval_miescatter_magnradial~~CallsGraph proc~tem_eval_miescatter_magnradial tem_eval_miescatter_magnradial proc~cart2polar cart2polar proc~tem_eval_miescatter_magnradial->proc~cart2polar proc~hankel2_n hankel2_n proc~tem_eval_miescatter_magnradial->proc~hankel2_n

Called by

proc~~tem_eval_miescatter_magnradial~~CalledByGraph proc~tem_eval_miescatter_magnradial tem_eval_miescatter_magnradial proc~tem_eval_miescatter_magnx tem_eval_miescatter_magnx proc~tem_eval_miescatter_magnx->proc~tem_eval_miescatter_magnradial proc~tem_eval_miescatter_magny tem_eval_miescatter_magny proc~tem_eval_miescatter_magny->proc~tem_eval_miescatter_magnradial proc~tem_spacetime_for_coord tem_spacetime_for_coord proc~tem_spacetime_for_coord->proc~tem_eval_miescatter_magnx proc~tem_spacetime_for_coord->proc~tem_eval_miescatter_magny proc~tem_spatial_for_coord tem_spatial_for_coord proc~tem_spatial_for_coord->proc~tem_eval_miescatter_magnx proc~tem_spatial_for_coord->proc~tem_eval_miescatter_magny interface~tem_spacetime_for tem_spacetime_for interface~tem_spacetime_for->proc~tem_spacetime_for_coord interface~tem_spatial_for tem_spatial_for interface~tem_spatial_for->proc~tem_spatial_for_coord proc~tem_spacetime_for_stcoord tem_spacetime_for_stcoord proc~tem_spacetime_for_stcoord->proc~tem_spacetime_for_coord proc~tem_spacetime_scalar_for_index tem_spacetime_scalar_for_index proc~tem_spacetime_scalar_for_index->proc~tem_spacetime_for_coord proc~tem_spatial_scalar_for_index tem_spatial_scalar_for_index proc~tem_spatial_scalar_for_index->proc~tem_spatial_for_coord

Source Code

  function tem_eval_miescatter_magnradial(me, coord, time, n) result(res)
    ! ---------------------------------------------------------------------------
    !> The function to evaluate
    type(tem_miescatter_field_type), intent(in) :: me
    !> Number of points to evaluate the function for.
    integer, intent(in) :: n
    !> barycentric Ids of an elements.
    !! 1st index goes over number of elements and
    !! 2nd index goes over x,y,z coordinates
    real(kind=rk), intent( in ) :: coord(n,3)
    !> The time to evaluate the function at.
    real(kind=rk), intent( in )  :: time
    !> return value of the function
    real(kind=rk) :: res(n)
    ! ---------------------------------------------------------------------------
    integer :: iPoint, iCoeff
    ! Polar coordinate vector in x-y plane.
    ! First entry: radius \n
    ! Second entry: angle
    real(kind=rk) :: polar(2)
    complex(kind=rk) :: tmp
    ! ---------------------------------------------------------------------------

    do iPoint = 1, n

      ! Convert to polar coordinates (relative to the center of the
      ! cylinder.
      polar = cart2polar( coord(iPoint,1)-me%miescatter%center(1), &
        &                 coord(iPoint,2)-me%miescatter%center(2)  )

      ! Inside the cylinder
      if(polar(1).le.me%miescatter%radius) then

        tmp = (0.0_rk, 0.0_rk)
        do iCoeff = -me%mieexpansion%nCoeffs , me%mieexpansion%nCoeffs
          tmp = tmp                                                           &
            & + me%mieexpansion%c_tot(iCoeff)                                 &
            & * bessel_jn(iCoeff, me%miescatter%wavenumber_cylinder*polar(1)) &
            & * (0.0_rk,1.0_rk) * iCoeff                                      &
            & * exp( (0.0_rk,1.0_rk) * iCoeff * polar(2) )
        end do
        ! Twiddle by phase and get the real part, first we check if distance from
        ! center is above a certain threshold (avoid division by zero).
        if(abs(polar(1)) > epsilon(1.0_rk) ) then
          res(iPoint) = real( tmp                                        &
            ! JZ: the paper has an additional factor of -1 here. However, I think
            ! that this factor is wrong at this place. So, I commented out
            ! the factor.
            !     & * (-1.0_rk) &
            & * (exp( (0.0_rk,1.0_rk)*me%miescatter%omega*time))         &
            & * (0.0_rk,1.0_rk)                                          &
            & / (me%miescatter%omega*me%miescatter%permeability_cylinder &
            &    *polar(1))                                              &
            & )
        else
          res(iPoint) = 0.0_rk
        end if

      ! Outside the cylinder
      else

        ! Add up incoming and scattered field
        tmp = (0.0_rk, 0.0_rk)
        do iCoeff = -me%mieexpansion%nCoeffs , me%mieexpansion%nCoeffs
          tmp = tmp +                                                     &
            & ( (0.0_rk,1.0_rk)**(-iCoeff)                                &
            &   * bessel_jn(iCoeff,                                       &
            &               me%miescatter%wavenumber_background*polar(1)) &
            &   + me%mieexpansion%c_scat(iCoeff)                          &
            &   * hankel2_n(iCoeff,                                       &
            &               me%miescatter%wavenumber_background*polar(1)) &
            & )                                                           &
            & * (0.0_rk,1.0_rk) * iCoeff                                  &
            & * exp( (0.0_rk,1.0_rk) * iCoeff * polar(2) )
        end do
        ! Twiddle by phase and get the real part, first we check if distance from
        ! center is above a certain threshold (avoid division by zero).
        if(abs(polar(1)) > epsilon(1.0_rk) ) then
          res(iPoint) = real( tmp                                          &
            ! JZ: the paper has an additional factor of -1 here. However, I think
            ! that this factor is wrong at this place. So, I commented out
            ! the factor.
            !     & * (-1.0_rk) &
            & * (exp( (0.0_rk,1.0_rk) * me%miescatter%omega * time ))      &
            & * (0.0_rk,1.0_rk)                                            &
            & / (me%miescatter%omega*me%miescatter%permeability_background &
            &    *polar(1))                                                &
            & )
        else
          res(iPoint) = 0.0_rk
        end if

      end if

    end do

  end function tem_eval_miescatter_magnradial