2018年12月31日 星期一

利用辛普森 理則 (Simpson's Rule) #計算 x=0 to 2 , 以曲線 y= (1+x^3) ^(1/3) 沿 z 軸旋轉一周的體積 --Python語言


# 利用辛普森 理則 (Simpson's Rule)
#計算 x=0 to 2 , 以曲線 y= (1+x^3) ^(1/3) 沿 z 軸旋轉一周的體積



from math import *

def F(x):
    temp=pow(pow((1+pow(x,3)),1.0/3),2)
    return temp

a=0.0
b=2.0
n=10
sum1=0.0
sum2=0.0
sn=0.0

m=int(n/2)

h=(b-a)/n
for i in range (1,2*m):
    x=a+i*h
    if(i%2==0):
        sum2=sum2+F(x)
    else:
        sum1=sum1+F(x)
    sn= pi *(h/3.0)*(F(a)+F(b)+2.0*sum2+4.0*sum1)
    print("S{%d}={%f}\n" % (n,sn))


===============================
輸出結果

======== RESTART: H:/2018-09勤益科大數值分析/數值分析/PYTHON/EX4-5.py ===========
S{10}={1.957852}

S{10}={2.394418}

S{10}={3.348841}

S{10}={3.900651}

S{10}={5.230509}

S{10}={6.048317}

S{10}={8.068276}

S{10}={9.308714}

S{10}={12.325078}

>>>

====================
辛普森理則的python另一寫法
====================


# Python code for simpson's 1 / 3 rule

import math

# Function to calculate f(x)
def func( x ):
    return (math.pi* math.pow((1+math.pow(x,3)),(2/3)))

# Function for approximate integral
def simpsons_( ll, ul, n ):

# Calculating the value of h
h = ( ul - ll )/n

# List for storing value of x and f(x)
x = list()
fx = list()

# Calcuting values of x and f(x)
i = 0
while i<= n:
x.append(ll + i * h)
fx.append(func(x[i]))
i += 1

# Calculating result
res = 0
i = 0
while i<= n:
if i == 0 or i == n:
res+= fx[i]
elif i % 2 != 0:
res+= 4 * fx[i]
else:
res+= 2 * fx[i]
i+= 1
res = res * (h / 3)
return res

# Driver code
lower_limit = 0.0 # Lower limit
upper_limit = 2.0 # Upper limit
n = 10 # Number of interval
print("%.6f"% simpsons_(lower_limit, upper_limit, n))


 RESTART: C:/Users/User/AppData/Local/Programs/Python/Python36-32/EX4-5-2.py 
12.325078
>>>

沒有留言:

張貼留言

Messaging API作為替代方案

  LINE超好用功能要沒了!LINE Notify明年3月底終止服務,有什麼替代方案? LINE Notify將於2025年3月31日結束服務,官方建議改用Messaging API作為替代方案。 //CHANNEL_ACCESS_TOKEN = 'Messaging ...