kill(all)$

load : packages
============

load(lrats)$

pde : separation variables heat equation.(problem :given :ic and bc),version 1.0
peter.vlasschaert@gmail.com,22/07/2019

see: website
=> pde : separation variables for principal of the heat equation.(not general solution),version 1.0
?  Why : -β^2 for separation constant.
?  General principal

given : problem. (unsteady probem : du/dt ≠0)
============

pde :  D*diff(u,x,2) = diff(u,x,t), u = u(x,t)
==============================
ic :  u(x,0) = f(x)            'initial condition'
bc : 1e)  u(0,t)  = 0   'boundary condition'
        2e) u(l,t)    = 0
==============================

assume(D>0);

p1:D*diff(u(x,t),x,2) = diff(u(x,t),t);

p2:u(x,t)  = X(x)*T(t);

p3:'diff(X(x),x,2)/X(x)='diff(T(t),t,1)/(D*T(t));

two : odes

assume(β>0);

p3a:(lhs(p3) =- β^2)*X(x);
p3b:(rhs(p3) = -β^2)*D*T(t);

solve : two odes with desolve.
=======================

X(x) : ?
======

p4a:desolve(p3a,X(x));

p4a1:lhs(p4a) =C1*sin(x*β)+C2*cos(x*β);

T(t) : ?
======

p4b:desolve(p3b,T(t));

a) boundary condition  ' into X(x).
=========================

1e) bc :  x=0  u(0,t)  => X(0) =0
--------------------------------------

p5a:subst([x=0],p4a1);

given : problem 'above'

p5a1:rhs(p5a) =0;

result : X(x)
=========

p5a11:subst([p5a1],p4a1);

1e) bc :  x=l  u(l,t)  => X(l) =0
--------------------------------------
use : 'p5a1'

p5a2:subst([x=l],p5a1);

given : problem 'above'

p5a3:rhs(p5a2) = 0;

assume : C1 ≠0;

assume(C1>0,l>0);

p5a4:solve(p5a3,β);

p5a5:solve(sin(β*l)=0,β);

wxplot2d([sin(x)], [x,-2*%pi,2*%pi],
 [plot_format, xmaxima])$

become zero 'sin(x)' :  x=0, x=+/-3.14, x= +/- 2*3.14.
see : figure  x=0 ,  n*π = x   => here problem :  x = β*l
=> found   'p5a4,p5a5'  n.π = 0  => β =0, l always l ≠ 0

p5a6:args(p5a3);

p5a7:β*l = n*π;

solve : β  from p5a7
---------------------------

p5a8:solve(p5a7,β);

final result : X(x)
============
use : p5a11

p5a9:subst(p5a8,p5a11);

=>  we find depend, 'n' for X(x)

b) initial condition  ' into T(t).
=========================
@ t= 0 

final result : T(t)
===========
use : p4b

p5b:lratsubst(p5a8,p4b);

p5b1:subst([T(0)=C3],p5b);

find  u(x,t)  =>
===========
use : 'p2'

p6:[p5b1,p5a9];
p6a:lratsubst(p6,p2);

p6a: we found product of two constant = C4

p6b:lratsubst([C1*C3=C4],p6a);

 ic :  t=0  u(x,0)  => u(x,0) = f(x)
-----------------------------------------

p7a:[t=0,u(x,0)=f(x)];

p7b:lratsubst(p7a,p6b);

How to find : C4
============
use : 'p7b';

p8a:lratsubst([C4=1],rhs(p7b));

principal : orthogonality relations
========================

1e) p7b multiplied p8a
==================

p8b:p7b*p8a;

2e) integration x:0 to l
=================

assume(n>0);

p8c:integrate(p8b,x,0,l);

before : fourier series n = 1,2, ...∞

p8d:subst([n=1],p8c);

build test function: Q(n)
==================
use :p8c

p8e:Q(n):=''p8c;

rem:make list with  Q(n)
==================

p8f:makelist(Q(n), n, 1, 10);

first element list 'p8f'
----------------------------

p8g:p8f[1];

solve : C4
========
use : 'p8f' and 'p8c'

general solution : C4
================
use :first element list

p9:lhs(p8c) = rhs(p8g);

p9a:solve(p9,C4);

find : 'use result' , 'p9a'  = > one term fourier series
============================-==========
use : 'p6b'

p10:subst(p9a,rhs(p6b));

General solution : u(x,t)
==================

u(x,t) :We make now general solution for n=1,2,...,∞,  ( fourier series)
----------------------------------------------------------------------
use :'p2'

;p10a:lhs(p2)=sum(p10,n,1,∞);

orthogonality   (relations)
==========

1e) multiplied 'p7b' by cos((%pi*n*x)/l)
=============================

p11:p7b*cos((%pi*n*x)/l);

2e) integration x:0 to l
=================

p11a:integrate(p11,x,0,l);

test n=1 ?
========

p11a1:subst([n=1],p11a);

test n=10 ?
=========

p11a2:subst([n=10],p11a);

build function : do it your self.

general rule: orthogonality relations
=========
integrate ;0 to l:{ f1,f2,f3,f4}
===============================
f1:sin(f(n)) * sin(f(m))  , n=m  => f1≠0,p9
f2:cos(f(n)) *cos(f(m)) , n=m => f2 ≠0
f3:cos(f(n)) *sin(f(m)) , n=m => f3 =0,p11a
f4:sin(f(n)) *cos(f(m)) , n=m => f4=0,p11a
================================
see above : 'p11a,p9' 

rem : representation , general function H(x)
================================
see :theory fourier series
=> piecewise continue and other conditions, see internet.

p12:H(x) = a[(0)]/2  + sum(a[(n)]*cos((n*π*x)/l) +b[n]*sin((n*π*x)/l),n,1,∞);

general rules : from above for 'orthogonality
==================================

a[n]: values

n=0
===

p13a:a[0]=1/(l/2)*integrate(H(x),x,0,l);

n>0  , n=1,2,....,∞
-----------------------

p13b:a[n]=1/(l/2) *integrate(H(x)*cos((%pi*n*x)/l),x,0,l);

p13c:b[n]=1/(l/2) *integrate(H(x)*sin((%pi*n*x)/l),x,0,l);

rem:n=0,into 'p13b'

p14:subst([n=0],p13b);

check : 'is' command
-----------------------------
p14  ?  p13a

before : check , use 'ratsimp,expand,ev,...'
best : ev( expression,nouns).

is(rhs(p14)=rhs(p13a));