kill(all)$

print( " quadratic elements'finite element'")$

TT:T(x):=a1+a2*x+a3*x^2$
T11:T(0)=T1;
T22:T(h/2)=T2;
T33:T(h)=T3;

T1;

eq1:algsys([T11,T22,T33],[a1,a2,a3]);

p1:eq1[1][1];
p2:eq1[1][2];
p3:eq1[1][3];

p4:subst([p1,p2,p3],T(x));

p4;

p5:ratsimp(p4);

p6:expand(p5);

p6;

print("N1(x),N2(x),N3(x) 'shape functions'")$
print("N1 = N1(x) ")$
print("N2 = N2(x) ")$
print("N3 = N3(x) ")$

N1:coeff(p6,T1);
N2:coeff(p6,T2);
N3:coeff(p6,T3);

print("N1'(x),N2'(x),N3'(x) 'first derivatives shape functions'")$
print("N1_D = N1'(x) ")$
print("N2_D = N2'(x) ")$
print("N3_D = N3'(x) ")$

N1_D:diff(N1,x,1);
N2_D:diff(N2,x,1);
N3_D:diff(N3,x,1);