How to solve a system of 8 equations having integrals!

4 views (last 30 days)
d=1.18; E=27000000; L=121.53;
I=(pi*d^4)/64;
x1=16.99; x2=34.16; x3=51.33; x4=68.49; x5=85.66; x6=102.83; x7=115.61; x8=120.12;
F1=49; F2=49; F3=49; F4=49; F5=49; F6=49; F7=49; F8=11;
K1=100; K2=100; K3=100; K4=100; K5=100; K6=100; K7=100; K8=100;
These are my moment equation, that are further used in the deflection (d1, d2, ...., d8) calculations.
M1 = @(x)(F1*x-K1*d1*x);
M2 = @(x)(F1*x+F2*(x-x1)-K1*d1*x-K2*d2*(x-x1));
M3 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2));
M4 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3));
M5 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)+F5*(x-x4)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3)-K5*d5*(x-x4));
M6 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)+F5*(x-x4)+F6*(x-x5)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3)-K5*d5*(x-x4)-K6*d6*(x-x5));
M7 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)+F5*(x-x4)+F6*(x-x5)+F7*(x-x6)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3)-K5*d5*(x-x4)-K6*d6*(x-x5)-K7*d7*(x-x6));
M8 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)+F5*(x-x4)+F6*(x-x5)+F7*(x-x6)+F8*(x-x7)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3)-K5*d5*(x-x4)-K6*d6*(x-x5)-K7*d7*(x-x6)-K8*d8*(x-x7));
Below are my eight equations with defined boundaries for integrals. And it contains eight variables i.e. d1, d2, ...., d8 as well. I am fine with both numerical or analytical solution.
eq1 = d1 ==(1/(E*I))*(integral(@(x)(M1(x).*x),0,x1)+integral(@(x)(M2(x).*x),x1,x2)+integral(@(x)(M3(x).*x),x2,x3)+integral(@(x)(M4(x).*x),x3,x4)+integral(@(x)(M5(x).*x),x4,x5)+integral(@(x)(M6(x).*x),x5,x6)+integral(@(x)(M7(x).*x),x6,x7)+integral(@(x)(M8(x).*x),x7,x8));
eq2 = d2 ==(1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*(x-x1)),x1,x2)+integral(@(x)(M3(x).*(x-x1)),x2,x3)+integral(@(x)(M4(x).*(x-x1)),x3,x4)+integral(@(x)(M5(x).*(x-x1)),x4,x5)+integral(@(x)(M6(x).*(x-x1)),x5,x6)+integral(@(x)(M7(x).*(x-x1)),x6,x7)+integral(@(x)(M8(x).*(x-x1)),x7,x8));
eq3 = d3 ==(1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*(x-x2)),x2,x3)+integral(@(x)(M4(x).*(x-x2)),x3,x4)+integral(@(x)(M5(x).*(x-x2)),x4,x5)+integral(@(x)(M6(x).*(x-x2)),x5,x6)+integral(@(x)(M7(x).*(x-x2)),x6,x7)+integral(@(x)(M8(x).*(x-x2)),x7,x8));
eq4 = d4 ==(1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*(x-x3)),x3,x4)+integral(@(x)(M5(x).*(x-x3)),x4,x5)+integral(@(x)(M6(x).*(x-x3)),x5,x6)+integral(@(x)(M7(x).*(x-x3)),x6,x7)+integral(@(x)(M8(x).*(x-x2)),x7,x8));
eq5 = d5 ==(1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*0),x3,x4)+integral(@(x)(M5(x).*(x-x4)),x4,x5)+integral(@(x)(M6(x).*(x-x4)),x5,x6)+integral(@(x)(M7(x).*(x-x4)),x6,x7)+integral(@(x)(M8(x).*(x-x4)),x7,x8));
eq6 = d6 ==(1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*0),x3,x4)+integral(@(x)(M5(x).*0),x4,x5)+integral(@(x)(M6(x).*(x-x5)),x5,x6)+integral(@(x)(M7(x).*(x-x5)),x6,x7)+integral(@(x)(M8(x).*(x-x5)),x7,x8));
eq7 = d7 ==(1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*0),x3,x4)+integral(@(x)(M5(x).*0),x4,x5)+integral(@(x)(M6(x).*0),x5,x6)+integral(@(x)(M7(x).*(x-x6)),x6,x7)+integral(@(x)(M8(x).*(x-x6)),x7,x8));
eq8 = d8 ==(1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*0),x3,x4)+integral(@(x)(M5(x).*0),x4,x5)+integral(@(x)(M6(x).*0),x5,x6)+integral(@(x)(M7(x).*0),x6,x7)+integral(@(x)(M8(x).*(x-x7)),x7,x8));

Accepted Answer

Torsten
Torsten on 29 Nov 2022
d0 = ones(1,8);
d = fsolve(@fun,d0)
Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient.
d = 1×8
0.5346 0.5133 0.4823 0.4784 0.3431 0.2181 0.0778 0.0067
norm(fun(d))
ans = 1.4844e-15
function res = fun(D)
d=1.18;
E=27000000;
L=121.53;
I=(pi*d^4)/64;
x1=16.99; x2=34.16; x3=51.33; x4=68.49; x5=85.66; x6=102.83; x7=115.61; x8=120.12;
F1=49; F2=49; F3=49; F4=49; F5=49; F6=49; F7=49; F8=11;
K1=100; K2=100; K3=100; K4=100; K5=100; K6=100; K7=100; K8=100;
d1 = D(1);
d2 = D(2);
d3 = D(3);
d4 = D(4);
d5 = D(5);
d6 = D(6);
d7 = D(7);
d8 = D(8);
M1 = @(x)(F1*x-K1*d1*x);
M2 = @(x)(F1*x+F2*(x-x1)-K1*d1*x-K2*d2*(x-x1));
M3 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2));
M4 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3));
M5 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)+F5*(x-x4)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3)-K5*d5*(x-x4));
M6 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)+F5*(x-x4)+F6*(x-x5)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3)-K5*d5*(x-x4)-K6*d6*(x-x5));
M7 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)+F5*(x-x4)+F6*(x-x5)+F7*(x-x6)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3)-K5*d5*(x-x4)-K6*d6*(x-x5)-K7*d7*(x-x6));
M8 = @(x)(F1*x+F2*(x-x1)+F3*(x-x2)+F4*(x-x3)+F5*(x-x4)+F6*(x-x5)+F7*(x-x6)+F8*(x-x7)-K1*d1*x-K2*d2*(x-x1)-K3*d3*(x-x2)-K4*d4*(x-x3)-K5*d5*(x-x4)-K6*d6*(x-x5)-K7*d7*(x-x6)-K8*d8*(x-x7));
res(1) = d1 -((1/(E*I))*(integral(@(x)(M1(x).*x),0,x1)+integral(@(x)(M2(x).*x),x1,x2)+integral(@(x)(M3(x).*x),x2,x3)+integral(@(x)(M4(x).*x),x3,x4)+integral(@(x)(M5(x).*x),x4,x5)+integral(@(x)(M6(x).*x),x5,x6)+integral(@(x)(M7(x).*x),x6,x7)+integral(@(x)(M8(x).*x),x7,x8)));
res(2) = d2 -((1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*(x-x1)),x1,x2)+integral(@(x)(M3(x).*(x-x1)),x2,x3)+integral(@(x)(M4(x).*(x-x1)),x3,x4)+integral(@(x)(M5(x).*(x-x1)),x4,x5)+integral(@(x)(M6(x).*(x-x1)),x5,x6)+integral(@(x)(M7(x).*(x-x1)),x6,x7)+integral(@(x)(M8(x).*(x-x1)),x7,x8)));
res(3) = d3 -((1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*(x-x2)),x2,x3)+integral(@(x)(M4(x).*(x-x2)),x3,x4)+integral(@(x)(M5(x).*(x-x2)),x4,x5)+integral(@(x)(M6(x).*(x-x2)),x5,x6)+integral(@(x)(M7(x).*(x-x2)),x6,x7)+integral(@(x)(M8(x).*(x-x2)),x7,x8)));
res(4) = d4 -((1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*(x-x3)),x3,x4)+integral(@(x)(M5(x).*(x-x3)),x4,x5)+integral(@(x)(M6(x).*(x-x3)),x5,x6)+integral(@(x)(M7(x).*(x-x3)),x6,x7)+integral(@(x)(M8(x).*(x-x2)),x7,x8)));
res(5) = d5-((1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*0),x3,x4)+integral(@(x)(M5(x).*(x-x4)),x4,x5)+integral(@(x)(M6(x).*(x-x4)),x5,x6)+integral(@(x)(M7(x).*(x-x4)),x6,x7)+integral(@(x)(M8(x).*(x-x4)),x7,x8)));
res(6) = d6 -((1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*0),x3,x4)+integral(@(x)(M5(x).*0),x4,x5)+integral(@(x)(M6(x).*(x-x5)),x5,x6)+integral(@(x)(M7(x).*(x-x5)),x6,x7)+integral(@(x)(M8(x).*(x-x5)),x7,x8)));
res(7) = d7 -((1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*0),x3,x4)+integral(@(x)(M5(x).*0),x4,x5)+integral(@(x)(M6(x).*0),x5,x6)+integral(@(x)(M7(x).*(x-x6)),x6,x7)+integral(@(x)(M8(x).*(x-x6)),x7,x8)));
res(8) = d8 -((1/(E*I))*(integral(@(x)(M1(x).*0),0,x1)+integral(@(x)(M2(x).*0),x1,x2)+integral(@(x)(M3(x).*0),x2,x3)+integral(@(x)(M4(x).*0),x3,x4)+integral(@(x)(M5(x).*0),x4,x5)+integral(@(x)(M6(x).*0),x5,x6)+integral(@(x)(M7(x).*0),x6,x7)+integral(@(x)(M8(x).*(x-x7)),x7,x8)));
end
  3 Comments

Sign in to comment.

More Answers (1)

Walter Roberson
Walter Roberson on 29 Nov 2022
Do not use == with integral() . integral() is for numeric integration, and if you use numeric integration returning a numeric value then the right side of the == would be numeric, and the == would be for asking whether d1 is exactly the same value, bit-for-bit identical (except for different nans being considered the same, and negative 0 being treated the same as regular 0)
If you want to process this system numerically using fsolve() then you should change the NAME == EXPRESSION into NAME - (EXPRESSION) and then you would have the function return [eq1, eq2, eq3, eq4, eq5, eq6, eq7, eq8]
You could alternately switch to using the Symbolic Toolbox and use int() or vpaintegral() instead of integral, and try to solve() or vpasolve() the vector of equations (which can be written with == in that case.)
  4 Comments

Sign in to comment.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by