I  encountered a non-integrable integral in the process of solving the following process, . How to achieve its numerical solution? Such as in a looping   code：

#######
pa[i] := pa[i-1]-(Int(subs(t = tau, Lpa[i-1]+Na1[i-1]-Na2[i-1]), tau = 0 .. t)); 

pw[i] := pw[i-1]-(Int(subs(t = tau, Lpw[i-1]+Nw1[i-1]-Nw2[i-1]), tau = 0 .. t)); u[i] := u[i-1]-(Int(subs(t = tau, Lu[i-1]+Nu1[i-1]+Nu2[i-1]), tau = 0 .. t))；

######
Detailed code see annexBC2.mw

plz help！ How to draw the three-dimensional graphics (p(x,t))? When I run,it can not run.I do not know where the problem lies.The code is as follows：

restart:
with(PDEtools); with(student); KN := 3;
C2 := 1/.3; C1 := 0.6e-2/(.3); C3 := 4.3/(.3); beta := 0.43e-1;
ADM1 := proc (n) options operator, arrow; convert(subs(lambda = 0, value((Diff(F(Sum(lambda^i*U[i], i = 0 .. n)), `$`(lambda, n)))/factorial(n))), diff) end proc; A0[0] := F(U[0]);
for n to KN do A0[n] := ADM1(n) end do;
for n from 0 to KN do A[n] := unapply(simplify(convert(C1*(diff(subs({seq(U[i] = p[i](x, t), i = 0 .. KN)}, expand(subs(F(U[0]) = U[0]*exp(U[0]), A0[n]))), t)), diff)+C2*convert(subs({seq(U[i] = Diff(p[i](x, t), x), i = 0 .. KN)}, expand(subs(F(U[0]) = U[0]^2, A0[n]))), diff)), x, t) end do;
p[0] := proc (x, t) options operator, arrow; .2*sin((1/2)*Pi*x)*exp((-1/4)*t*Pi*Pi/C1) end proc;
p[1] := proc (x, t) options operator, arrow; (-(int(subs(t = s, C3*(diff(p[0](x, t), x))), s = 0 .. t))+int(subs(t = s, A[0](x, t)), s = 0 .. t))/beta end proc;
p[2] := proc (x, t) options operator, arrow; -(int(subs(t = s, C3*(diff(p[1](x, t), x))), s = 0 .. t))+int(subs(t = s, A[1](x, t)), s = 0 .. t) end proc;
p := unapply(subs(x = Zeta, t = tau, p[0](x, t)+p[1](x, t)+p[2](x, t)), Zeta, tau);
with(plots);
plot3d(p(x,t),x=0..1,t=0..1)

Plz help me！ How to accelerate the calculation speed the following the program？

restart;
with(PDEtools); declare((u, W)(x, t)); KN := 10;
AFP := proc (C1, C2, C3, C4, H, KN, N) local ADM1, n, lambda, F, i, A0, A, U, W, u, L, R, NL, w, PDE1, IC1, d, Eq1, Eq2, Eq3, LI, trL, tr1, trN, Apr, AprK, trSol, Sol, AD; declare((u, W)(x, t)); ADM1 := proc (n) options operator, arrow; convert(subs(lambda = 0, value((Diff(F(Sum(lambda^i*U[i], i = 0 .. n)), `$`(lambda, n)))/factorial(n))), diff) end proc; A0[0] := F(U[0]); for n to KN do A0[n] := ADM1(n) end do; for n from 0 to KN do A[n] := convert(C1*(diff(subs({seq(U[i] = W[i](x, t), i = 0 .. KN)}, expand(subs(F(U[0]) = U[0]*exp(U[0]), A0[n]))), x)), diff)+convert(C2*subs({seq(U[i] = Diff(W[i](x, t), x), i = 0 .. KN)}, expand(subs(F(U[0]) = U[0]^2, A0[n]))), diff) end do; L := proc (w) options operator, arrow; diff(w(x, t), t) end proc; R := proc (w) options operator, arrow; C3*(diff(w(x, t), x)) end proc; NL := proc (w) options operator, arrow; C1*(diff(w(x, t)*exp(w(x, t)), t))-C2*(diff(w(x, t), x))^2 end proc; PDE1 := proc (w) options operator, arrow; L(w)-R(w) = -NL(w) end proc; IC1 := u(x, 0) = sum(2*(int(sin((d+1/2)*Pi*x/H), x = 0 .. H))*exp(-C4*(d+1/2)^2*Pi^2*t/H^2)*sin((d+1/2)*Pi*x/H)/H, d = 0 .. N); LI := proc (w) options operator, arrow; Int(w(x, t), t = 0 .. t) end proc; tr1 := u-rhs(IC1); Eq1 := LI(lhs(PDE1(u))) = LI(rhs(PDE1(u))); Eq2 := simplify(subs(lhs(Eq1) = tr1, Eq1)); trL := u = add(u[j](x, t), j = 0 .. KN); trN := LI(NL(u)) = Int(Sum(A[i], i = 0 .. KN), t = 0 .. t); Eq3 := subs(trL, lhs(Eq2)) = subs(trN, rhs(Eq2)); Apr[0] := u[0](x, t) = rhs(IC1); AprK := u[k+1](x, t) = -(Int(AD[k], t = 0 .. t)); for i from 0 to KN do Apr[i+1] := value(subs({seq(Apr[m], m = 0 .. i)}, subs({seq(W[m] = u[m], m = 0 .. i)}, subs(k = i, AD[i] = A[i], AprK)))) end do; trSol := {seq(Apr[i], i = 0 .. KN)}; value(subs(trSol, trL)) end proc;

Download AFP.mw

Please help me.I don't know how to achieve the following iteration relation by maple code.

hi,i am studying the maple most recent.But when calculating function integral，I ran into trouble.I hope to get your help.Here is the code I wrote, but it runs a very long time. How to effectively reduce the integration time?

restart;
with(student);
assume(n::integer);
Fourierf := proc (sigma, a, b, N) local A, A0, B, T, S, Ff; T := b-a; A0 := int(sigma, t = a .. b); A := int(sigma*sin(n*Pi*t/T), t = a .. b); B := int(sigma*cos(n*Pi*t/T), t = a .. b); S := sum(A*sin(n*Pi*t/T)+B*cos(n*Pi*t/T), n = 1 .. N)+(1/2)*A0; Ff := unapply(S, t) end proc;

f := proc (t) options operator, arrow; piecewise(t < .13*2.6 and 0 <= t, 100*t/(.13*2.6), .13*2.6 <= t and t < 2.6, 100, 2.6 <= t and t < 2.6*1.1, 0) end proc;

sigma := f(t);
a := 0;
b := 1.1*2.6;
s1 := unapply((Fourierf(sigma, a, b, 500))(t)/uw0, t);

s2 := unapply((Fourierf(sigma, a, b, 500))(t)/ua0, t);
A1 := (2*n+1)^2*Pi^2*(C3+1+sqrt(4*C1*C2*C3+C3^2-2*C3+1))/(8*C1*C2-8);
A2 := (2*n+1)^2*Pi^2*(C3+1-sqrt(4*C1*C2*C3+C3^2-2*C3+1))/(8*C1*C2-8);
g := -C2*Cww*(diff(s1(x), `$`(x, 2)))+Caa*(diff(s2(x), `$`(x, 2))+(n+1/2)^2*Pi^2*(diff(s2(x), x)));
f1 := -(1/2)*(n+1/2)^2*Pi^2*sqrt(4*C1*C2*C3+C3^2-2*C3+1)+C2*Cww*((D@@1)(s1))(0)-Caa*((D@@1)(s2))(0)+(n+1/2)^2*Pi^2*(C2-(1/2)*C3+1/2);

CN := ((2*(int(exp(-A1*x)*g, x = 0 .. t)-f1))*exp(A1*t)-(2*(int(exp(-A2*x)*g, x = 0 .. t)-f1))*exp(A2*t))/((n+1/2)^3*Pi^3*sqrt(4*C1*C2*C3+C3^2-2*C3+1));
ua := sum(CN*sin((n+1/2)*Pi*z), n = 0 .. 100);