Skip to main content
1 of 3

Orthogonal Collocation Using Jacobi Polynomials

I'm trying to solve a PDE(diffusion-reaction in a spherical catalyst pellet) using Jacobi Orthogonal Collocation method. But at the stage of solving the system of ODEs(using NDSolve) resulting from the collocation, Mathematica does not return a result. This is the eqn:(dC/dt)=(d^2C/dr^2)+(1/r)(dC/dr)+C/(C+1)^2 Please find the code below:

Np = 10;
sol = Solve[JacobiP[Np - 2, 0, 0, 2 x - 1.] == 0, x]
tbl = Select[Table[x /. sol[[i]], {i, 1, Length[sol]}], # > 0 &];
x[Np] = 1.0;
x[1] = 10^-20; For[i = 2, i < Np, {x[i] = tbl[[i - 1]], i++}];
Q = Chop[Table[Table[x[j]^(i - 1), {i, 1, Np}], {j, 1, Np}]];
IQ = Inverse[Q] // Quiet;
Ci = Chop[Table[Table[(i - 1) x[j]^(i - 2), {i, 1, Np}], {j, 1, Np}]];
A = Ci.IQ;
Di = Chop[
Table[Table[(i - 1) (i - 2) x[j]^(i - 3), {i, 1, Np}], {j, 1, 
 Np}]];
    B = Di.IQ;
B // MatrixForm;
Z = IdentityMatrix[8];
Table[C[i][0] == 10, {i, 2, 9}]
Table[C'[i][t] == (B.Table[C[i][t], {i, 1, 10}])[[i]], {i, 2, 9}] +Z.Table[C[i][t]/(1 + C[i][t])^2, {i, 2, 9}];
 Flatten[Table[C'[i][t] == (B.Table[C[i][t], {i, 1, 10}])[[i]], {i, 2, 9}] +Z.Table[C[i][t]/(1 + C[i][t])^2, {i, 2, 9}]];
sol = NDSolve[Flatten[{Table[C'[i][t] == (B.Table[C[i][t], {i, 1, 10}])[[i]], {i, 2, 9}] +Z.Table[C[i][t]/(1 + C[i][t])^2, {i, 2, 9}], C'[0][t] == 0,C[1][t] == 1, Table[C[i][0] == 10, {i, 2, 9}]}],Table[C[i][t], {i, 1, 10}],{t, 0, 10}];