kill(all)$

maxima : fvm : simple formulation of 4 basic rules , 1D  fvm (W,w,p,e,E) ,CV =[w,p,e],version1.0
peter.vlasschaert@gmail.com,16/11/2022

ratprint:false$ load(lrats)$

Rule 1: Flux across the Control volumes.
----------------------------------------------------

Face common to two adjacent control volumes ,the flux must be same
for two control volumes.

simple : fluid. (no diffusion,steady,Φ(x)=1)
-----------------

see :new:12/11/2022:maxima : :convection -diffusion relation peclet cell number (max :limit),version 1.0

→ conservation law  = continuity equation (Φ(x)=1)

→ start : advection equation

p1:diff(m(x),x) = 0;

p1a:integrate(p1,x,w(i),e(i));

p1a1:part(p1a,1,2,1);

p1a2:p1a+p1a1;

p1aa:m(x)=ρ(x)*u(x)*Φ(x);

transport propeties Φ(x) =1 , conservation law

p1b:subst([Φ(x)=1],p1aa);

p1b1:subst([x=e(i)],p1b);
p1b2:subst([x=w(j)],p1b);

p1b3:[p1b1,p1b2];

flow : pos-direction (u = velocity ,ρ = density)
------------------------------------------------------------

→ first control volume : i=1

p1bb1:subst([i=1],p1a2);

→ second control volume: i=2 ( = adjacent control volume on the rhs)

p1bb2:subst([i=2],p1a2);

result :

p1bbb:lhs(p1bb1)  = rhs(p1bb2);

p1d:subst([i=1,j=2],p1b3);

p1bbbb:subst(p1d,p1bbb);

density uniform in flow between adjacent control volumes : ρ(e(1)) =ρ(w(2))

p1f:subst([ρ(e(1)) =ρ(w(2))],p1bbbb)/ρ(w(2));

result : cell out (1) = cell in (2)

Rule 2 : Positive Coefficients .
--------------------------------------

proof :

see : new:12/11/2022:maxima : :convection -diffusion relation peclet cell number (max :limit),version 1.0

p2a:a[p]*Φ[p] = a[W]*Φ[w] + a[p]*Φ[w] ;

P2b: [a[p]=sign(1) ,a[W]=sign(1) ,a[E]=sign(1) ];

rem :a[p]=2*pos
rem : all coef . "same sign ",(= take in account the flow direction)

Rule 3:Negative linearization of the source term.
----------------------------------------------------------------

p3a:S =S[c]+S[p]*Φ[p];

source : term must be included ,p2a

p2aa:lhs(p2a) = rhs(p2a)+b;

b must not depent on Φ[p]

p3b:lhs(p2a) = rhs(p2a) +S;
p3b1:subst([p3a],p3b);

p3b2:part(rhs(p3a),1);

p3b3:factor(p3b1-p3b2);

a[p] include part linear source term : negative slope of S

p3b4 :reverse(a[p]=ratcoeff(lhs(p3b3),Φ[p]));

p3b5:lratsubst([p3b4],lhs(p3b3))=rhs(p3b3);

→ find :b
use :p2aa,p3b5

p3b6:p2aa-p3b5;

p3b6a:solve(p3b6,b)[1];

--------------------------------------------------------

p3b6aa1:subst([reverse(p3b6a)],p3b5);

p3b6aa2:reverse(p3b4);

---------------------------------------------------------

rem : a[p] =a[W]+a[E] = sum(a[nb],nb,1,2).
---------------------------------------------------------
nc cells  → sum(a[nb],nb,1,nc) ,nc =2 ,1D 

when : S=0
---------------

see :12/11/2022:maxima : :convection -diffusion relation peclet cell number (max :limit),version 1.0
→ p34cb

Rule 4: Sum of Neighbouring Coefficients.=nb
-------------------------------------------------------------

p4:a[p] =sum(a[nb],nb,1,2) -S[p];

p4a:[a[1]=a[W],a[2]=a[E]];
p4a1:subst(p4a,p4);

no source : S=0 -> {S[c] =0,S[p]=0}

p4a2:subst([S[p]=0],p4a1);

see : Rule 3,rem.
---------------------
rem : a[p] with source is different from this without source (=S)

------------------------------------------------------------------------------------------------------------------------------------------