Fortran Newton法求解非线性方程组.rar_fortran 非线性方程组_newton_牛顿迭代_牛顿迭代法_非线性

module M_GAUSS !高斯列主元消去法模块 contains subroutine LINEQ(A,B,X,N) !高斯列主元消去法 implicit real*8(A-Z) integer::I,K,N integer::ID_MAX !主元素标号 real*8::A(N,N),B(N),X(N) real*8::AUP(N,N),BUP(N) !A,B为增广矩阵 real*8::AB(N,N+1) real*8::VTEMP1(N+1),VTEMP2(N+1) AB(1:N,1:N)=A AB(:,N+1)=B !列主元消去法的核心部分 do K=1,N-1 ELMAX=DABS(AB(K,K)) ID_MAX=K !这段程序不是为了赋值最大元素给ELMAX,而是为了找出最大元素对应的标号 do I=K+1,N if(dabs(AB(I,K))>ELMAX) then ELMAX=AB(I,K) ID_MAX=I end if end do !交换两行元素，其他不变 VTEMP1=AB(K,:) VTEMP2=AB(ID_MAX,:) AB(K,:)=VTEMP2 AB(ID_MAX,:)=VTEMP1 !以上交换两行元素，然后开始消元 do I=K+1,N TEMP=AB(I,K)/AB(K,K) AB(I,:)=AB(I,:)-TEMP*AB(K,:) end do end do AUP(:,:)=AB(1:N,1:N) BUP(:)=AB(:,N+1) !调用上三角方程组的回代方法 call UP_TRI(AUP,BUP,X,N) end subroutine LINEQ subroutine UP_TRI(A,B,X,N) !上三角方程组的回代方法 implicit real*8(A-Z) integer::I,J,N real*8::A(N,N),B(N),X(N) X(N)=B(N)/A(N,N) !回代部分 do I=N-1,1,-1 X(I)=B(I) do J=I+1,N X(I)=X(I)-A(I,J)*X(J) end do X(I)=X(I)/A(I,I) end do end subroutine UP_TRI end module M_GAUSS module M_NEWTON !关于方程的模块 contains subroutine SOLVE() !计算非线性方程组 use M_GAUSS implicit real*8(A-Z) integer::I,ITMAX=50 integer::N=2 real*8::X(2),F(2),DX(2) real*8::DF(2,2) !ITMAX为最大允许迭代次数 !N为方程组维数 !DF为偏导数矩阵 open(unit=11,file="jieguo.txt") write(11,100) 100 format(/,T6,"Newton法计算非线性方程组迭代序列：",/) X=(/1.6,1.2/) TOL=1D-9 do I=1,ITMAX call FUNC(F,X) call JAC(DF,X) call LINEQ(DF,-F,DX,N) X=X+DX write(11,200)I,X 200 format(I5,2F16.10) !判断计算精度，满足误差要求时推出循环 DX2=DSQRT(DX(1)**2+DX(2)**2) if(DX2<TOL) exit end do end subroutine SOLVE subroutine FUNC(F,X) !方程函数 implicit real*8(A-Z) real*8::X(2),F(2) F(1) = X(1)**2 + X(2)**2 - 4 F(2) = X(1)**2 - X(2)**2 - 1 end subroutine FUNC subroutine JAC(DF,X) !偏导数矩阵 implicit real*8(A-Z) real*8::X(2),DF(2,2) DF(1,1) = 2*X(1) DF(1,2) = 2*X(2) DF(2,1) = 2*X(1) DF(2,2 )= -2*X(2) end subroutine JAC end module M_NEWTON program newton use M_NEWTON !调用Newton法中的函数 call SOLVE() end program newton