[math-fun] "new" kind of fourier series
From: Neil Bickford <techie314@gmail.com> Unfortunately (if I've interpreted the above correctly), it looks like expanding the D-bit version over individual bits of T just results in a product of a whole bunch of phase-shifted sine waves, scaled by some amount*.
This can't represent arbitrary functions (or even as many as any individual Fourier transform can), because it has zeroes all over the place.**
*Since C_{A,0...}*Sin[i*x1*a1]*(rest of the sum...) + C_{A,1...}*Cos[i*x1*a1]*(rest of the sum...) = Sqrt[c1^2+c2^2]*Sin[i*x1*a1 + ArcTan[c2/c1]]*(rest of the sum...).
** Unless the series expresses a constant, in which case it's, um, well, not extraordinarily interesting.
--Neil Bickford
--WDS responds: The fact that sin(x) and cos(x) have many zeros is no obstacle to SUM_k a[k] * sin(x*k) + b[k] * cos(x*k) representing any 2*pi-periodic nice-enough function. This sum is, in fact, a standard 1-dimensional fourier series, which DOES represent any 2*pi-periodic nice-enough function. Therefore I fail to understand your "This can't represent arbitrary functions (or even as many as any individual Fourier transform can), because it has zeroes all over the place... unless a constant." If you expand exp(i*x) = cos(x)+i*sin(x) and exp(i*x+i*y) = [cos(x)+i*sin(x)]*[cos(y)+i*sin(y)] = cos(x)*cos(y)+i*sin(x)*cos(y)+...-sin(x)*sin(y) etc, then the usual fourier series in the form USUAL FOURIER = SUM_A C_A * exp(i*A.X) summed over integer vectors A -- which will represent any nice-enough complex-valued function of D variables that is periodic in a D-dimensional box with sidelengths 2*pi -- is converted via the expansion process to a form like my second kind of fourier series. In fact, it is the same as my second kind of fourier series, EXCEPT that the coefficients of the new series are determined (via the expansion process) from a much small number of coeffs of the old series, so the "extra information" in all those extra coefficients is illusory. But if we now permit the extra coeffs to be altered independently so that that extra information really is there, then we get a wider class of series -- which was my point. And it still seems to me to be a good point, I'm not seeing anything wrong about it. This "new" kind of fourier series seems likely to me to be more useful than the old plane-wave kind for many applications. Which is more useful? Probably depends on the application.
Ah - it looks like I both forgot about the D-dimensional sum, and employed a fallacious implicit induction, which doesn't work for dimensions above 1, making my statements completely incorrect. Whoops! However, it actually looks like any real-valued series in Warren's form can be converted to an equivalent real-valued usual Fourier series, if you include both positive and negative frequencies (technically, positive and negative for every coordinate except the first - more on that later) with the same number of coefficients as the new kind of Fourier series. One quick way of seeing this, experimentally, is just to use any local computer algebra system (in this case, Mathematica, sorry) and expand the internal summand for some finite number of dimensions, in this case 3: (We've switched Warren's notation such that a bit of 0 = COS and a bit of 1 = SIN (mostly for mathematical convenience, so that we don't have to include a whole bunch of factors of -i later)) Expanding over the 3-bit word (we're expressing the sin/cos choice using a phase shift) In[54]:= Sum[a[b1,b2,b3](Cos[x-Pi/2 b1]*Cos[y-Pi/2 b2]*Cos[z-Pi/2 b3]), {b3,0,1},{b2,0,1},{b1,0,1}] Out[54]= a[0,0,0] Cos[x] Cos[y] Cos[z]+ a[1,0,0] Cos[y] Cos[z] Sin[x]+ a[0,1,0] Cos[x] Cos[z] Sin[y]+ a[1,1,0] Cos[z] Sin[x] Sin[y]+ a[0,0,1] Cos[x] Cos[y] Sin[z]+ a[1,0,1] Cos[y] Sin[x] Sin[z]+ a[0,1,1] Cos[x] Sin[y] Sin[z]+ a[1,1,1] Sin[x] Sin[y] Sin[z] Multiplying out and collecting terms, In[55]:= Collect[TrigReduce[%],Sin[var_]|Cos[var_]] Out[55]= (1/4 a[0,0,0]-1/4 a[0,1,1]+1/4 a[1,0,1]+1/4 a[1,1,0]) Cos[x-y-z]+ (1/4 a[0,0,0]+1/4 a[0,1,1]+1/4 a[1,0,1]-1/4 a[1,1,0]) Cos[x+y-z]+ (1/4 a[0,0,0]+1/4 a[0,1,1]-1/4 a[1,0,1]+1/4 a[1,1,0]) Cos[x-y+z]+ (1/4 a[0,0,0]-1/4 a[0,1,1]-1/4 a[1,0,1]-1/4 a[1,1,0]) Cos[x+y+z]+ (-(1/4) a[0,0,1]-1/4 a[0,1,0]+1/4 a[1,0,0]-1/4 a[1,1,1]) Sin[x-y-z]+ (-(1/4) a[0,0,1]+1/4 a[0,1,0]+1/4 a[1,0,0]+1/4 a[1,1,1]) Sin[x+y-z]+ (1/4 a[0,0,1]-1/4 a[0,1,0]+1/4 a[1,0,0]+1/4 a[1,1,1]) Sin[x-y+z]+ (1/4 a[0,0,1]+1/4 a[0,1,0]+1/4 a[1,0,0]-1/4 a[1,1,1]) Sin[x+y+z] which looks like a few terms of a Fourier series! But you can go further than that, and actually write out a full derivation. (At this point, the thread gets really long. Feel free to skip this if you wish - I've also included testable code samples at the bottom that you can also use.) Disclaimer: This proof is not exactly short, probably not neat, and certainly not elegant. It's mostly just the most direct way of proving the statement possible - in this case, by expanding out terms and collecting them. I'd deeply appreciate a different proof of equivalence, possibly by taking the Fourier transform of the new series and noting that the only peaks occur at Dirac deltas bounded by a box, observing that the function is periodic and low-frequency, using a singular value decomposition-like approach, or even something competely different. One final note: I've changed the notation a bit, not just to confuse everybody, but simply because Sum[ c[f1,f2,...]*e^(i*x1*f1)*e^(i*x2*f2)..., {0,0,...}<={f1,f2,...}<={a1,a2,...} ] seems more easily readable than SUM_A B_A * exp(i*X*A), even though the two mean the same thing. On to the proof! Expressing each Cos[fi*xi-Pi/2*bi] as a sum of two exponentials, Sum[c[f1,f2,f3,...,b1,b2,b3,...]* (e^(i*(f1*x1 - Pi/2*b1))+e^(-i*(f1*x1 - Pi/2*b1)))/2* (e^(i*(f2*x2 - Pi/2*b2))+e^(-i*(f2*x2 - Pi/2*b2)))/2* (e^(i*(f3*x3 - Pi/2*b3))+e^(-i*(f3*x3 - Pi/2*b3)))/2*..., 0<={f1,f2,f3...}<=A, {b1,b2,b3...bd} all in {0,1}] = Sum[a[f...,b...](exp(i*(f1*x1 - Pi/2*b1))+exp(i*(-f1*x1 + Pi/2*b1)))/2* (exp(i*(f1*x1 - Pi/2*b1))+exp(i*(-f1*x1 + Pi/2*b1)))/2 (exp(i*(f1*x1 - Pi/2*b1))+exp(i*(-f1*x1 + Pi/2*b1)))/2..., 0<={f...}<=A, {b...} all in {0,1}] Distributing, multiplying out, and expressing the resulting binary explosion using D new indices, each either -1 or 1, to denote signs: Sum[a[f...,b...](1/2^D)* e^(i*(s1*f1*x1 + s2*f2*x2 + s3*f3*x3... - Pi/2(s1*b1 + s2*b2 + s3*b3...))) 0<={f...}<=A, {b...} all in {0,1}, {s...} all in {-1,1}] Expanding over s1 and combining {1,s2,s3...} with {-1,-s2,-s3...}, Sum[a[f...,b...](2/2^D)*(1/2)*( e^( i*(f1*x1 + s2*f2*x2 + s3*f3*x3... - Pi/2(b1 + s2*b2 + s3*b3...)))+ e^(-i*(f1*x1 + s2*f2*x2 + s3*f3*x3... - Pi/2(b1 + s2*b2 + s3*b3...)))) 0<={f...}<=A, {b...} all in {0,1}, {s2...} all in {-1,1}] = Sum[a[f...,b...](2/2^D)*(( Cos[f1*x1 + s2*f2*x2 + s3*f3*x3... - Pi/2(b1 + s2*b2 + s3*b3...)], 0<={f...}<=A, {b...} all in {0,1}, {s2...} all in {-1,1}] We now only need to combine the coefficients of terms with the same frequency- space coordinates (that is, terms with a factor of Cos[f1*x1+s2*f2*x2...]), which is as simple as saying Sum[Sum[(2/2^D)*a[f...,b...]* Cos[f1*x1 + s2*f2*x2 + s3*f3*x3... - Pi/2(b1 + s2*b2 + s3*b3...)], {b...} all in {0,1}], 0<={f...}<=A, {s2...} all in {-1,1}]. Which, if you want a more obvious Fourier series, is equivalent to Sum[(2/2^D)*(-1)^Floor[(b1 + s2*b2 + s3*b3...)/2] * (Sum[a[f...,b...], {b...} all in {0,1} such that Total[b] is even] * Cos[f1*x1 + s2*f2*x2 + s3*f3*x3 ...] + Sum[a[f...,b...], {b...} all in {0,1} such that Total[b] is odd] * Sin[f1*x1 + s2*f2*x2 + s3*f3*x3 ...]), 0<={f...}<=A, {s2...} all in {-1,1}]. And now, for code! We're assuming the user has provided the array in an order such that coeffs[[f1+1, f2+1, f3+1, ..., fD+1, t+1]] is C_{F,T}. (Mathematica uses an array offset of 1 - sorry) So, for instance, let's suppose we have a two-dimensional function, with each frequency ranging between 0 and 2. Then a random coefficient array is: In[59]:= RandomReal[{-1,1},{2+1,2+1,2^2}] Out[59]= {{{-0.556713,-0.935173,0.198351,-0.560439}, {-0.673702,0.447996,0.0991223,0.586851}, {-0.10958,-0.0223377,0.904763,0.737812}}, {{-0.115752,-0.248104,-0.858183,0.654516}, {-0.57084,-0.300345,0.904632,0.294495}, {-0.296196,0.97974,0.350523,0.956529}}, {{-0.0282067,-0.70006,0.686257,0.557952}, {0.104394,-0.0305375,0.592332,0.728952}, {-0.00136298,0.479633,-0.822651,0.695625}}} in which case the coefficient cooresponding to the {1,2,{1,0}} harmonic is: In[66]:= Extract[%59,Append[{1,2}+1,FromDigits[{1,0},2]+1]] Out[66]= 0.350523 where our part index is In[67]:= Append[{1,2}+1,FromDigits[{1,0},2]+1] Out[67]= {2,3,3} Here's a function that constructs a (cos and sin-swapped) sum of the new type, given the coefficients. This code, unfortunately, would probably qualify for the Unintentionally Obfuscated Mathematica competition: newfourier[coeffs_]:=Block[{numDimensions,boxSize,maxT,iteratives}, (numDimensions=Depth[coeffs]-2; boxSize=Dimensions[coeffs]; (* gets the values of D *) maxT=Last[boxSize]; (* gets the number of possible values for T *) boxSize=Drop[boxSize,-1]; iteratives=Range/@boxSize; (* prepares the n-dimensional iteration*) Sum[Extract[coeffs,Append[f,t]]* (*c[f...,t...] *) Times@@Table[ (*form a length-D product...*) If[BitGet[t-1,numDimensions-1-bi]==0,Cos,Sin][(f[[bi+1]]-1)*x[bi]], {bi,0,numDimensions-1}], {f,Tuples[iteratives]},{t,1,maxT}] ) ] Finally, here's full code to convert Fourier series of Warren's new type to the old type of Fourier series, maintaining the number of coefficients: fourierconvert[coeffs_]:=Block[ {numDimensions,boxSize,maxT,iteratives, signs,vars}, ( (*More or less the same initialization as newfourier, except for the initialization of an array of signs and variables.*) numDimensions=Depth[coeffs]-2; If[numDimensions<=1, (Print["Already in Fourier form!"]; Return[newfourier[coeffs]])]; boxSize=Dimensions[coeffs]; maxT=Last[boxSize]; boxSize=Drop[boxSize,-1]; iteratives=Range/@boxSize; signs=ConstantArray[{-1,1},numDimensions-1]; vars=Table[x[i-1],{i,1,numDimensions}]; (2/2^numDimensions)* Sum[Extract[coeffs,Append[f,t]]* Cos[Prepend[s,1].((f-1)*vars)- Pi/2(Prepend[s,1].IntegerDigits[t-1,2,numDimensions])], {f,Tuples[iteratives]},{t,1,maxT},{s,Tuples[signs]}] )] and just a quick kinda-sorta-verification that they're identical: In[123]:= {newfourier[%59],fourierconvert[%59]}/.{x[0]->x,x[1]->y} Out[123]= { -0.556713 -0.115752 Cos[x] -0.0282067 Cos[2 x] -0.673702 Cos[y] -0.57084 Cos[x] Cos[y] +0.104394 Cos[2 x] Cos[y] -0.10958 Cos[2 y] -0.296196 Cos[x] Cos[2 y] -0.00136298 Cos[2 x] Cos[2 y] -0.858183 Sin[x] +0.904632 Cos[y] Sin[x] +0.350523 Cos[2 y] Sin[x] +0.686257 Sin[2 x] +0.592332 Cos[y] Sin[2 x] -0.822651 Cos[2 y] Sin[2 x] +0.447996 Sin[y] -0.300345 Cos[x] Sin[y] -0.0305375 Cos[2 x] Sin[y] +0.294495 Sin[x] Sin[y] +0.728952 Sin[2 x] Sin[y] -0.0223377 Sin[2 y] +0.97974 Cos[x] Sin[2 y] +0.479633 Cos[2 x] Sin[2 y] +0.956529 Sin[x] Sin[2 y] +0.695625 Sin[2 x] Sin[2 y] , 1/2 (-1.11343-0.231504 Cos[x] -0.0564135 Cos[2 x] +0.660334 Cos[x-2 y] +0.694262 Cos[2 x-2 y] -0.276345 Cos[x-y] +0.833346 Cos[2 x-y] -1.3474 Cos[y] -0.21916 Cos[2 y] -0.865335 Cos[x+y] -0.624559 Cos[2 x+y] -1.25273 Cos[x+2 y] -0.696988 Cos[2 x+2 y] -1.71637 Sin[x] +1.37251 Sin[2 x] -0.629218 Sin[x-2 y] -1.30228 Sin[2 x-2 y] +1.20498 Sin[x-y] +0.62287 Sin[2 x-y] +0.895992 Sin[y] -0.0446754 Sin[2 y] +0.604287 Sin[x+y] +0.561795 Sin[2 x+y] +1.33026 Sin[x+2 y] -0.343018 Sin[2 x+2 y])} In[128]:= newfourier[%59]-fourierconvert[%59]//Simplify Out[128]= 5.55112*10^-17 Cos[x[0]-2 x[1]] +5.55112*10^-17 Cos[2 x[0]+x[1]] +1.11022*10^-16 Sin[x[0]-x[1]] +1.11022*10^-16 Sin[2 (x[0]-x[1])] -5.55112*10^-17 Sin[2 x[0]-x[1]] -5.55112*10^-17 Sin[x[0]+x[1]] -5.55112*10^-17 Sin[2 x[0]+x[1]] (*removing miniscule powers of 10 created by floating point noise*) In[130]:= Chop[%128] Out[130]= 0 and there you go! Alternatively, see http://neilbickford.com/assets/fourier-plottest.png . Hope this works (this message has now exceeded 2^8 lines, sorry), --Neil Bickford
participants (2)
-
Neil Bickford -
Warren D Smith