kill(all)$

maxima :apply general BC on a stiff system ( 'finite element'),version 1.0
peter.vlasschaert@gmail.com,20/12/2020

load packages;option
=================

load(lrats)$ ratprint:false $

general boundary condition
======================

powerdisp:true$

Menu :Generate Matrix

1e) Stiffness Matrix ( = SM)

number of rows       :ro
number of columms :co

ro:3;
co:3;

K: genmatrix(k, ro, co);

ne:length(K);

LHS : 

UU1:makelist(u(i),i,1,length(K));
UU2:transpose(UU1);
UU3:K.UU2;

RHS :

FF1:makelist(f(i),i,1,length(K));
FF2:transpose(FF1);

LHS = RHS

M:UU3=FF2;

We get 'ne' :equations
==================

ME:makelist(UU3[i,1]=FF2[i,1],i,1,length(K));

suppose: u(1) is given  : u(i) , i=2,3 unknown
=================

BC : Boundary Condition 

first equation            :boundary equations

p1a:u(1)=α(1);

second equation,.....: general equations ≠ boundary equations

p1b:subst([p1a],ME);

p1b[1]:reverse(p1a);

p1c:p1b;

p1d:makelist(ratcoef(lhs(p1c[i]),rhs(p1a))*rhs(p1a),i,1,length(K));

p1e:p1c-p1d;

p1fa:makelist(p1e[i],i,2,length(K));

p1fb:makelist(u(i),i,2,length(K));

final equations 
===========

LHS :

p1ga:coefmatrix(p1fa,p1fb);

p1gb:transpose(p1fb);

RHS :

rhs(p1fa[1]);

p1gc:makelist(rhs(p1fa[i]),i,1,length(p1ga));
p1gc1:transpose(p1gc);

LHS = RHS

simp;

simp:false;

p1gd:p1ga.p1gb=p1gc1;/* no evaluation */

simp:true;

p1ge:p1gd; /* evaluation */

How construct : New Stifness Matrix ( =NSM)
============================
delete : first row        (when u(1) = α(1))
delete : first column   (when u(1) = α(1))

first step     : orginal SM

p2a:K;

second step : new  NSM

p2b:submatrix(1,p2a,1);

p2a1:row(p2a,1);
p2a2:col(p2a,1);

solve : system:'p1gd' or 'p1ge'
=======================

use :'p1fa', see :above

qa:p1fa;

use :'p1fb'

qb:p1fb;

find :u(i) , i =2,3  i≠1 ( = given,see above)
-------------------------------------------------------

sol:flatten(solve(qa,qb));

sol1:sol[1];
sol2:sol[2];

example : delete (α(2),α(3) is given)

KA: genmatrix(ka, 4, 4);

v1:[2,3];

→  α(2) , dim(KA1) = dim(KA)-1

KA1:submatrix(v1[1],KA,v1[1]);

→  α(3) , dim(KA2) = dim(KA)-1

KA2:submatrix(v1[2]-1,KA1,v1[2]-1);

example : delete (α(2),α(3),α(4) is given)

KB:KA;

v2:[2,3,4];/*list low to high numbers */

→  α(2) , dim(KB) = dim(KB1)-1

KB1:submatrix(v2[1],KB,v2[1]);

→  α(4) , dim(KB2) = dim(KB1)-1

KB2:submatrix(v2[2]-1,KB1,v2[2]-1);

→  α(3) , dim(KB3) = dim(KB2)-1

KB3:submatrix(v2[3]-(length(v2)-1),
    KB2,v2[3]-(length(v2)-1))[1];