Hi all

I have following program to solve time delay system...

the solution is good by choosing r:=11... but for r greater than 11(e.g r=20) the system take about4-5 minute to do and then it says:"" K2:=simplify(inverse(K1)):Error, (in minor) object too large""

and for plotting says:"" Warning, unable to evaluate 1 of the 2 functions to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct""

i don't know what is the problem?

any one can help me???

best wishes

> restart:

> with(plots):

> with(linalg):

> with(LinearAlgebra):

> L:=1:

> r:=20:

> tau:=0.3:

> #definition of exact solution

>

> T1 := piecewise(t<0,0,t>=0 and t<0.3,1+t^2,t>=0.3 and t<0.6,(691/1000)+(109/100)*t+(7/10)*t^2+(1/3)*t^3,0.6<=t and t<0.9,(4409/5000)+(209/500)*t+(69/50)*t^2+(2/15)*t^3+(1/12)*t^4,0.9<=t and t<1,(1500917/2000000)+(35107/40000)*t+(1617/2000)*t^2+(87/200)*t^3+(1/120)*t^4+(1/60)*t^5):

> plot(T1,t=0..1,numpoints=10000,discont = true):

> g:=t^2:

> a[0]:=evalf((1/L)*Int(g,t=0..L)):

> for i1 from 1 to r do

> a[i1]:=evalf((2/L)*Int(g*cos(2*i1*Pi*t/L),t=0..L)):

> od:

> for j1 from 1 to r do

> b[j1]:=evalf((2/L)*Int(g*sin(2*j1*Pi*t/L),t=0..L)):

> od:

> X00:=matrix([[a[0]]]);

> X10:=matrix(1,r,0):

> for j from 1 to r do

> X10[1,j]:=a[j]:

> od:

>

> X20:=matrix(1,r,0):

> for j from 1 to r do

> X20[1,j]:=b[j]:

> od:

X00 := [0.3333333333]

> X00:=blockmatrix(1,3,[X00,X10,X20]):

>

>

>

> Z:=linalg[matrix](2*r+1,2*r+1):

> Z[1,1]:=tau:

> for iz from 2 to r+1 do

> Z[iz,1]:=(L/(2*(iz-1)*Pi))*sin(2*(iz-1)*Pi*tau/L):

> od:

> for iz from r+2 to 2*r+1 do

> Z[iz,1]:=(L/(2*(iz-1-r)*Pi))*(1-cos(2*(iz-1-r)*Pi*tau/L)):

> od:

> for jz from 2 to 2*r+1 do

> for iz from 1 to 2*r+1 do

> Z[iz,jz]:=0;

> od:

> od:

> Dtau00:=matrix([[1]]):

> Dtau01:=matrix(1,r,0):

> Dtau02:=matrix(1,r,0):

> Dtau10:=matrix(r,1,0):

> Dtau20:=matrix(r,1,0):

>

> Dtau1:=linalg[matrix](r,r):

> for i from 1 to r do

> for j from 1 to r do

> Dtau1[i,i]:=cos(2*(i)*Pi*tau/L):

> if i<>j then Dtau1[i,j]:=0 fi:

> od:

> od:

>

> Dtau2:=linalg[matrix](r,r):

> for i from 1 to r do

> for j from 1 to r do

> Dtau2[i,i]:=sin(2*(i)*Pi*tau/L):

> if i<>j then Dtau2[i,j]:=0 fi:

> od:

> od:

>

> Dtau3:=linalg[matrix](r,r):

> for i from 1 to r do

> for j from 1 to r do

> Dtau3[i,i]:=-sin(2*(i)*Pi*tau/L):

> if i<>j then Dtau3[i,j]:=0 fi:

> od:

> od:

>

> Dtau4:=linalg[matrix](r,r):

> for i from 1 to r do

> for j from 1 to r do

> Dtau4[i,i]:=cos(2*(i)*Pi*tau/L):

> if i<>j then Dtau4[i,j]:=0 fi:

> od:

> od:

>

> Dtau:=blockmatrix(3,3,[Dtau00,Dtau01,Dtau02,Dtau10,Dtau1,Dtau2,Dtau20,Dtau3,Dtau4]):

>

> P00:=matrix([[L/2]]):

> P01:=matrix(1,r,0):

>

> P02:=matrix(1,r,0):

> for j from 1 to r do

> P02[1,j]:=-L/(j*Pi):

> od:

>

> P10:=matrix(r,1,0):

>

> P20:=matrix(r,1,0):

> for i from 1 to r do

> P20[i,1]:=L/(2*i*Pi):

> od:

>

>

> P1:=linalg[matrix](r,r):

> for i from 1 to r do

> for j from 1 to r do

> P1[i,j]:=0

> od;

> od;

> P2:=linalg[matrix](r,r):

> for i from 1 to r do

> for j from 1 to r do

> P2[i,i]:=L/(2*i*Pi):

> if i<>j then P2[i,j]:=0 fi:

> od:

> od:

>

> P3:=linalg[matrix](r,r):

> for i from 1 to r do

> for j from 1 to r do

> P3[i,i]:=-L/(2*i*Pi):

> if i<>j then P3[i,j]:=0 fi:

> od:

> od:

>

> P4:=linalg[matrix](r,r):

> for i from 1 to r do

> for j from 1 to r do

> P4[i,j]:=0:

> od:

> od:

>

> P:=blockmatrix(3,3,[P00,P01,P02,P10,P1,P2,P20,P3,P4]):

> I1:=Matrix(2*r+1,shape=identity):

> K1:=simplify(evalm(I1-Dtau&*P+Dtau&*Z)):

> K2:=simplify(inverse(K1)):

Error, (in minor) object too large

>

> X0:=matrix(1,2*r+1,0):

> X0[1,1]:=1:

> for j from 2 to 2*r+1 do

> X0[1,j]:=0:

> od:

>

> X:=simplify(evalm(evalm((X0+X00))&*K2)):

>

> for h from 1 to r do

> f1(h):=cos(2*h*Pi*t/L):

> od:

> for k from 1 to r do

> f2(k):=sin(2*k*Pi*t/L):

> od:

> XP:=X[1,1]+evalf(sum(X[1,l+1]*f1(l)+X[1,r+l+1]*f2(l),l=1..r )):

> plot({XP,T1},t=0..1);

Warning, unable to evaluate 1 of the 2 functions to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct