From Wikiversity - Reading time: 4 minIn the python programming language a complex square root is available for floating point numbers with a precision of 15.
# python code:
>>> (5.27 + 3.36*1j) ** 0.5
(2.4+0.7000000000000001j)
>>>
>>> (-18j) ** .5
(3-2.9999999999999996j)
>>>
Function cmath.sqrt() provides clean output:
# python code:
>>> import cmath
>>> cmath.sqrt (5.27 + 3.36*1j)
(2.4+0.7j)
>>> cmath.sqrt (-18j)
(3-3j)
>>>
If it is desired to calculate complex square root with precision greater than that available for python's floating point math,
the following code using python's decimal module will do the job. The following code also shows how complex square root
is calculated.
# python code:
import decimal
dD = decimal.Decimal
dgt = decimal.getcontext()
Precision = dgt.prec = 100 # Adjust as necessary.
def ComplexSquareRoot (v1, v2 = None) :
'''
p,q = ComplexSquareRoot (a,b) or
p,q = ComplexSquareRoot ((a,b))
a,b are the real and imaginary parts of complex number W = (a+bi)
p,q are the real and imaginary parts of complex number w = (p+qi)
where W = w ** 2
This function preserves +/- 0 as calculated by python function cmath.sqrt().
'''
thisName = 'ComplexSquareRoot (v1, v2 = None) :'
if v2 == None : a,b = v1
else : a,b = v1,v2
dgt.prec += 3
a,b = [ dD(str(v)) for v in (a,b) ]
if b == 0 :
if a == 0 : p,q = dD(0),b
else :
root = abs(a).sqrt()
if a > 0: p,q = root, b
else : p,q = dD(0), root.copy_sign(b)
elif a == 0 :
root = abs(b/2).sqrt()
if b > 0 : p,q = root, root
else : p,q = root, -root
else :
Wmod = (a**2 + b**2).sqrt()
wmod = Wmod.sqrt()
cosWφ = a/Wmod
# 2
# cos2A = 2cos A - 1
coswφ = ((cosWφ + 1) / 2).sqrt()
p = wmod * coswφ
q = b /(2*p)
dgt.prec -= 3
return [ s.normalize() for s in (p,q) ]
# python code:
import cmath
sqrt = cmath.sqrt
for (a,b) in (
(0.,0.), (0.,-0.), (-0.,0.), (-0.,-0.),
(4,0.), (4,-0.), (-4,0.), (-4,-0.),
(0.,50), (0.,-50), (-0.,50), (-0.,-50),
(-394200411798404114010884279663511687236816, 157994206778295991363266285626991662856270),
):
result1 = ComplexSquareRoot (a,b)
result2 = sqrt (complex(a,b))
str1 = 'result1, result2' ; print (str1, eval(str1))
result1, result2 ([Decimal('0'), Decimal('0')], 0j)
result1, result2 ([Decimal('0'), Decimal('-0')], -0j)
result1, result2 ([Decimal('0'), Decimal('0')], 0j)
result1, result2 ([Decimal('0'), Decimal('-0')], -0j)
result1, result2 ([Decimal('2'), Decimal('0')], (2+0j))
result1, result2 ([Decimal('2'), Decimal('-0')], (2-0j))
result1, result2 ([Decimal('0'), Decimal('2')], 2j)
result1, result2 ([Decimal('0'), Decimal('-2')], -2j)
result1, result2 ([Decimal('5'), Decimal('5')], (5+5j))
result1, result2 ([Decimal('5'), Decimal('-5')], (5-5j))
result1, result2 ([Decimal('5'), Decimal('5')], (5+5j))
result1, result2 ([Decimal('5'), Decimal('-5')], (5-5j))
# Function ComplexSquareRoot() preserves precision of large numbers:
result1, result2 ([Decimal('123456789012345678935'), Decimal('639876543210987654321')], (1.2345678901234568e+20+6.398765432109876e+20j))