Я пытаюсь написать программу трапециевидного правила с использованием maple

Question

Я пытаюсь написать программу трапецеидального правила с использованием клена, я отправлю свой код, но проблема в том, что когда я прочитаю эту статью, я обнаружил, что сложно запрограммировать все числовые примеры, может ли кто-нибудь мне помочь?

статья по этому адресу: http://www.cs.technion.ac.il/~asidi/Sidi_Journal_Papers/P097_JCAM.sinpq.trans.pdf код моего Maple в:

with(linalg); with(LinearAlgebra); with(plots); with(MTM); 
a := 0; b := 1; p := 20; 
k := 1; R := vector(p); 
f := proc (x) options operator, arrow; x*(1-x)/(1+x) end proc; 
Uex := int(f(x), x = a .. b); 
Uex := evalf(Uex); 
for N to p do 
    h := (b-a)/N; 
    U[N] := (f(a)+f(b))*(1/2); 
    for j to N-1 do 
     U[N] := U[N]+f(a+j*(b-a)/N) 
    end do; 
    U[N] := h.U[N] 
end do; 
for N to p do 
    R[N] := -evalf(log(abs(Uex-U[N]))/log(10)); 
    print(N, evalf(Uex-U[N])) 
end do; 
plot(f); listplot(R) 

Answer 1

Попробуйте это:

a:= 0: b:= 1: p:= 20: 
R:= Vector(p): H:= Vector(p): 
f:= x-> x*(1-x)/(1+x); 
Uex:= int(f(x), x= a..b); 
for n to p do 
    N:= 2^n; 
    h:= (b-a)/N; 
    H[n]:= h; 
    R[n]:= evalhf(abs(Uex-h*((f(a)+f(b))/2 + add(f(a+j*h), j= 1..N-1))))  
end do: 
unassign('h'); 
err = Statistics:-PowerFit(H,R,h);