calculate_metric_terms Subroutine

private pure subroutine calculate_metric_terms(gradstot, gradrptot, gradthtot, dpsidrho, dqdrp, qval, rhoc, bmag, force_sym, gds2, gds21, gds22, gds23, gds24, gds24_noq)

Arguments

Type IntentOptional Attributes Name
real, intent(in), dimension(-ntgrid:, :) :: gradstot
real, intent(in), dimension(-ntgrid:, :) :: gradrptot
real, intent(in), dimension(-ntgrid:, :) :: gradthtot
real, intent(in) :: dpsidrho
real, intent(in) :: dqdrp
real, intent(in) :: qval
real, intent(in) :: rhoc
real, intent(in), dimension(-ntgrid:) :: bmag
logical, intent(in) :: force_sym
real, intent(out), dimension(-ntgrid:) :: gds2
real, intent(out), dimension(-ntgrid:) :: gds21
real, intent(out), dimension(-ntgrid:) :: gds22
real, intent(out), dimension(-ntgrid:) :: gds23
real, intent(out), dimension(-ntgrid:) :: gds24
real, intent(out), dimension(-ntgrid:) :: gds24_noq

Contents


Source Code

  pure subroutine calculate_metric_terms(gradstot, gradrptot, gradthtot, dpsidrho, dqdrp, &
       qval, rhoc, bmag, force_sym, gds2, gds21, gds22, gds23, gds24, gds24_noq)
    real, dimension(-ntgrid:, :), intent(in) :: gradstot, gradrptot, gradthtot
    real, intent(in) :: dpsidrho, dqdrp, qval, rhoc
    real, dimension(-ntgrid:), intent(in) :: bmag
    logical, intent(in) :: force_sym
    real, dimension(-ntgrid:), intent(out) :: gds2, gds21, gds22, gds23, gds24, gds24_noq
    real, dimension(:), allocatable :: gdsdum1, gdsdum2, gdsdum3, gdsdum4, gdsdum5
    allocate(gdsdum1(-ntgrid:ntgrid), gdsdum2(-ntgrid:ntgrid), gdsdum3(-ntgrid:ntgrid))
    call dottgridf(gradstot, gradstot, gdsdum1)
    call dottgridf(gradstot, gradrptot, gdsdum2)
    call dottgridf(gradrptot, gradrptot, gdsdum3)
    gds2  = gdsdum1            *dpsidrho**2
    gds21 = gdsdum2*dqdrp      *dpsidrho**2
    gds22 = gdsdum3*dqdrp*dqdrp*dpsidrho**2

    ! v_E . grad theta = i*ky*phi*(rho_{r}*v_{t,r}/a**2)*(gds23 + gsd24*theta0)
    ! calculation of gds23 and gds24 assumes rp=rho, which is not necessarily true
    ! for numerical equilibria
    allocate(gdsdum4(-ntgrid:ntgrid), gdsdum5(-ntgrid:ntgrid))
    call dottgridf(gradrptot, gradthtot, gdsdum4) !grad rho dot grad theta
    call dottgridf(gradstot, gradthtot, gdsdum5) !grad alpha dot grad theta
    gds23 = (gdsdum4 * gds2 - gdsdum5 * gdsdum2 * dpsidrho**2) / bmag**2
    gds24_noq = (gdsdum4 * gdsdum2 - gdsdum5 * gdsdum3) * (dpsidrho / bmag)**2
    gds24 = gds24_noq * dqdrp
    gds24_noq = gds24_noq * qval / rhoc

    if (force_sym) then
       call sym(gds2, 0, ntgrid)
       call sym(gds21, 1, ntgrid)
       call sym(gds22, 0, ntgrid)
       ! not sure about the following yet -- MAB
       ! call sym(gds23, 0, ntgrid); call sym(gds24, 1, ntgrid)
    end if
  end subroutine calculate_metric_terms