kill(all)$

fluid : divergence of tensor ,version 1.0,part 1
peter.vlasschaert@gmail.com,16/05/2019

tensor :  scalar  ( tensor of  rank = 0)   , 0 order tensor
                 vector ( tensor of   rank = 1)   ,1 order tensor
                 second order tensor ( tensor of rank = 2)*
==================================================
example : scalar  : temperature ,pressure.    (number)
                   vector  : velocity                                 (3 components)
                  *tensor : stress tensor                       (9 components)
==================================================

components
===========
---------------------------------------------------------------------------------------------------
number : π                                                                                     (1)
vector    : [u,v,w] =[ux,uy,uz]                                                      (3)
tensor    : [σxx,σxy,σxz,σyx,σyy,σyz,σzx,σzy,σzz]             (9)
-----------------------------------------------------------------------------------------------------
x,y,z : index (other use  x:1,y:2,z:3), see simple example below 'm2'
first       index  : direction                            : x
second index : on every  point on area  : xy,xz (vectors inside area)
-----------------------------------------------------------------------------------------------------

show when : T  ( = tensor second order)

div(T)  ?
=======

nabla( m) operator.   (cartesian coordinates)
==============

unit vector
----------------
(ex , ey ,ez)  => (e1,e2,e3)      1:x,2:y,3:z

p1:'diff(m,x)*e1+ 'diff(m,y)*e2 +'diff(m,z)*e3;

p2:grad(m) := ''p1;

simple : example for grad operator.
=============================
m1: x^2*y^5+z

p2a:'grad(m1) = subst([m1=x^2*y^5+z],grad(m1));

1e)  option : diff  , good for this example
-------------------------------------------------------------

p2b1: lhs(p2a) = ev(rhs(p2a),diff);

2e)  option : nouns  , 'general'
---------------------------------------------

p2b2: lhs(p2a) = ev(rhs(p2a),nouns);

rule : tensor , scalar product, commutative (a*b =b*a)
============================================ 

p3:[e1*e1=1,e2*e2 =1,e3*e3 = 1,e1*e2 =0,e1*e3=0,e2*e3=0];

simple : example for tensor scalar product.
m2: e1* T12* e1*e1    = >  show this : T12*e1
m3: e1 *T22 *e2*e2    =>   show this : 0

?  We want control evaluation of a 'expression'
----------------------------------------------------------------------

1e) step :  control evaluation.  (option: general true for 'simp')
===================================================
set option false simp

p4a:simp :false;

2e) step : input of the example ' m2+m3'
=================================
rem : no evaluation

a) m2+m3
---------------

p4b:e1* T12* e1*e1+e1 *T22 *e2*e2 ;

rem : sum , number of terms, 'length'

p4ba:length(p4b);

rem : length 'term'
combination : length,part.

p4ba1:length(part(p4b,1));

p4ba2:length(part(p4b,2));

b) m2,m3
--------------

p4b1:args(p4b);

p4b2:length(p4b1);

p4b11:p4b1[1];
p4b22:p4b1[2];

p4b1a:length(p4b11);

how to use part : 'maxima'
---------------------------------------

p4cc11:part(p4b11,2);

3e) How to write expression on the two first :ei*ej
=========================================
i=1,2,3
i=1,2,3
=>  m2 : (e1*e1)*T12*e1
=>  m3 : (e1*e2) *T22*e2

rulet1(f) :=block( [simp :false],
 l1: [e1^2=1,e2^2=1,e3^2=1,e1*e2=0,e1*e3=0,e2*e3=0],
q:length(f),
t1:part(f,1) *part(f,3),
t2:ev(subst(l1,t1),simp) ,
t3:t2*part(f,2)*part(f,q) ,
t4:ev(t3,simp)          
);

evaluation : m2
=============

p5a:rulet1(p4b11);

evaluation : m3
=============

p5b:rulet1(p4b22);

4e) How to evaluate sum of terms m2+m3
===================================

map(rulet1,p4b1);

rulet1sum(f) := map(rulet1,f);

rulet1sum(p4b1);

rulet2sum(f):= block([], t1:rulet1sum(f));

p6:rulet2sum(p4b1);

p6a:p6[1];

rulet3sum(f):= block([], t1:rulet1sum(f),
ev(sum (t1[i], i, 1, length(t1)),simp));

rulet3sum(p4b1);

rulet4sum(f):= block([], t1:rulet1sum(f),simp:true,
sum (t1[i], i, 1, length(t1)));

rulet4sum(p4b1);

info : about simp
===============

properties(simp);

simp;

simp:false;

bad : example ( simp:false)
=======================

rulet5sum(f):= block([], t1:rulet1sum(f),
sum (t1[i], i, 1, length(t1)));

doesn't   work :  because simp false

apply: subscript must be an integer; found: g441
#0: rulet5sum(f=[e1*T12*e1*e1,e1*T22*e2*e2])
 -- an error. To debug this try: debugmode(true);

good : example ( simp:false)
========================

rulet6sum(f):= block([], t1:rulet1sum(f),
ev(sum (t1[i], i, 1, length(t1)),simp));

rulet6sum(p4b1);