function [bv,bv2] = bevaluate_unimorph(N,Np,X,h) h3 = 1/(h^3); bval = zeros(Np,N+3); bval2 = zeros(Np,N+3); one = 1; for ell=1:Np x = X(ell,1); for kk = 1:N+3 k = kk-2; ind1 = k-2; ind2 = k-1; ind3 = k; ind4 = k+1; ind5 = k+2; if (ind1*h<=x)&(x