kill(all)$

maxima : enter : read ,solve reaction equation   , version 1.4
peter.vlasschaert@gmail.com,16/11/2019

rem : yellow band must be go away
===========================

definition : lcmlistresult(npb,npa,vv)
------------------------------------------------
vv   =>  [x1=%r1,x2=%r1/2,y1=%r1/2,y2=%r1]
npb => products                 (b = before)
npa => reaction products  (a = after    )

lcmlistresult(npb,npa,vv):=block([],
    s:npb+npa,
    v1:makelist(rhs(vv[i]),i,1,s),
    v2:makelist(denom(v1[i]),i,1,s),
    v3:lcm(v2),
    return(v3))$    

make  varlist :  {x1,x2,...,y1,y2,...}
==========

varlistba(xx,npb,yy,npa):=block([],
tb:makelist(concat(xx,i),i,1,npb),
ta:makelist(concat(yy,i),i,1,npa),
tt:[tb,ta],
tt:flatten(tt),
return(tt)
)$

definition
=======

inprodre(ll,xx,jj):=block([lza:[] ,ll1:[] ],
for i:1 thru length(ll)/2 do(
lza:endcons([ll[2*i-1],ll[2*i]],lza)),
nf1:length(lza),
for i:1 thru nf1+1 do(
ll1:endcons(0,ll1) ),
for i:1 thru nf1 do(
    ll1[i]:lza[i]
),
ll1[nf1+1]:concat(xx,jj),
 for i:1 thru nf1 do(
  if  ll1[i] [2] = 1 then ll1[i] [2] :z
),return(ll1)
)$   

rem :  xx = x  =>  reaction
                   y =>  reaction products
ll : [Na,1,O,1,H,1]
-----------------------------------------------------------------
optimized : concat inside , three argument input
new :definition

formula(m):=block([],
  mx:delete(last(m) ,m),
    mx:flatten(mx),
    mx:delete(z,mx),
    mp:concat(mx[1],mx[2]),
        for i:3 thru length(mx) do(
    mp:concat(mp,mx[i])
      ) ,
return(mp)
)$

remove numbers from the list
======================

removenuml( ii,mrl ) := block([ls1:[]],
   for j:1 thru length(mq2[ii]) /2  do(
    ls1:endcons(mrl[ii][2*j-1],ls1),
    mq3[ii] : ls1          
    ),
    ls1:[] ,
    return(mq3[ii])    
)$

find elements in the : reaction
=======================

findel(nw,mrel):=block([ls2:[]],
for i:1 thru nw do(
  ls2:endcons(mrel[i],ls2)
    ),
ls2:flatten(ls2),
el1:setify(ls2) ,
el2: listify(el1),     
return(el2))$

reaction :
=======
x2*H2SO4+x1*NaOH=y1*Na2SO4+y2*H2O
(before reaction       )    (after reaction           )
-----------------------------------------------------------------
npb = 2  ( number of products before : reaction)
npa = 2  ( number of products after    : reaction)
-----------------------------------------------------------------

print(" number of products : integer")$

products : before reaction
====================

npb : read( "number products (before reaction)   :'shift' & 'enter'       = ")$

 => enter product  ( before reaction )
products : {NaOH,H2SO4}

for i:1 thru npb do(
qbf[i]:read( "number products (before reaction)   :'shift' & 'enter'       = ")
        )$

rem : 'shift & enter' for every 'product'

convert :computation

for i:1 thru  npb do(
    qb[i]:inprodre(qbf[i],x,i)
)$

for i:1 thru npb do(
    qqb[i] :delete(last(qb[i]) ,qb[i])
)$

products :after reaction
====================

npa : read( "number products (after reaction)   :'shift' & 'enter'       = ")$

 => enter product  ( before reaction )
products : {Na2SO4,H2O}

for i:1 thru npa do(
qaf[i]:read( "number products (before reaction)   :'shift' & 'enter'       = ")
        )$

convert :computation

for i:1 thru  npa do(
    qa[i]:inprodre(qaf[i],y,i)
)$

for i:1 thru npa do(
    qqa[i] :delete(last(qa[i]) ,qa[i])
)$

elements in the reaction :
===================
before same = after same

for i:1 thru npb do(
mq2[i] : flatten(qb[i]),
mq2[i] :delete(last(mq2[i]),mq2[i])   
)$

mq2[2];

remove numbers from the list
-----------------------------------------

for i:1 thru npb do(
    removenuml( i,mq2 )
)$

element list ( = before the reaction)
===========================
same : unique.

el:findel(npb,mq3);

1e) before  reaction : display
---------------------------------------

br:makelist([formula((qb[j])),last(qb[j])],j,1,npb);
br1:[sum(br[j][2]*br[j][1],j,1,npb)]$

2e)after  reaction : display
------------------------------------

ar:makelist([formula(qa[j]),last(qa[j])],j,1,npa);
ar1:[sum(ar[j][2]*ar[j][1],j,1,npa)]$

combination (1e,2e)
================

g1:br1[1] =ar1[1];

to find : coef {x1,x2,.....} & {y1,y2,....} :reaction : g1
======================================

lhs 
==

mb[i]:=makelist([qqb[i][j][1],qqb[i][j][2]*concat(x,i)],j,1,length(qqb[i]))$
lb:[]$
for i:1  thru npb do(
lb:endcons(mb[i],lb)
    )$
lb:flatten(lb)$

lhs :  equation  (list)
---------------------------

for j:1 thru length(el) do(
    sbb:[],         
for i:1  thru length(lb) do(
    if el[j] = lb[i] then  
        (
        sbb:endcons(lb[i+1],sbb)
        )                   
)  ,
    mfb[j]:sbb    
)$

lhs => equation, [H,Na,O,S],sum
=========================

for i:1 thru  length(el) do(
    mfb[i]:apply("+",mfb[i])
    )$

rhs
===

ma[i]:=makelist([qqa[i][j][1],qqa[i][j][2]*concat(y,i)],j,1,length(qqa[i]))$
la:[]$
for i:1  thru npa do(
la:endcons(ma[i],la)
    )$
la:flatten(la)$

rhs :  equation  (list)
---------------------------

for j:1 thru length(el) do(
    saa:[],
for i:1  thru length(la) do(
    if el[j] = la[i] then  
        (
        saa:endcons(la[i+1],saa)
        )                   
)  ,
    mfa[j]:saa
    
)$

rhs => equation, [H,Na,O,S],sum
=====================

for i:1 thru  length(el) do(
    mfa[i]:apply("+",mfa[i])
    )$

final : element balancce: [H,Na,O,S]
--------------------------------------------------
replace :  z = 1

for i:1 thru length(el) do(
    fer[i]:subst([z=1],mfb[i]=mfa[i])
)$

output : equations : [H,Na,O,S]
==============

for i:1 thru length(el) do(
  print(  fer[i] )
)$

find coef : {x1,x2,...} & {y1,y2,...}
========================

list : [x1,x2,....,y1,y2,]
----------------------------

use varlistba: 'see above '
====================

qq:varlistba(x,npb,y,npa);

qq1:length(qq)$

list : equations  fer[i], i=1,2,..,el
-------------------------------------------
=>   [H,Na,O,S]

ss:makelist(ev(fer[i]),i,1,length(el))$

c1: coefmatrix(ss,tt);

c2 : matrix
========

c2:genmatrix(s, qq1+1, qq1+1)$

1e)  coef matrix 
-----------------------

for i:2  thru qq1+1 
   do(
    for j:2 thru qq1+1 
    do(
        c2[i,j]:c1[i-1,j-1]
    )
)$

2e)  col    : [x1,x2,y1,y2]   
        rows : [H,Na,O,S]
----------------------------------

for i:2 thru  qq1+1 do (
   c2[i,1]:el[i-1]
    )$
for j:2 thru qq1+1 do(
     c2[1,j]:qq[j-1] 
    )$
c2[1,1] :x$

final matrix :
=========

output : matrix : [H,Na,O,S]
=====================

c2:c2;

find :coefmatrix
============
use : submatrix

c3:submatrix(1,c2,1)$

solve =
==================================
x2*H2SO4+x1*NaOH=y1*Na2SO4+y2*H2O
==================================

ratprint :false$
frsm:solve(ss,qq)$

frsm:flatten(frsm);

website:

definition : 'see above'

appendix : LCM ( = least common multiplier)
peter.vlasschaert@gmail.com,14/11/2019

pq:lcmlistresult(npb,npa,frsm);

frm:[%r1=1,%r2=1,%r3=1];

frsm1:subst(frm,frsm);

subst:{x1,x2,....,y1,y2,...} in the reaction g1
================================

g2:expand(pq*subst(frsm1,g1));