How to use multi-part function to find Fourier series coefficients with this program?

Hi I am writing a script to find Fourier series coefficients. Code is given below. In the code, I have function f=x*x defined on interval -L to L. However if I have some function defined, say

f(x) = -1 -L<=x<0
     = 1   0<=x<=L

then, how I can use such function to find Fourier series coefficients? Any ideas?

close all;
clc;

syms x k L n
evalin(symengine,'assume(k,Type::Integer)'); %k is treated as integer
f=x*x; %f is the function to find Fourier series

a = int(f*cos(k*pi*x/L)/L,x,-L,L); %Fourier series coefficient a_n
b =  int(f*sin(k*pi*x/L)/L,x,-L,L); %Fourier series coefficient b_n
afun=@(f,x,k,L) int(f*cos(k*pi*x/L)/L,x,-L,L);
bfun = @(f,x,k,L) int(f*sin(k*pi*x/L)/L,x,-L,L);
fs = vpa(afun(f,x,0,L)/2 + symsum(afun(f,x,k,L)*cos(k*pi*x/L) + ... bfun(f,x,k,L)*sin(k*pi*x/L),k,1,3),5); %Partial sum of Fourier series coefficients

fs=subs(fs,L,1); %making L=1
fs1=symfun(fs,x); %defining symbolic function fs1 for plotting

% disp(a);
% disp(b);
% disp(fs);

figure();
ezplot(fs1,[-5*pi,5*pi]);

Answer 1

[Star Strider]

You are not defining the integral limits correctly.

For example:

syms a b k L w0 x
krnl = symfun(a*cos(-j*k*w0*x) + j*b*sin(-j*k*w0*x), x) % Transform kernel
F1 = int(-1*krnl, x, -L, 0);            % -L -> 0
F2 = int(+1*krnl, x,  0, L);            % 0 -> +L
FT = F1 + F2;                           % Add integrals
FT = subs(FT, L, 1)                     % L = 1
FT = simplify(collect(FT))

produces (in R2014a):

FT =

(4*b*sinh((k*w0)/2)^2)/(k*w0)

MATLAB (correctly) substituted sinh for the sin of an imaginary argument, so this is equivalent to:

FT =

(4*b*sin((j*k*w0)/2)^2)/(k*w0)

The individual ‘b’ coefficients are functions of k. The ‘a’ coefficients are identically zero.

The frequency argument, w0, is the fundamental frequency of the sampled signal. I refer you to the Wikipedia article on Fourier series for the details.

[Christopher Creutzig]

If you would like to get rid of the sinh (to get a more standard representation), try rewrite(FT,'sin').