Need help with this code- Chevyshev rational function on matlab

Question

sagar pokhrel le 23 Déc 2014

0
Lien

Utiliser le lien direct vers cette question

https://fr.mathworks.com/matlabcentral/answers/167906-need-help-with-this-code-chevyshev-rational-function-on-matlab

Commenté : John D'Errico le 24 Déc 2014

I am trying to approximate some data. Matlab couldn't approximate polynomial so I approximated using another software which can generate matlab file. This is the actual function.

%F(x,y)=(a+dT1(x )+eT1(y )+hT2(x)+iT2(y)+lT3(x )+mT3(y )+pT4(x)+qT4(y )+tT5(x )+uT5(y )... +abT6(x)+acT6(y)+afT7(x)+agT7(y )+ajT8(x )+akT8(y )+anT9(x)+aoT9(y)+arT10(x)+... asT10(y))/((1+bT1(x)+cT1(y)+fT2(x)+gT2(y )+jT3(x)+kT3(y )+nT4(x)+oT4(y )+... rT5(x )+sT5(y)+vT6(x)+aaT6(y)+adT7(x)+aeT7(y)+ahT8(x)+aiT8(y)+alT9(x)+amT9(y 1)+... apT10(x)+aqT10(y ) )+atT11(x )+auT11(y ))

% where Tn(x) = 2xTn-1(x) - Tn-2(x), n≥2 is chebysheb rational function

I do understand in this code is I need to provide input x and y and matrix c is the coefficient of above function. I get lost in function z = evalcratl(order, logx, logy, x, y, p,... s0, s1, s2, s3, s4, s5, s6, s7, s8, s9). How is this function working. Meaning of Evalcratl and tcnt .

Any help on this would be appreciated.

--------------------------------------------------------------------------------

 [rowx colx]=size(xa);
if(rowx~=1 & colx~=1)
  error('x must be scalar or 1D array');
  return;
end
[rowy coly]=size(ya);
if(rowy~=1 & coly~=1)
  error('y must be scalar or 1D array');
  return;
end
c=[
  0.7823326556648382,
  -0.003073888389053987,
  6.311956873270135E-16,
  -0.0003235701036352635,
  0,
  0.1184529306672700,
  0.7000618410087535,
  -0.06699883001305451,
  0.6882427553087074,
  -0.001289845952369463,
  7.834476942891433E-16,
  0.0003263428678678636,
  0,
  -0.01433571833230648,
  -0.2173515431820155,
  0.01552171410016730,
  -0.2211421352264687,
  -0.0001289154504598550,
  6.559509737139252E-16,
  0.0004529420876541368,
  0,
  -0.01218726621626016,
  -0.09562258140958195,
  0.003063157214620869,
  -0.08732570079084124,
  0.0002432826238712180,
  4.144106621110242E-16,
  0.0003501722099832917,
  0,
  -0.008018886525623193,
  -0.04998847489715176,
  0.003499439465963643,
  -0.04330683933474295,
  6.716713718394032E-05,
  1.191551227290313E-16,
  8.914863235537785E-05,
 0,
  0.001270750287501435,
  -0.0001942757598005615,
  ];
lenx=length(xa);
leny=length(ya);
for(j=1:leny)
  for(i=1:lenx)
    x=xa(i);
    y=ya(j);
    z=evalcratl(38,0,0,x,y,c,...
      0.000000000000000,5.000000000000000,...
      0.000000000000000,0.000000000000000,...
      0.000000000000000,5.000000000000000,...
      0.000000000000000,0.000000000000000,...
      0.1962500000000000,0.8042500000000000);
     za(i,j)=z;
     end
   end
%--------------------------------------------------------------
function z = evalcratl(order, logx, logy, x, y, p,...
s0, s1, s2, s3, s4, s5, s6, s7, s8, s9)
%--------------------------------------------------------------
tx=[];
ty=[];
if(logx==0)
  x=(x-s0)/s1;
else
  x=(log(x)-s2)/s3;
end
if(logy==0)
  y=(y-s4)/s5;
else
  y=(log(y)-s6)/s7;
end
if (order==6)
  tcnt=3;
elseif (order==10)
  tcnt=4;
elseif (order==14)
  tcnt=5;
elseif (order==18)
  tcnt=6;
elseif (order==22)
  tcnt=7;
elseif (order==26)
  tcnt=8;
elseif (order==30)
  tcnt=9;
elseif (order==34)
  tcnt=10;
elseif (order==38)
  tcnt=11;
elseif (order==42)
  tcnt=12;
else
  return;
end
if(tcnt>7)
  if(x<-1.0)
    x=-1.0;
  end
  if(x>1.0)
    x= 1.0;
  end
  if(y<-1.0)
    y=-1.0;
  end
  if(y>1.0)
    y= 1.0;
  end
end
tx(1)=1.0;
ty(1)=1.0;
tx(2)=x;
ty(2)=y;
for j=2:1:tcnt-1
  tx(j+1)=2*x*tx(j)-tx(j-1);
  ty(j+1)=2*y*ty(j)-ty(j-1);
end
m=2;
num=p(1);
den=1.0+p(2)*tx(m)+p(3)*ty(m);
for(j=3:4:order-1)
  num=num+p(j+1)*tx(m);
  num=num+p(j+2)*ty(m);
  m=m+1;
  den=den+p(j+3)*tx(m);
  den=den+p(j+4)*ty(m);
end
if (den==0)
  z=0;
else
  z=(num/den)*s8+s9;
end
return;

Thank You in advance

3 commentaires
Afficher 1 commentaire plus ancienMasquer 1 commentaire plus ancien

sagar pokhrel le 23 Déc 2014

I tried to approximate using matlab but it could approximate polynomial only in X^5 and Y^5. I was looking for a approximation that can go to x^10 and Y^10. Is there other way I can approximate to this limit.

Matlab could run code after the inputs are given. Since I know the equation and coefficient I want to use the first equation but I dont know how this code is working.

John D'Errico le 24 Déc 2014

Again, those limits were more about your understanding of the numerical analysis involved than about the true capabilities of a numerical tool.

As far as an explanation of the code you have been given, I don't see much purpose to it. Sorry, but reading undocumented code that does something you don't understand is a losing proposition. Either trust the provider or don't use it.

Connectez-vous pour commenter.

Connectez-vous pour répondre à cette question.