Re: Empirical othogonal functions and pure functions
- To: mathgroup at smc.vnet.net
- Subject: [mg84836] Re: Empirical othogonal functions and pure functions
- From: dh <dh at metrohm.ch>
- Date: Wed, 16 Jan 2008 03:28:03 -0500 (EST)
- References: <fmhr82$c60$1@smc.vnet.net>
Hi,there are several issues here:
-Function has the attribute HoldAll, therefore you must wrap its
argument into Evaluate:Table[Function[x,Evaluate[x^i]],{i,0,4}] However,
you do not need this anyway.
-the input to Orthogonalize are expressions, not function.
-the inner product function takes 2 expressions and returns a number
- Orthogonalize returns a linear combinations of the input expressions
Here is an example. As X we take Chebyshev points, therefore we should
get Chebyshev polynomials (up to norm):
X=Solve[ChebyshevT[4,x]==0,x][[All,1,2]];
EmpirScal[f_,g_]:=Dot[Table[f/.x->X[[i]],{i,Length@X}],Table[g/.x->X[[i]],{i,Length@X}]];
Orthogonalize[funs,EmpirScal]//Expand//Simplify
hope this helps, Daniel
mante wrote:
> Hello!
>
> I need to othogonalize polynomials relatively to a fixed sampling
> mesh {x1, ............., xp}, ie. to the *discrete* scalar product :
> <P,Q>=Sum[P(xi) Q(xi),{i,p}].
> This scalar product corresponds to the function :
> EmpirScal[f_, g_, X_] := Dot[Table[f[X[[i]]], {i, Length@X}],
> Table[g[X[[i]]], {i, Length@X}]];
>
> Let's try with
> rr = Table[x^i, {i, 0, 2}]
> using first the "canonical" orthogonalization procedure :
> Orthogonalize[rr, Integrate[#1*#2, {x, -1, 1}] & ]
> ->{1/Sqrt[2], Sqrt[3/2]*x, (3/2)*Sqrt[5/2]*(-(1/3) + x^2)}
> Fine!
>
> Now, putting for instance X = Range[3] we can note that :
> Orthogonalize[rr, EmpirScal[#1, #2, X] & ]
> does not work, because EmpirScal doesn't works itself :
> In[94]:= EmpirScal[rr[[2]], rr[[3]]]
> Out[94]= x[1]*(x^2)[1] + x[2]*(x^2)[2] + x[3]*(x^2)[3] + x[4]*(x^2)[4]
> + x[5]*(x^2)[5]
>
> It seems that in this case x^k is not interpreted as a function, but as
> the name of some function ; I don't understand why.
>
> So, I tried to work with pure functions, superseding rr by tt :
>
> In[53]:= tt = Table[Function[x, x^i], {i, 0, 2}]
> Out[53]= {Function[x, x^i], Function[x, x^i], Function[x, x^i]}
> This is not what I expected...
>
> Then, I used an old module of mine, converting a polynomiak to a pure
> function :
>
> PolyToPure[pp_] := Module[{x, cc, pol, degre},
> pol = Expand[pp];
> x = First@Variables[N@pol];
> degre = Check[Exponent[pol, x], 0];
> cc = If[degre > 0, CoefficientList[pol, x], {pol}];
> Function[x,
> Apply[Plus,
> Prepend[Table[cc[[i]] x^(i - 1), {i, 2, Length[cc]}], cc[[1]]]]]
> ];
>
> In[122]:= PolyToPure[(1 - x)^3][z]
> Out[122]= 1 - 3 z + 3 z^2 - z^3
>
> Now, we can compute scalar products :
>
> In[7]:= rr = Table[x^i, {i, 0, 2}]
> ss = Map[PolyToPure, rr];
> EmpirScal[ss[[2]], ss[[3]], X]
>
> Out[7]= {1, x, x^2}
> Out[9]= 36
>
>
> But the final result is again unusable :
>
> In[48]:= res = Orthogonalize[ss, EmpirScal[#1, #2, X] & ];
> res[[1]]
>
> Out[49]= (1/Sqrt[3])Function[x$, Plus @@ Prepend[ Table[cc$894[[i]]
> x$^(i - 1), {i, 2, Length[cc$894]}],
> cc$894[[1]]]]
>
> Could you enlighten me?
>
> Thanks,
> Claude Manté
>
>