clear all;
clc;

% 1 trial function
syms c1 x N1 u
N1 = x*(1-x);
u = c1*N1;
du_dx = diff(u);
f = du_dx^2-u^2+2*u*x^2;
I = 0.5*int(f,0,1);
del_I = diff(I, c1);
c1_1 = solve(del_I);

% 2 trial functions
syms c1 c2 x N1 N2 u
N1 = x*(1-x);
N2 = x^2*(1-x);
u = c1*N1+c2*N2;
du_dx = diff(u);
f = du_dx^2-u^2+2*u*x^2;
I = 0.5*int(f,0,1);
del_I1 = diff(I, c1);
del_I2 = diff(I, c2);
[c1_2 c2_2] = solve(del_I1,del_I2);

% 3 trial functions
syms c1 c2 c3 x N1 N2 N3 u
N1 = x*(1-x);
N2 = x^2*(1-x);
N3 = x^3*(1-x);
u = c1*N1+c2*N2+c3*N3;
du_dx = diff(u);
f = du_dx^2-u^2+2*u*x^2;
I = 0.5*int(f,0,1);
del_I1 = diff(I, c1);
del_I2 = diff(I, c2);
del_I3 = diff(I, c3);
[c1_3 c2_3 c3_3] = solve(del_I1,del_I2, del_I3);