! This is a basic starting shell for your programs MODULE const INTEGER, PARAMETER :: dp = KIND(1.0d0), sp = KIND(1.0e0) REAL(KIND=dp), PARAMETER :: dpi = 3.141592653589793d0 END MODULE const ! PROGRAM prog USE const IMPLICIT NONE INTEGER :: n REAL(kind=dp) :: x,y,left,right,error REAL(kind=dp) :: binexp WRITE(*,*) 'Input n,x,y' READ(*,*) n,x,y left = (x+y)**n right = binexp(n,x,y) error = (right-left)/left WRITE(*,*) 'left,right = ',left,right WRITE(*,*) 'error = ',error END PROGRAM prog !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! FUNCTION binexp(n,x,y) USE const IMPLICIT NONE INTEGER,INTENT(IN) :: n REAL(kind=dp),INTENT(IN) :: x,y REAL(kind=dp) :: binexp REAL(kind=dp) :: bincoeff INTEGER :: i binexp = 0 DO i = 0,n binexp = binexp + bincoeff(n,i) * x**i * y**(n-i) ENDDO RETURN END FUNCTION binexp !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! FUNCTION bincoeff(n,k) USE const IMPLICIT NONE INTEGER,INTENT(IN) :: n,k REAL(kind=dp) :: bincoeff REAL(kind=dp) :: factfun bincoeff = factfun(n)/(factfun(k)*factfun(n-k)) RETURN END FUNCTION bincoeff !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! FUNCTION factfun(n) USE const IMPLICIT NONE INTEGER,INTENT(IN) :: n REAL(kind=dp) :: factfun INTEGER :: i factfun = 1 DO i = 1,n factfun = factfun*i ENDDO END FUNCTION factfun ! Input n,x,y ! 6,1.3,2.4 ! left,right = 2565.7264090000008 2565.7264090000003 ! error = -1.7723922133369754E-016 ! ! Input n,x,y ! 20,1.0000001,-1.0000002 ! left,right = 9.9999996726842304E-141 5.2493120961116801E-011 ! error = 5.2493122679299485E+129