Hi all
I am going to implement exercise below in fortran:
i wrote code below :
program main
implicit none
integer :: k
real :: a_1, a_2, a_3, Q, R, discriminant, theta, s, x_1
complex :: z_2, z_3, omega
real , dimension(3) :: x_k
real, parameter :: pi = 3.1415
print *, "Please enter coefficients for cubic equation:"
print *, "Cubic equation -> x^3 + a_1 * x^2 + a_2 * x + a_3 = 0"
read *, a_1, a_2, a_3
Q = (a_1 ** 2 - 3 * a_2) / 9
R = ((2 * a_1 ** 3) - 9 * a_1 * a_2 + 27 * a_3) / 54
discriminant = Q ** 3 - R ** 2
if (discriminant >= 0) then
theta = acos(R / sqrt(Q ** 3))
do k = 1, 3
x_k(k) = (-2 * sqrt(Q)) * cos((theta + 2 * pi * (k - 1)) / 3) - (a_1 / 3)
end do
print *, "Real roots are:", x_k
else
s = -1 * sign(R, R) * (abs(R) + sqrt(R ** 2 - Q ** 3))
omega = (-1, 1.7320) / 2
x_1 = s + Q / s - a_1 / 3
z_2 = (omega * s) + (omega ** 2 * (Q / s)) - a_1 / 3
z_3 = (omega ** 2 * s) + (omega * Q / s) - a_1 / 3
print *, "Real root is:", x_1, "Complex roots are:", z_2, z_3
end if
end program main
For a_1 = 1, a_2 = 2, a_3 = 3 ,
it should print
x_1 = -1.27568
z_2 = 0.13784 + i * 1.52731
z_3 = 0.13784 - i * 1.52731
But it prints
x_1 = -3.12811947
z_2 = ( 1.06406796 , -2.74305439 )
z_3 = ( 1.06392860 , 2.74305439 )