# 例題 4-8 利用 梯形理則 計算 雙重積分 f(x,y)= 4xy 在[0 , 1-x^2] dy 與 [0,1] dx 的定積分

# 例題 4-8
# 利用 梯形 (trapezoidal Rule ) 計算 雙重積分
# f(x,y)= 4xy  在[0 , 1-x^2] dy  與 [0,1] dx 的定積分

# Based on Trapezoidal Rule to  compute the double integral.
def F(x,y):
    return (4*x*y + 0.0 )
def C(x):
    return (0.0)
def D(x):
    return (1-pow(x,2) + 0.0)

def gy(n):
    sum=0.0;
    for i in range (0,n+1):
        for j in range (1,n):
            sum=sum+F(x[i],y[i][j])

        gg[i]=(hy[i]/2.0)*(F(x[i],y[i][0])+F(x[i],y[i][n])+2*sum)
        sum=0.0
       
#======================================
x=[0.0 for i in range(30)]
y=[ [0.0 for i in range(30)] for j in range (30)]
#print(y)
gg= [0.0 for i in range(30) ]
hy= [0.0 for i in range(30) ]

sum1=0.0
n=10
a=0.0
b=1.0
hx=(b-a)/n

for i in range (0 ,n+1) :
    x[i]=a+i*hx
    #print('x[%2d]={%.6f}' %(i,x[i]))
    hy[i]=(D(x[i])-C(x[i]))/n
    #print('hy[%2d]={%.6f}' %(i,hy[i]))
    for j in range (0,n+1):
        y[i][j]=C(x[i]) + j*hy[i] + 0.0
        #print('y[%2d][%2d]={%.6f}' %(i,j,y[i][j]),end='' )
    #print()
gy(n)

print('i={%2d} gg[%2d]={%.6f}' %(0,0,gg[0]))
for i in range (1,n):
    sum1=sum1+gg[i]
    print('i={%2d} gg[%2d]={%.6f}' %(i,i,gg[i]))

print('i={%2d} gg[%2d]={%.6f}' %(n,n,gg[n]))
ts=(hx/2.0)*(gg[0]+gg[n]+2*sum1)
print("T{%d}={%.6f}\n" %(n,ts))


#=============================
# n= 10 --> n=20
#=============================
sum1=0.0
n=29
a=0.0
b=1.0
hx=(b-a)/n

for i in range (0 ,n+1) :
    x[i]=a+i*hx
    #print('x[%2d]={%.6f}' %(i,x[i]))
    hy[i]=(D(x[i])-C(x[i]))/n
    #print('hy[%2d]={%.6f}' %(i,hy[i]))
    for j in range (0,n+1):
        y[i][j]=C(x[i]) + j*hy[i] + 0.0
        #print('y[%2d][%2d]={%.6f}' %(i,j,y[i][j]),end='' )
    #print()
gy(n)

print('i={%2d} gg[%2d]={%.6f}' %(0,0,gg[0]))
for i in range (1,n):
    sum1=sum1+gg[i]
    print('i={%2d} gg[%2d]={%.6f}' %(i,i,gg[i]))

print('i={%2d} gg[%2d]={%.6f}' %(n,n,gg[n]))
ts=(hx/2.0)*(gg[0]+gg[n]+2*sum1)
print("T{%d}={%.6f}\n" %(n,ts))

========= RESTART: F:/2018-09勤益科大數值分析/數值分析/PYTHON/EX4-8-1.py ==========
i={ 0} gg[ 0]={0.000000}
i={ 1} gg[ 1]={0.196020}
i={ 2} gg[ 2]={0.368640}
i={ 3} gg[ 3]={0.496860}
i={ 4} gg[ 4]={0.564480}
i={ 5} gg[ 5]={0.562500}
i={ 6} gg[ 6]={0.491520}
i={ 7} gg[ 7]={0.364140}
i={ 8} gg[ 8]={0.207360}
i={ 9} gg[ 9]={0.064980}
i={10} gg[10]={0.000000}
T{10}={0.331650}

i={ 0} gg[ 0]={0.000000}
i={ 1} gg[ 1]={0.068802}
i={ 2} gg[ 2]={0.136622}
i={ 3} gg[ 3]={0.202492}
i={ 4} gg[ 4]={0.265465}
i={ 5} gg[ 5]={0.324631}
i={ 6} gg[ 6]={0.379126}
i={ 7} gg[ 7]={0.428143}
i={ 8} gg[ 8]={0.470947}
i={ 9} gg[ 9]={0.506885}
i={10} gg[10]={0.535398}
i={11} gg[11]={0.556029}
i={12} gg[12]={0.568443}
i={13} gg[13]={0.572429}
i={14} gg[14]={0.567920}
i={15} gg[15]={0.555000}
i={16} gg[16]={0.533915}
i={17} gg[17]={0.505088}
i={18} gg[18]={0.469130}
i={19} gg[19]={0.426851}
i={20} gg[20]={0.379269}
i={21} gg[21]={0.327627}
i={22} gg[22]={0.273400}
i={23} gg[23]={0.218312}
i={24} gg[24]={0.164340}
i={25} gg[25]={0.113733}
i={26} gg[26]={0.069021}
i={27} gg[27]={0.033025}
i={28} gg[28]={0.008870}
i={29} gg[29]={0.000000}
T{29}={0.333135}

>>>