Calculate the derivative of the flux label rho with respect to the internal flux label rp (the normalised poloidal flux)
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real, | intent(in) | :: | rp | |||
real, | intent(in) | :: | dr |
real function drho_drp(rp, dr)
implicit none
real, intent(in) :: rp, dr
real :: rp1, rp2, rho1, rho2, phirpmax
rp1 = rp * (1 - dr) ; rp2 = rp * (1 + dr)
if (irho == 1) then
phirpmax = phi(rpmax)
rho1 = sqrt(abs(phi(rp1) / phirpmax)) ; rho2 = sqrt(abs(phi(rp2) / phirpmax))
else if (irho == 2) then
rho1 = 0.5 * diameter(rp1) ; rho2 = 0.5 * diameter(rp2)
else if (irho == 3) then
rho1 = psifun(rp1) ; rho2 = psifun(rp2)
else if (irho == 4) then
rho1 = rhofun(rp1) ; rho2 = rhofun(rp2)
end if
drho_drp = (rho2 - rho1) / (rp2 - rp1)
end function drho_drp