Prev: simulink image process
Next: what distance metric is this?
From: kp pk on 27 Jun 2010 09:30
I am trying to generate orthonormal polynomials of the form
1 - P(1)
a+bx - P(2)
c+dx+ex^2 - P(3)

and so on..

the way in which i am generating above polynomials is that
integral(P(i)*P(j)) x varying from 0 to 1 = 0 if i is not equal to j
= 1 if i is equal to j
the code i am using to do this is

n=20;-specifies the order of polynomials i will generate till
d=(n*(n+1)/2)-1;
syms x;
A = sym(zeros(1, d));
A(1,1)=1;
for i=2:n
for k=1:i
A(i,k)=sym(sprintf('m%d',i,k));
end
end /* -in this part,i just created a matrix with the coefficients of the 'x' in the polynomials as the symbolic entries
so the matrix is
1
m2m1 m2m2
m3m1 m3m2 m3m3..and so on.
observe that in order to get m3m1 and m3m2 and m3m3 i must first find out the values of m2m1 and m2m2 and substitute them into P(2) */
P(1)=sym(1);
for i=2:n
P(i)=sym(0);
for k=1:i
P(i)=P(i)+A(i,k)*(x^(k-1));
end-here i have created the polynomials
for m=1:i
temp(m)=0;
temp(m)=int(P(m)*P(i),x,0,1); -here i have created the expressions
end
S(i)=solve(temp
/* this is where i get a problem,suppose i am solving for m2m1 then i need to write
S=solve(temp(1),temp(2));
and then use S.m2m1 and S.m2m2 in order to find out their values and then substitute those values back into the polynomial.now obviously if i write this in the loop,it will fail for the next iteration viz when i=3.
so is there a way in which i could write S=solve(temp) ,i.e.solve the entire array of equations temp without writing out temp(1),temp(2) and so on.
or if this isnt possible do i have to manually solve the set of equations without using a inbuilt function?
*/

thanks
From: John D'Errico on 27 Jun 2010 09:46
"kp pk" <krishna123321(a)gmail.com> wrote in message <i07jpf$53d$1(a)fred.mathworks.com>...
> I am trying to generate orthonormal polynomials of the form
> 1 - P(1)
> a+bx - P(2)
> c+dx+ex^2 - P(3)
>
> and so on..
>
> the way in which i am generating above polynomials is that
> integral(P(i)*P(j)) x varying from 0 to 1 = 0 if i is not equal to j
> = 1 if i is equal to j
> the code i am using to do this is

Why write it yourself like this?

Why not just recognize that these are Legendre
polynomials, shifted to be orthogonal over the
interval [0,1] instead of [-1,1].

Given that, there are simple methods to generate
Legendre polynomials.

Why expend all of this effort to do something that
is trivially done with other methods?

John
From: kp pk on 27 Jun 2010 10:01
"John D'Errico" <woodchips(a)rochester.rr.com> wrote in message <i07kms$3ns$1(a)fred.mathworks.com>...
> "kp pk" <krishna123321(a)gmail.com> wrote in message <i07jpf$53d$1(a)fred.mathworks.com>...
> > I am trying to generate orthonormal polynomials of the form
> > 1 - P(1)
> > a+bx - P(2)
> > c+dx+ex^2 - P(3)
> >
> > and so on..
> >
> > the way in which i am generating above polynomials is that
> > integral(P(i)*P(j)) x varying from 0 to 1 = 0 if i is not equal to j
> > = 1 if i is equal to j
> > the code i am using to do this is
>
> Why write it yourself like this?
>
> Why not just recognize that these are Legendre
> polynomials, shifted to be orthogonal over the
> interval [0,1] instead of [-1,1].
>
> Given that, there are simple methods to generate
> Legendre polynomials.
>
> Why expend all of this effort to do something that
> is trivially done with other methods?
>
> John

well,the only shifted legendre polynomials over [0.1] i could find were lists of polynomials .i couldnt find any code in MATLAB to generate them.i would be grateful if you could give me a good link though.
thanks
From: kp pk on 27 Jun 2010 10:05
"John D'Errico" <woodchips(a)rochester.rr.com> wrote in message <i07kms$3ns$1(a)fred.mathworks.com>...
> "kp pk" <krishna123321(a)gmail.com> wrote in message <i07jpf$53d$1(a)fred.mathworks.com>...
> > I am trying to generate orthonormal polynomials of the form
> > 1 - P(1)
> > a+bx - P(2)
> > c+dx+ex^2 - P(3)
> >
> > and so on..
> >
> > the way in which i am generating above polynomials is that
> > integral(P(i)*P(j)) x varying from 0 to 1 = 0 if i is not equal to j
> > = 1 if i is equal to j
> > the code i am using to do this is
>
> Why write it yourself like this?
>
> Why not just recognize that these are Legendre
> polynomials, shifted to be orthogonal over the
> interval [0,1] instead of [-1,1].
>
> Given that, there are simple methods to generate
> Legendre polynomials.
>
> Why expend all of this effort to do something that
> is trivially done with other methods?
>
> John


the only shifted legendres over [0,1] i could find were lists of polynomials,i couldnt find any code to generate them.i would be grateful if you could give me a good link though.
thanks.
From: John D'Errico on 27 Jun 2010 11:23
"kp pk" <krishna123321(a)gmail.com> wrote in message <i07lr1$fpg$1(a)fred.mathworks.com>...
> "John D'Errico" <woodchips(a)rochester.rr.com> wrote in message <i07kms$3ns$1(a)fred.mathworks.com>...
> > "kp pk" <krishna123321(a)gmail.com> wrote in message <i07jpf$53d$1(a)fred.mathworks.com>...
> > > I am trying to generate orthonormal polynomials of the form
> > > 1 - P(1)
> > > a+bx - P(2)
> > > c+dx+ex^2 - P(3)
> > >
> > > and so on..
> > >
> > > the way in which i am generating above polynomials is that
> > > integral(P(i)*P(j)) x varying from 0 to 1 = 0 if i is not equal to j
> > > = 1 if i is equal to j
> > > the code i am using to do this is
> >
> > Why write it yourself like this?
> >
> > Why not just recognize that these are Legendre
> > polynomials, shifted to be orthogonal over the
> > interval [0,1] instead of [-1,1].
> >
> > Given that, there are simple methods to generate
> > Legendre polynomials.
> >
> > Why expend all of this effort to do something that
> > is trivially done with other methods?
> >
> > John
>
>
> the only shifted legendres over [0,1] i could find were lists of polynomials,i couldnt find any code to generate them.i would be grateful if you could give me a good link though.
> thanks.

Download my sympoly tools from the FEX.

http://www.mathworks.com/matlabcentral/fileexchange/9577

In there, you will find orthpoly. For example,

for i = 0:5
% standard legendre polynomials
P(i + 1) = orthpoly(i,'leg');
end

P =
Sympoly array has size = [1 11]
Sympoly array element [1 1]
1
Sympoly array element [1 2]
x
Sympoly array element [1 3]
-0.5 + 1.5*x^2
Sympoly array element [1 4]
-1.5*x + 2.5*x^3
Sympoly array element [1 5]
0.375 - 3.75*x^2 + 4.375*x^4
Sympoly array element [1 6]
1.875*x - 8.75*x^3 + 7.875*x^5

These are orthogonal over the interval [-1,1], although
they are not orthonormal. First, change the interval.

sympoly y
Q = subs(Q,'x',2*y - 1)
Q =
Sympoly array has size = [1 6]

Sympoly array element [1 1]
1
Sympoly array element [1 2]
-1 + 2*y
Sympoly array element [1 3]
1 - 6*y + 6*y^2
Sympoly array element [1 4]
-1 + 12*y - 30*y^2 + 20*y^3
Sympoly array element [1 5]
1 - 20*y + 90*y^2 - 140*y^3 + 70*y^4
Sympoly array element [1 6]
-1 + 30*y - 210*y^2 + 560*y^3 - 630*y^4 + 252*y^5


These polynomials are orthogonal over [0,1]. It is
easy to convince yourself of that. For example...

defint(Q(2)*Q(6),'y',[0 1])
ans =
0

But they are not orthonormal.

defint(Q(6)*Q(6),'y',[0 1])
ans =
0.090909

So make them so.

for i = 0:5
Q(i+1) = Q(i+1)/sqrt(defint(Q(i+1)*Q(i+1),'y',[0 1]));
end

Are the polynomials now orthonormal? Yes.

defint(Q(6)*Q(6),'y',[0 1])
ans =
1

HTH,
John
 |  Next  |  Last
Pages: 1 2
Prev: simulink image process
Next: what distance metric is this?