MATLAB: For the force vector part ri and stiffness matrix ki,error is getting displayed while plotting the plots for displacement graph analytically and by approximation methods.What seems to be the probelm,if found ,please explain

entiregalerkin

syms pi p q E A l x
n=3;
k=(n:n);
r=(1:n);
a=(1:n);
N=sym(1:n);
pi=eval(pi);
p=q*x/E*A*l;
for i=1:n  
      syms x
      N(i)=(x.^2-2*x).^i;
      %N(i)=sin((2*i-1)/2*pi*x);
      %N(i)=piecewise(x>=(i-1)/n & x<1/n+(i-1)/n,n*(x-(i-1)/n),x>=i/n & x<=(i+1)/n & x<1,-n*x+i+1,0);
      Y1=(N(i));
      Y2=diff(N(i),x);
      Y3=diff(N(i),x,2);
      x=linspace(0,1);
      plot(x,subs(Y1))
      hold on
      plot(x,subs(Y2))
      hold on
      plot(x,subs(Y3))
      hold off
      legend('N','N''','N''''');
      pause(1)
end
for i=1:n
      r(i)=-int(N(i)*p,x=linspace(0,1))*l ;   %Galerkin 
      %r(i)=-int(diff(N(i),x,2)*p,x=linspace(0,1))/l ;     %Least squares
      for j=1:n
          %k(i,j)=1/l*int(N(i)*diff(N(i),x,2),x=linspace(0,1));      %Galerkin
           %k(i,j)=1/l^3*int(diff (N(j),x,2)*diff(N(i),x,2),x=linspace(0,1));     %Least Squares
           k(i,j)=-1/l*int(diff(N(i),x)*diff(N(j),x),x=linspace(0,1));     %Displacement formulation
      end
end
display(k)
display(r)
a=linsolve(k,r);
U=q*l^2*(x/2-x^3/6)/E*A;
u=o;
for i=1:n
      u=u+a(i)*N(i);
end
u=simplify(u);
plot((U/(q*l^2/E*A)),(u/(q*l^2/E*A),x=linspace(0,1))
X1=0:0.9999;
Z3=(U/(q*l^2/E/A)-u/(q*l^2/E/A));
plot(X1,Z3)

Best Answer

This is as close as you can get:
syms pi p q E A l x
n=3;
k = sym(n:n);
r = sym(1:n);
a=(1:n);
N=sym(1:n);
pi=eval(pi);
p=q*x/E*A*l;
for i=1:n  
        syms x
        N(i)=(x.^2-2*x).^i;
        %N(i)=sin((2*i-1)/2*pi*x);
        %N(i)=piecewise(x>=(i-1)/n & x<1/n+(i-1)/n,n*(x-(i-1)/n),x>=i/n & x<=(i+1)/n & x<1,-n*x+i+1,0);
        Y1=(N(i));
        Y2=diff(N(i),x);
        Y3=diff(N(i),x,2);
        x=linspace(0,1);
        plot(x,subs(Y1))
        hold on
        plot(x,subs(Y2))
        hold on
        plot(x,subs(Y3))
        hold off
        legend('N','N''','N''''');
        pause(1)
end
syms x
for i=1:n
        r(i)=-int(N(i)*p,x,0, 1)*l ;   %Galerkin 
        %r(i)=-int(diff(N(i),x,2)*p,x=linspace(0,1))/l ;     %Least squares
        for j=1:n
            %k(i,j)=1/l*int(N(i)*diff(N(i),x,2),x=linspace(0,1));      %Galerkin
             %k(i,j)=1/l^3*int(diff (N(j),x,2)*diff(N(i),x,2),x=linspace(0,1));     %Least Squares
             k(i,j)=-1/l*int(diff(N(i),x)*diff(N(j),x),x, 0, 1);     %Displacement formulation
        end
end
display(k)
display(r)
a=linsolve(k,r.');
U=q*l^2*(x/2-x^3/6)/E*A;
u = sym(0);
for i=1:n
        u=u+a(i)*N(i);
end
u=simplify(u);
X = linspace(0,1);
expr1 = subs((U/(q*l^2/E*A)), x, X);
expr2 = subs((u/(q*l^2/E*A)), x, X);
plot(expr1, expr2)
X1=0:0.9999;
Z3 = subs((U/(q*l^2/E/A)-u/(q*l^2/E/A)), x, X1);
plot(X1,Z3)
The two plot() at the end will fail. The first one will fail because expr2 involves the unresolved symbolic variable l (lower-case L). The second one will fail because Z3 involves the unresolved symbolic variables A and l (lower-case L)
Related Question