Cube root code example 1

cbrt_1.f90

module cbrt_mod1

implicit none
private

public :: cbrt, sp, r

interface cbrt
  module procedure cbrt_sp_sp
  module procedure cbrt_csp_csp
end interface

integer, parameter :: sp = kind(1.)
complex(sp), parameter :: r = cmplx(-1,sqrt(3._sp))/2._sp

contains

  function cbrt_sp_sp(x) result(cbrt)
    real(sp), intent(in) :: x
    real(sp) :: cbrt

    if (x >= 0.0_sp) then
      cbrt = x**(1._sp/3)
    else 
      cbrt = -(-x)**(1._sp/3)
    end if
  end function

  function cbrt_csp_csp(z) result(cbrt)
    complex(sp), intent(in) :: z
    complex(sp) :: cbrt
    cbrt = z**(1._sp/3)
  end function

end module

program test_cbrt

  use cbrt_mod1, only: cbrt,sp, r
  implicit none

  complex(sp) :: z
  character(40) :: fmtr, fmtc

  fmtr = '(A,F7.4)'
  fmtc = '(A,F7.4,SP,F7.4,"j")'

  write(*,'(A)')  "real to real"
  
  write(*,fmtr) "cbrt( 8._sp) = ", cbrt(8._sp)
  write(*,fmtr) "cbrt(-8._sp) = ", cbrt(-8._sp)

  write(*,'(/,A)')  "complex to complex"

  z = cmplx(-8._sp)
  write(*,fmtc) "z = ", z
  write(*,fmtc) "r = ", r
  write(*,fmtc) "cbrt(z)    = ", cbrt(z)
  write(*,fmtc) "cbrt(z)*r  = ",cbrt(z)*r
  write(*,fmtc) "cbrt(z)*r' = ",cbrt(z)*conjg(r)

end program