代码语言：fortran

所属分类：

代码描述：fortran编写一个盛金公式解一元三次方程代码

代码标签： 盛金 公式 一元 三次 方程

下面为部分代码预览，完整代码请点击下载或在bfwstudio webide中打开

Module Shengjin_mod
  Implicit None

contains

  Function Cubic_equation(Co) Result (X)
!盛金公式求解一元三次方程
!默认浮点为8字节,以保证精度
    Integer , Parameter :: P = Selected_Real_Kind(12)
    Real (Kind=P), Parameter :: Eps = 1.0_P
    Real (Kind=P), Intent (In) :: Co(4)
    Real (Kind=P) :: A, B, C, D, A0, B0, C0, D0, Y1, Y2
    Complex (Kind=P) :: X(3)
    A = Co(1)
    B = Co(2)
    C = Co(3)
    D = Co(4)
    X = Huge(A)/100.0_P
!非三次方程
    If (Abs(A) < Eps) Return
    A0 = B*B - 3.0_P*A*C
    B0 = B*C - 9.0_P*A*D
    C0 = C*C - 3.0_P*B*D
    D0 = B0*B0 - 4.0_P*A0*C0

    If (Abs(A0) < Eps .And. Abs(B0) < Eps) Then
!三重实根
      X(1:3) = -B/(3.0_P*A)
    Else
      If (D0 > 0.0_P) Then
        Y1 = (A0*B+1.5_P*A*(-B0+Sqrt(D0)))
        If (Y1 < 0.0_P) Then
          Y1 = -Abs(Y1)**(1.0_P/3.0_P)
        Else
          Y1 = Y1**(1.0_P/3.0_P)
        End If
        Y2 = (A0*B+1.5_P*A*(-B0-Sqrt(D0)))
        If (Y2 < 0.0_P) Then
          Y2 = -Abs(Y2)**(1.0_P/3.0_P)
        Else
          Y2 = Y2**(1.0_P/3.0_P)
        End If
!一个实根，一对共轭复根
        X(1) = -(B+Y1+Y2)/(3.0_P*A)
        X(2) = Cmplx(-2.0_P*B+Y1+Y2, Sqrt(3.0_P)*(Y1-Y2))/(6.0_.........完整代码请登录后点击上方下载按钮下载查看

网友评论0

取消

评论提交中...

提交评论