Incorrect dimensions for raising a matrix to a power. Check that the matrix is square and the power is a scalar. To perform elementwise matrix powers, use '.^'.

I don't really get what is wrong with this code:

clear all; close all; clc;

L = 200;

%length of beam in cm

E = 52000;

%modulus of eleasticity in kN/cm^2

I = 32000;

%moment of inertia in cm^4

w = 0.4;

%distributed load in kN/cm

x = linspace(0,200,20);

%discretizing the domain

dx = L/(20-1);

%dx = domain length/(no. of nodes-1) for

i = 1:length(x);

%deflection of the beam

y(i)= (-w*(x(i)^5-(3*L^2*x(i).^3)+(2*L^3*x(i)^2)))/(120*E*I);

%slope analytical

s(i)= (-w*(5*x(i)^4-(9*L^2*x(i)^2)+(4*L^3*x(i))))/(120*E*I);

%bending moment analytical

M(i)= (-w*(20*x(i)^3-(18*L^2*x(i))+(4*L^3)))/120;

%numerical slope

sn(i)= ((-w*((x(i)+dx)^5-(3*L^2*(x(i)+dx)^3)+(2*L^3*(x(i)+dx)^2)))/(120*E*I)- (-w*((x(i)-dx)^5-(3*L^2*(x(i)-dx)^3)+(2*L^3*(x(i)-dx)^2)))/(120*E*I))/(2*dx);

%numerical bending moment

Mn(i) = (E*I/dx^2)*(sn(i)- (2*y(i))+2*(-w*((x(i)-dx)^5-(3*L^2*(x(i)-dx)^3)+(2*L^3*(x(i)-dx)^2)))/(120*E*I));

figure(1)

plot(x,y,x,s,'*',x,sn)

grid on

figure(2)

plot(x,M,x,Mn)

grid on

Errors:

Error using ^ (line 51)

Incorrect dimensions for raising a matrix to a power. Check that the matrix is square and the power is a scalar. To perform elementwise matrix powers, use

'.^'.

Error in untitled (line 16)

y(i)= (-w*(x(i)^5-(3*L^2*x(i).^3)+(2*L^3*x(i)^2)))/(120*E*I);

Answer 1

[Steven Lord]

x^5 is conceptually equivalent to x*x*x*x*x using matrix multiplication.

x.^5 is conceptually equivalent to x.*x.*x.*x.*x using element-wise multiplication.

If x is a 1-by-n vector, x*x is not defined (the inner matrix dimensions don't match) and so neither is x^5.

If x is a 1-by-n vector, both x.*x and x.^5 are defined.

In one place on that line you're using .^ on the vector x and in two places you're using ^.

Answer 2

[VBBV]

The element wise power for a scalar should not affect much on your calculation unless vectors are multiplied or raised to power

 clear all; close all; clc; 
 L = 200; 
 %length of beam in cm 
 E = 52000; 
 %modulus of eleasticity in kN/cm^2 
 I = 32000; 
 %moment of inertia in cm^4 
 w = 0.4; 
 %distributed load in kN/cm 
 x = linspace(0,200,20); 
 %discretizing the domain 
 dx = L/(20-1); 
 %dx = domain length/(no. of nodes-1) 
 for i = 1:length(x)
     %deflection of the beam 
     y(i)= (-w*(x(i)^5-(3*L^2*x(i)^3)+(2*L^3*x(i)^2)))/(120*E*I); % This should not affect 
     %slope analytical 
     s(i)= (-w*(5*x(i)^4-(9*L^2*x(i)^2)+(4*L^3*x(i))))/(120*E*I); 
     %bending moment analytical 
     M(i)= (-w*(20*x(i)^3-(18*L^2*x(i))+(4*L^3)))/120; 
     %numerical slope 
     sn(i)= ((-w*((x(i)+dx)^5-(3*L^2*(x(i)+dx)^3)+(2*L^3*(x(i)+dx)^2)))/(120*E*I)- (-w*((x(i)-dx)^5-(3*L^2*(x(i)-dx)^3)+(2*L^3*(x(i)-dx)^2)))/(120*E*I))/(2*dx); 
     %numerical bending moment 
     Mn(i) = (E*I/dx^2)*(sn(i)- (2*y(i))+2*(-w*((x(i)-dx)^5-(3*L^2*(x(i)-dx)^3)+(2*L^3*(x(i)-dx)^2)))/(120*E*I)); 
 end
 figure(1) 
 plot(x,y,x,s,'*',x,sn) 
 grid on 

 figure(2) 
 plot(x,M,x,Mn) 
 grid on

 y(i)= (-w*(x(i)^5-(3*L^2*x(i).^3)+(2*L^3*x(i)^2)))/(120*E*I);  % whether you use .^ or ^ bcoz x(i) is scalar