I have a matrix equation system and I need help to solve it
1 次查看(过去 30 天)
显示 更早的评论
mahdi
2023-5-24
So I have these 3, 7*1 variable matrices in those system of equations and C{1,1} and C{1,2} are 7*7 matrices and they are constant, also S_1 to S_15 are 1*7 constant matrices. Problem is, it gives me U=0, W=0 and PHI=0 and it doesnt calculate the variables in those U, W and PHI matrices.
syms U [7 1]
U;
syms W [7 1]
W;
syms PHI [7 1]
PHI;
syms U W PHI
eq1 = (S_1*C{1,2}*U)+(S_2*C{1,1}*U)-(S_3*U)+(S_4*C{1,2}*W)+(S_5*C{1,1}*W) == 0;
eq2 = (S_6*C{1,2}*PHI)+(S_7*C{1,1}*PHI)+(S_8*PHI)-(S_9*C{1,1}*W) == 0;
eq3 = (S_10*C{1,2}*W)+(S_11*C{1,1}*W)+(S_12*C{1,1}*PHI)+(S_13*PHI)+(S_14*C{1,1}*U)...
-(S_15*U) == 0;
[U, W, PHI] = solve([eq1, eq2, eq3],[U, W, PHI])
44 个评论
Dyuman Joshi
2023-5-24
You are overwriting the variables U, W and PHI. Remove the line where you define them again.
Note that this doesn't guarantee obtaining a solution as there are 21 unknowns and 3 equations.
syms U [7 1]
syms W [7 1]
syms PHI [7 1]
U
W
PHI
%Remove this line
%This line overwrites U, W and PHI from 7x1 syms variable
%to 1x1 syms variable
syms U W PHI
U
W
PHI
mahdi
2023-5-24
Well, thanks for that point you mentioned. Actually for example S_1 = ((2*A_11).*X1.^2) and same for the other S and also X1 is a 1*7 vector and those W_1 and U_1 and PHI_1 must be the answer for first element of that X1 vector:
X1 =
0 0.0335 0.1250 0.2500 0.3750 0.4665 0.5000
A_11 is just scalar.
There are also some boundary conditions which I have not yet implemented them yet.
mahdi
2023-5-25
I update my question here: C{1,1} and C{1,2} are 7*7 matrices with numeric elements.
U_g1 = 0; U_g2 = 0; U_g3 = 0;
B1 = [-0.0493 -0.0490 -0.0455 -0.0348 -0.0188 -0.0051 0];
W0 = B1; Nu = 0.285;
A_11 = 64; A_55 = 37; D_11 = 2;
syms U [7 1]
syms W [7 1]
syms PHI [7 1]
U;
W;
PHI;
for z=[0,0.03349,0.1250,0.2500,0.3750,0.4665,0.5000]
for v=[-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0]
for q=[0,-0.0016,-0.0060,-0.0111,-0.0145,-0.0156,-0.0157]
eq1 = (((2*A_11).*z.^2)*C{1,2}*U)+(((2*A_11).*z)*C{1,1}*U)-((2*A_11)*U)+...
(((A_11.*v).*z)*C{1,2}*W)+((((A_11-(Nu*A_11)).*q).*z)*C{1,1}*W) == 0;
eq2 = ((D_11.*z.^2)*C{1,2}*PHI)+((D_11.*z)*C{1,1}*PHI)+(((D_11+A_55).*z)*PHI)...
-((A_55.*z.^2)*C{1,1}*W) == 0;
eq3 = ((((2*A_55).*z)+((2*A_11.*U_g2).*z)+((A_11.*v.^2).*z)+...
(2*Nu*A_11.*U_g3))*C{1,2}*W)+(((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z)+(A_11.*q.^2)-(((A_11.*W0).*q).*z))*C{1,1}*W)+...
(((2*A_55).*z)*C{1,1}*PHI)+((2*A_55)*PHI)+(((2*Nu*A_11)-((2*A_11.*W0).*z))*C{1,1}*U)...
-((2*Nu*A_11.*W0)*U) == 0;
[U, W, PHI] = solve([eq1, eq2, eq3],[U, W, PHI])
end
end
end
Undefined variable 'C'.
The error message I get: (U, W and PHI are 7*1 matrices, why does it say inconsistent output with only 3 variables?)
Error using sym/solve (line 279)
Inconsistent output with 3 variables for input argument with 21 variables.
Error in Untitled (line 24)
[U, W, PHI] = solve([eq1, eq2, eq3],[U, W, PHI])
Dyuman Joshi
2023-5-25
"Note that this doesn't guarantee obtaining a solution as there are 21 unknowns and 3 equations."
I said this in an above comment.
How can you solve for 21 unknowns when you only have 3 equations?
mahdi
2023-5-25
Thank you for your suggestion, I did make Y = zeros(7,1) and substitute it with 0's of right hand side, but it still gives me the same error as before
Torsten
2023-5-25
编辑:Torsten
2023-5-25
Works for me - although the solution looks boring:
U_g1 = 0; U_g2 = 0; U_g3 = 0;
B1 = [-0.0493 -0.0490 -0.0455 -0.0348 -0.0188 -0.0051 0];
W0 = B1; Nu = 0.285;
A_11 = 64; A_55 = 37; D_11 = 2;
syms U [7 1]
syms W [7 1]
syms PHI [7 1]
for z=[0.03349]
for v=[-0.0493]
for q=[-0.0016]
C11 = rand(7);
C12 = rand(7);
eq1 = (((2*A_11).*z.^2)*C11*U)+(((2*A_11).*z)*C11*U)-((2*A_11)*U)+...
(((A_11.*v).*z)*C12*W)+((((A_11-(Nu*A_11)).*q).*z)*C11*W) == zeros(7,1);
eq2 = ((D_11.*z.^2)*C12*PHI)+((D_11.*z)*C11*PHI)+(((D_11+A_55).*z)*PHI)...
-((A_55.*z.^2)*C11*W) == zeros(7,1);
eq3 = ((((2*A_55).*z)+((2*A_11.*U_g2).*z)+((A_11.*v.^2).*z)+...
(2*Nu*A_11.*U_g3))*C12*W)+(((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z)+(A_11.*q.^2)-(((A_11.*W0).*q).*z))*C11*W)+...
(((2*A_55).*z)*C11*PHI)+((2*A_55)*PHI)+(((2*Nu*A_11)-((2*A_11.*W0).*z))*C11*U)...
-((2*Nu*A_11.*W0)*U) == zeros(7,1);
UWPHI = solve([eq1, eq2, eq3],[U, W, PHI])
end
end
end
UWPHI = struct with fields:
U1: 0
U2: 0
U3: 0
U4: 0
U5: 0
U6: 0
U7: 0
W1: 0
W2: 0
W3: 0
W4: 0
W5: 0
W6: 0
W7: 0
PHI1: 0
PHI2: 0
PHI3: 0
PHI4: 0
PHI5: 0
PHI6: 0
PHI7: 0
mahdi
2023-5-25
编辑:mahdi
2023-5-25
When I run this with complete z, v and q vectors, it takes some time and it shows this:
UWPHI =
struct with fields:
U1: [1×1 sym]
U2: [1×1 sym]
U3: [1×1 sym]
U4: [1×1 sym]
U5: [1×1 sym]
U6: [1×1 sym]
U7: [1×1 sym]
W1: [1×1 sym]
W2: [1×1 sym]
W3: [1×1 sym]
W4: [1×1 sym]
W5: [1×1 sym]
W6: [1×1 sym]
W7: [1×1 sym]
PHI1: [1×1 sym]
PHI2: [1×1 sym]
PHI3: [1×1 sym]
PHI4: [1×1 sym]
PHI5: [1×1 sym]
PHI6: [1×1 sym]
PHI7: [1×1 sym]
mahdi
2023-5-25
So you are saying that it is better to save it in 3d cell array right?(How can I?)
And the solutions are all 0.(Why is it though...?)
Once again thank you for help and time.
Torsten
2023-5-25
编辑:Torsten
2023-5-25
The solutions are all 0 because you have a system of the form A*[U,V,PHI] = 0.
If you change the right-hand sides to something different from 0, you will get more interesting solutions.
U_g1 = 0; U_g2 = 0; U_g3 = 0;
B1 = [-0.0493 -0.0490 -0.0455 -0.0348 -0.0188 -0.0051 0];
W0 = B1; Nu = 0.285;
A_11 = 64; A_55 = 37; D_11 = 2;
syms U[7 1]
syms W[7 1]
syms PHI[7 1]
Z = [0,0.03349,0.1250,0.2500,0.3750,0.4665,0.5000];
V = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0];
Q = [0,-0.0016,-0.0060,-0.0111,-0.0145,-0.0156,-0.0157];
for iz = 1:numel(Z)
z = Z(iz);
for iv = 1:numel(V)
v = V(iv);
for iq = 1:numel(Q)
q = Q(iq);
eq1 = (((2*A_11).*z.^2)*C{1,2}*U)+(((2*A_11).*z)*C{1,1}*U)-((2*A_11)*U)+...
(((A_11.*v).*z)*C{1,2}*W)+((((A_11-(Nu*A_11)).*q).*z)*C{1,1}*W) == zeros(7,1);
eq2 = ((D_11.*z.^2)*C{1,2}*PHI)+((D_11.*z)*C{1,1}*PHI)+(((D_11+A_55).*z)*PHI)...
-((A_55.*z.^2)*C{1,1}*W) == zeros(7,1);
eq3 = ((((2*A_55).*z)+((2*A_11.*U_g2).*z)+((A_11.*v.^2).*z)+...
(2*Nu*A_11.*U_g3))*C{1,2}*W)+(((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z)+(A_11.*q.^2)-(((A_11.*W0).*q).*z))*C{1,1}*W)+...
(((2*A_55).*z)*C{1,1}*PHI)+((2*A_55)*PHI)+(((2*Nu*A_11)-((2*A_11.*W0).*z))*C{1,1}*U)...
-((2*Nu*A_11.*W0)*U) == zeros(7,1);
UWPHI{iz,iv,iq} = solve([eq1, eq2, eq3],[U, W, PHI])
end
end
end
Chunru
2023-5-26
==> There are 3*7 equations because each eq_i is made up of 7 equations. I think it should suffice to substitute the 0 of the right-hand side by a 7x1 vector of 0's and keep the output variables general (as "UWPHI", e.g.).
eq1 = (S_1*C{1,2}*U)+(S_2*C{1,1}*U)-(S_3*U)+(S_4*C{1,2}*W)+(S_5*C{1,1}*W) == 0;
This is a single equation instead of 7 eqations. Since S_1 is 1x7, C{1,2} is 7x7 and U is 7x1, S_1*C{1,2} is a constant 1x7 vector. Therefore S_1*C{1,2}*U)+(S_2*C{1,1}*U is a linear combination of 7 elements of U. Similarly for the other terms in the left side of eq1. Therefore this is a single equation of U and W (with 14 unknows).
mahdi
2023-5-26
编辑:mahdi
2023-5-26
Thanks for help guys.
More update: This is so much harder that I thought
U_g1 = 0; U_g2 = 0; U_g3 = 0;
syms U [7,1]
Q1 = sum(C{1,1}*U);
V1 = sum(C{1,2}*U);
syms W [7,1]
Q2 = sum(C{1,1}*W);
V2 = sum(C{1,2}*W);
syms PHI [7,1]
Q3 = sum(C{1,1}*PHI);
V3 = sum(C{1,2}*PHI);
z1 = [0,0.03349,0.1250,0.2500,0.3750,0.4665,0.5000];
z2 = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0];
z3 = [0,-0.0016,-0.0060,-0.0111,-0.0145,-0.0156,-0.0157];
B1 = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0];
eq1 = (((2*A_11).*z1.^2)*V1)+(((2*A_11).*z1)*Q1)-((2*A_11)*U)+...
(((A_11.*z2).*z1)*V2)+((((A_11-(Nu*A_11)).*z3).*z1)*Q2) == 0;
eq2 = ((D_11.*z1.^2)*V3)+((D_11.*z1)*Q3)+(((D_11+A_55).*z1)*PHI)...
-((A_55.*z1.^2)*Q1) == 0;
eq3 = (((((2*A_55).*z1)+((2*A_11.*U_g2).*z1)+((A_11.*z2.^2).*z1)+...
(2*Nu*A_11.*U_g3))*V2)+(((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z1)+(A_11.*z3.^2)+((2*A_11).*U_g3)+...
((2*A_11).*z1))*Q1)+(((2*A_55).*z1)*Q3)+((2*A_55)*PHI)) == ...
(((((2*Nu*A_11).*U_g2).*z1)+((A_11.*z3.^2).*z1)...
+((2*Nu*A_11).*U_g3)).*B1);
I must actually write a code that puts the first elements of z1, z2, z3 and B1 into the equations and also gets the U1, W1 and PHI1 from the Q1, V1 and Q2, V2 and Q3, V3 and puts them into the equations and then solve the equations for U1, W1 and PHI1 (3 variables and 3 equations I think) and goes for the second elements of z1, z2, z3 and B1 and also solve for U2, W2 and PHI2 and so on.
Torsten
2023-5-26
I don't understand what you mean by "get the U1, W1 and PHI1 from the Q1, V1 and Q2, V2 and Q3, V3".
You have defined 99 equations for 21 unknowns.
Please explain what you are trying to do.
C{1,1} = rand(7);
C{1,2} = rand(7);
syms U [7,1]
Q1 = sum(C{1,1}*U);
V1 = sum(C{1,2}*U);
syms W [7,1]
Q2 = sum(C{1,1}*W);
V2 = sum(C{1,2}*W);
syms PHI [7,1]
Q3 = sum(C{1,1}*PHI);
V3 = sum(C{1,2}*PHI);
z1 = [0,0.03349,0.1250,0.2500,0.3750,0.4665,0.5000];
z2 = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0];
z3 = [0,-0.0016,-0.0060,-0.0111,-0.0145,-0.0156,-0.0157];
B1 = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0];
U_g1 = 0; U_g2 = 0; U_g3 = 0;
W0 = B1; Nu = 0.285;
A_11 = 64; A_55 = 37; D_11 = 2;
eq1 = (((2*A_11).*z1.^2)*V1)+(((2*A_11).*z1)*Q1)-((2*A_11)*U)+...
(((A_11.*z2).*z1)*V2)+((((A_11-(Nu*A_11)).*z3).*z1)*Q2) == 0
eq1 =
eq2 = ((D_11.*z1.^2)*V3)+((D_11.*z1)*Q3)+(((D_11+A_55).*z1)*PHI)...
-((A_55.*z1.^2)*Q1) == 0
eq2 =
eq3 = (((((2*A_55).*z1)+((2*A_11.*U_g2).*z1)+((A_11.*z2.^2).*z1)+...
(2*Nu*A_11.*U_g3))*V2)+(((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z1)+(A_11.*z3.^2)+((2*A_11).*U_g3)+...
((2*A_11).*z1))*Q1)+(((2*A_55).*z1)*Q3)+((2*A_55)*PHI)) == ...
(((((2*Nu*A_11).*U_g2).*z1)+((A_11.*z3.^2).*z1)...
+((2*Nu*A_11).*U_g3)).*B1)
eq3 =
size(eq1)
ans = 1×2
7 7
size(eq2)
ans = 1×2
1 7
size(eq3)
ans = 1×2
7 7
Torsten
2023-5-26
Maybe you mean
C{1,1} = rand(7);
C{1,2} = rand(7);
syms U [7,1]
Q1 = sum(C{1,1}*U);
V1 = sum(C{1,2}*U);
syms W [7,1]
Q2 = sum(C{1,1}*W);
V2 = sum(C{1,2}*W);
syms PHI [7,1]
Q3 = sum(C{1,1}*PHI);
V3 = sum(C{1,2}*PHI);
z1 = [0,0.03349,0.1250,0.2500,0.3750,0.4665,0.5000].';
z2 = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0].';
z3 = [0,-0.0016,-0.0060,-0.0111,-0.0145,-0.0156,-0.0157].';
B1 = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0].';
U_g1 = 0; U_g2 = 0; U_g3 = 0;
W0 = B1; Nu = 0.285;
A_11 = 64; A_55 = 37; D_11 = 2;
eq1 = (((2*A_11).*z1.^2)*V1)+(((2*A_11).*z1)*Q1)-((2*A_11)*U)+...
(((A_11.*z2).*z1)*V2)+((((A_11-(Nu*A_11)).*z3).*z1)*Q2) == 0
eq1 =
eq2 = ((D_11.*z1.^2)*V3)+((D_11.*z1)*Q3)+(((D_11+A_55).*z1).*PHI)...
-((A_55.*z1.^2)*Q1) == 0
eq2 =
eq3 = (((((2*A_55).*z1)+((2*A_11.*U_g2).*z1)+((A_11.*z2.^2).*z1)+...
(2*Nu*A_11.*U_g3))*V2)+(((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z1)+(A_11.*z3.^2)+((2*A_11).*U_g3)+...
((2*A_11).*z1))*Q1)+(((2*A_55).*z1)*Q3)+((2*A_55)*PHI)) == ...
(((((2*Nu*A_11).*U_g2).*z1)+((A_11.*z3.^2).*z1)...
+((2*Nu*A_11).*U_g3)).*B1)
eq3 =
[A,b] = equationsToMatrix([eq1,eq2,eq3])
A =
b =
rank(A)
ans = 16
rank([A,b])
ans = 17
But as you can see, your linear system A*X = b of 21 equations in 21 unknowns cannot be solved since the right-hand side vector b is not in the column space of the matrix A.
mahdi
2023-5-26
编辑:mahdi
2023-5-26
Sorry for not being able to explain well, what I mean is this:
% To solve equations for U2 W2 PHI2
syms U2 W2 PHI2
U_g1 = 0; U_g2 = 0; U_g3 = 0;
eq1 = (((2*A_11)*0.03349^2)*U2)+(((2*A_11)*0.03349)*U2)-((2*A_11)*U2)+...
(((A_11*-0.0490)*0.03349)*W2)+((((A_11-(Nu*A_11))*-0.0016)*0.03349)*W2) == 0;
eq2 = ((D_11*0.03349^2)*PHI2)+((D_11*0.03349)*PHI2)+(((D_11+A_55)*0.03349)*PHI2)...
-((A_55*0.03349^2)*W2) == 0;
eq3 = (((((2*A_55)*0.03349)+((2*A_11*U_g2)*0.03349)+((A_11*-0.0490^2)*0.03349)+...
(2*Nu*A_11*U_g3))*W2)+(((2*A_55)+(2*A_11*U_g2)...
+((2*A_11*U_g1)*0.03349)+(A_11*-0.0016^2)+((2*A_11)*U_g3)+...
((2*A_11)*0.03349))*W2)+(((2*A_55)*0.03349)*PHI2)+((2*A_55)*PHI2)) == ...
(((((2*Nu*A_11)*U_g2)*0.03349)+((A_11*-0.0016^2)*0.03349)...
+((2*Nu*A_11)*U_g3))*-0.0490);
[U2, W2, PHI2] = vpasolve([eq1, eq2, eq3],[U2, W2, PHI2])
% To solve equations for U5 W5 PHI5
syms U5 W5 PHI5
U_g1 = 0; U_g2 = 0; U_g3 = 0;
eq1 = (((2*A_11)*0.3750^2)*U5)+(((2*A_11)*0.3750)*U5)-((2*A_11)*U5)+...
(((A_11*-0.0188)*0.3750)*W5)+((((A_11-(Nu*A_11))*-0.0145)*0.3750)*W5) == 0;
eq2 = ((D_11*0.3750^2)*PHI5)+((D_11*0.3750)*PHI5)+(((D_11+A_55)*0.3750)*PHI5)...
-((A_55*0.3750^2)*W5) == 0;
eq3 = (((((2*A_55)*0.3750)+((2*A_11*U_g2)*0.3750)+((A_11*-0.0188^2)*0.3750)+...
(2*Nu*A_11*U_g3))*W5)+(((2*A_55)+(2*A_11*U_g2)...
+((2*A_11*U_g1)*0.3750)+(A_11*-0.0145^2)+((2*A_11)*U_g3)+...
((2*A_11)*0.3750))*W5)+(((2*A_55)*0.3750)*PHI5)+((2*A_55)*PHI5)) == ...
(((((2*Nu*A_11)*U_g2)*0.3750)+((A_11*-0.0145^2)*0.3750)...
+((2*Nu*A_11)*U_g3))*-0.0188);
[U5, W5, PHI5] = vpasolve([eq1, eq2, eq3],[U5, W5, PHI5])
U2 =
-0.0000000000028296697285378618465161974969249
W2 =
0.0000000032533718692757014678911154291619
PHI2 =
0.000000000097278799663482288883395933977357
U5 =
-0.0000000058690397669296133169393354511255
W5 =
0.00000051981438465991820998979600647868
PHI5 =
0.00000017096290394124225864833987217252
Lets take the vector z1 as 7 nodes on a line z1 = [0,0.03349,0.1250,0.2500,0.3750,0.4665,0.5000], now there is a function w_0 = w_0*cos(pi*r/(2*R)) and we got O = diff(w_0,r); P = diff(w_0,r,2); then we got z3 = eval(subs(O,r,z1));
z2 = eval(subs(P,r,z1)); then we got B1 = eval(subs(P,r,z1)); which is the same as z2 but we need that cause if you could remember B1=W0 isnt gonna change at all in the last equation but the z2 and z3 are gonna be initial guesses and I must replace them with the new values of W's the equations give and also the same replacing goes for those U_g1 = 0; U_g2 = 0; U_g3 = 0; too but I have not made it to this replacing part yet, cause Im gonna need some error function, epsilon and loop of i=1:100 for that.
So here I am now I can just write these equations 7 times for each set of U's and W's and PHI's or is there gonna be a faster code that helps me out. Cause this time its 7 nodes, next time its gonna be for example 9,11,13,15 nodes and I will lose my hands to write these equations 15 times for each set of U's and W's and PHI's.
I hope I explained this well this time. (Really sry if I got you a headache)
mahdi
2023-5-26
编辑:mahdi
2023-5-26
Another explaination for this part for example:
syms U [7,1]
Q1 = sum(C{1,1}*U)
V1 = sum(C{1,2}*U)
Q1 =
(46802780229323767*U2)/1125899906842624 - (200410183417987123*U1)/3377699720527872 - (9007199254740995*U3)/1688849860263936 + (5*U4)/2251799813685248 + (18014398509482023*U5)/3377699720527872 - (93605560458647603*U6)/2251799813685248 + (801640733671948765*U7)/13510798882111488 + 8*3^(1/2)*U2 - 8*3^(1/2)*U6
V1 =
(626000348204499175*U1)/281474976710656 - (246196779629587011*U2)/70368744177664 + (144115188075855919*U3)/70368744177664 - (4640*U4)/3 + (72057594037927949*U5)/35184372088832 - (984787118518348255*U6)/281474976710656 + (156500087051124827*U7)/70368744177664
So solving the equations for U1, W1 and PHI1 Im gonna need:
- (200410183417987123*U1)/3377699720527872 from Q1 and (626000348204499175*U1)/281474976710656 from V1 and and z1 = [0], z2 = [-0.0493], z3 = [0] and B1 = [-0.0493] and - (200410183417987123*W1)/3377699720527872 from Q2 and (626000348204499175*W1)/281474976710656 from V2 and - (200410183417987123*PHI1)/3377699720527872 from Q3 and (626000348204499175*PHI1)/281474976710656 from V3 inside the equtions of this comment: https://www.mathworks.com/matlabcentral/answers/1972419-i-have-a-matrix-equation-system-and-i-need-help-to-solve-it#comment_2760489
Another notice is that you must change the remaining U, W, PHI in the equations of the comment to U1, W1 and PHI1.
Torsten
2023-5-26
编辑:Torsten
2023-5-26
C{1,1} = rand(7);
C{1,2} = rand(7);
syms U [7,1]
Q1 = (sum(C{1,1},1).*U.').';
V1 = (sum(C{1,2},1).*U.').';
syms W [7,1]
Q2 = (sum(C{1,1},1).*W.').';
V2 = (sum(C{1,2},1).*W.').';
syms PHI [7,1]
Q3 = (sum(C{1,1},1).*PHI.').';
V3 = (sum(C{1,2},1).*PHI.').';
z1 = [0,0.03349,0.1250,0.2500,0.3750,0.4665,0.5000].';
z2 = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0].';
z3 = [0,-0.0016,-0.0060,-0.0111,-0.0145,-0.0156,-0.0157].';
B1 = [-0.0493,-0.0490,-0.0455,-0.0348,-0.0188,-0.0051,0].';
U_g1 = 0; U_g2 = 0; U_g3 = 0;
W0 = B1; Nu = 0.285;
A_11 = 64; A_55 = 37; D_11 = 2;
eq1 = (((2*A_11).*z1.^2).*V1)+(((2*A_11).*z1).*Q1)-((2*A_11)*U)+...
(((A_11.*z2).*z1).*V2)+((((A_11-(Nu*A_11)).*z3).*z1).*Q2) == 0;
eq2 = ((D_11.*z1.^2).*V3)+((D_11.*z1).*Q3)+(((D_11+A_55).*z1).*PHI)...
-((A_55.*z1.^2).*Q1) == 0;
eq3 = (((((2*A_55).*z1)+((2*A_11.*U_g2).*z1)+((A_11.*z2.^2).*z1)+...
(2*Nu*A_11.*U_g3)).*V2)+(((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z1)+(A_11.*z3.^2)+((2*A_11).*U_g3)+...
((2*A_11).*z1)).*Q1)+(((2*A_55).*z1).*Q3)+((2*A_55).*PHI)) == ...
(((((2*Nu*A_11).*U_g2).*z1)+((A_11.*z3.^2).*z1)...
+((2*Nu*A_11).*U_g3)).*B1);
vars = [U;W;PHI];
[A,b] = equationsToMatrix([eq1,eq2,eq3],vars);
rank(A)
ans = 20
rank([A,b])
ans = 20
sol = A\b
Warning: Solution is not unique because the system is rank-deficient.
sol =
vpa(sol(2)) %U2
ans =
0.00000000011733770039154233650122075968525
vpa(sol(2+7)) %W2
ans =
vpa(sol(2+14)) %PHI2
ans =
0.000000000011399565221492271690978509835807
vpa(sol(5)) %U5
ans =
vpa(sol(5+7)) %W5
ans =
vpa(sol(5+14)) %PHI5
ans =
mahdi
2023-5-26
Is this the complete form that gives the solution from U_1 to U_7 and W_1 to W_7 and PHI_1 to PHI_7? If it is then thank you so much.
Walter Roberson
2023-5-26
z2 = eval(subs(P,r,z1));
Never eval() a symbolic expression. Mathworks does not document any meaning for eval() of a symbolic expression, and you will get errors and unexpected answers sometimes.
Torsten
2023-5-26
编辑:Torsten
2023-5-26
Is this the complete form that gives the solution from U_1 to U_7 and W_1 to W_7 and PHI_1 to PHI_7?
Should be so. But as you can see, rank(A) is only 20, not 21. I did not study in detail what implications this might have. It seems that W(1) is arbitrary since z1(1) = z3(1) = 0.
mahdi
2023-5-27
Ok, so whats your suggestion for that part? Cause rather than giving me numeric data, it gives me this:
z1 =
[0, -(pi*sin((1206735883187029*pi)/36028797018963968))/200, -(pi*(2 - 2^(1/2))^(1/2))/400, -(pi*2^(1/2))/400, -(pi*(2^(1/2) + 2)^(1/2))/400, -(pi*sin((4201915656573739*pi)/9007199254740992))/200, -pi/200]
mahdi
2023-5-27
It seems that W(1) is arbitrary since z1(1) = z3(1) = 0.
How can I implement W(1) = 0.005 into the equtions now? and does that fix the problem
syms U [7,1]
Q1 = (sum(C{1,1},1).*U.').';
V1 = (sum(C{1,2},1).*U.').';
syms W [7,1]
Q2 = (sum(C{1,1},1).*W.').';
V2 = (sum(C{1,2},1).*W.').';
syms PHI [7,1]
Q3 = (sum(C{1,1},1).*PHI.').';
V3 = (sum(C{1,2},1).*PHI.').';
z1 = X1.';
z2 = B.';
z3 = A.';
z4 = B1.';
U_g1 = 0; U_g2 = 0; U_g3 = 0;
W(1) = 0.005; %trying to fix that arbitary problem
eq1 = (((2*A_11).*z1.^2).*V1)+(((2*A_11).*z1).*Q1)-((2*A_11)*U)+...
(((A_11.*z2).*z1.^2).*V2)+((((A_11-(Nu*A_11)).*z3).*z1).*Q2) == 0;
eq2 = ((D_11.*z1.^2).*V3)+((D_11.*z1).*Q3)-((D_11+(A_55.*z1.^2)).*PHI)...
-((A_55.*z1.^2).*Q2) == 0;
eq3 = (((((2*A_55).*z1)+((2*A_11.*U_g2).*z1)+((A_11.*z2.^2).*z1)+...
(2*Nu*A_11.*U_g3)).*V2)+((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z1)+(A_11.*z3.^2)-(((A_11.*z4).*z3).*z1).*Q2)+((2*A_55.*z1).*Q3)...
+(2*A_55.*PHI)+((2*Nu*A_11-((2*A_11.*z4).*z1)).*Q1)-((2*Nu*A_11.*z4).*U)) == 0;
vars = [U;W;PHI];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
rank(Omega)
rank([Omega,b])
sol = Omega\b;
Cause this does not work, do I need a for loop for that?
Torsten
2023-5-27
编辑:Torsten
2023-5-27
Why do you always include code that I cannot execute - thus code that I cannot correct ?
Remove the line
W(1) = 0.005
and set
sol(7+1) = 0.005
This way, you have set W(1) = 0.005.
If you type
Omega(:,8)
in your code, you will see a column of zeros. This means that W(1) can be chosen arbitrarily since it is multiplied by 0 for all equations present.
mahdi
2023-5-29
So sorry for all the trouble I caused by not including the full code all the time, and thank you for your help.
mahdi
2023-5-30
编辑:Torsten
2023-5-30
I wanna solve these equations 100 times and use the solutions of W's and U's in the equations like below and I also need to compare the new solutions with the previous stage solutions with the help of the epsilon and stop the loop but I get errors, any help or suggestions would be appreciated. Considering C a 1*2 cell array with 2 7*7 matrices in it.
{}; % for example only
C{1,1} = rand(7);
C{1,2} = rand(7);
Nu = 0.285; A_11 = 64; A_55 = 37; D_11 = 2;
syms r
w_0 = 0.005; R = 0.5;
w_0 = w_0*cos(pi*r/(2*R));
O = diff(w_0,r);
P = diff(w_0,r,2);
N = 7; a = 0; b = 0.5;
X1 = a+(b-a)*(1-cos(((1:N)-1)*pi/(N-1)))/2;
A = double(subs(O,r,X1));
B = double(subs(P,r,X1));
B1 = double(subs(P,r,X1));
C_1 = cell2mat(C(1));
al = zeros(N);
al(2:N, 1:N) = C_1(2:N, 1:N);
bl = zeros(N);
bl(1:N-1, 1:N) = C_1(1:N-1, 1:N);
% Equations and their Solutions:
syms U [N,1]
Q1 = (sum(C{1,1},1).*U.').';
V1 = (sum(C{1,2},1).*U.').';
syms W [N,1]
Q2 = (sum(al).*W.').';
V2 = (sum(C{1,2},1).*W.').';
syms PHI [N,1]
Q3 = (sum(C{1,1},1).*PHI.').';
V3 = (sum(C{1,2},1).*PHI.').';
z1 = X1.';
z2 = B.'; % Initial guess for W_,rr
z3 = A.'; % Initial guess for W_,r
z4 = B1.'; % W_0,rr
U_g1 = 0; U_g2 = 0; U_g3 = 0; % Initial guesses for U_,rr and U_,r and U
eq1 = (((2*A_11).*z1.^2).*V1)+(((2*A_11).*z1).*Q1)-((2*A_11)*U)+...
(((A_11.*z2).*z1.^2).*V2)+((((A_11-(Nu*A_11)).*z3).*z1).*Q2) == 0;
eq2 = ((D_11.*z1.^2).*V3)+((D_11.*z1).*Q3)-((D_11+(A_55.*z1.^2)).*PHI)...
-((A_55.*z1.^2).*Q2) == 0;
eq3 = (((((2*A_55).*z1)+((2*A_11.*U_g2).*z1)+((A_11.*z2.^2).*z1)+...
(2*Nu*A_11.*U_g3)).*V2)+((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z1)+(A_11.*z3.^2)-(((A_11.*z4).*z3).*z1).*Q2)+((2*A_55.*z1).*Q3)...
+(2*A_55.*PHI)+((2*Nu*A_11-((2*A_11.*z4).*z1)).*Q1)-((2*Nu*A_11.*z4).*U)) == 0;
vars = [U;W;PHI];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
rank(Omega)
ans = 20
rank([Omega,b])
ans = 21
Omega(:,8) = -500; % for W(1)
rank(Omega)
ans = 21
sol = Omega\b;
JabejayiU = vpa(sol(1:N));
KhizW = vpa(sol(1+N:N+N));
TaghirzaviyePHI = vpa(sol(1+(2*N):N+(2*N)));
k = 100;
for i=1:k
U_g3 = JabejayiU.'; % must be new value in every loop
A = double(subs(O,r,KhizW.')); % new value in every loop
B = double(subs(P,r,KhizW.')); % new value in every loop
z2 = B.';
z3 = A.';
eq1 = (((2*A_11).*z1.^2).*V1)+(((2*A_11).*z1).*Q1)-((2*A_11)*U)+...
(((A_11.*z2).*z1.^2).*V2)+((((A_11-(Nu*A_11)).*z3).*z1).*Q2) == 0;
eq2 = ((D_11.*z1.^2).*V3)+((D_11.*z1).*Q3)-((D_11+(A_55.*z1.^2)).*PHI)...
-((A_55.*z1.^2).*Q2) == 0;
eq3 = (((((2*A_55).*z1)+((2*A_11.*U_g2).*z1)+((A_11.*z2.^2).*z1)+...
(2*Nu*A_11.*U_g3)).*V2)+((2*A_55)+(2*A_11.*U_g2)...
+((2*A_11.*U_g1).*z1)+(A_11.*z3.^2)-(((A_11.*z4).*z3).*z1).*Q2)+((2*A_55.*z1).*Q3)...
+(2*A_55.*PHI)+((2*Nu*A_11-((2*A_11.*z4).*z1)).*Q1)-((2*Nu*A_11.*z4).*U)) == 0;
vars = [U;W;PHI];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
rank(Omega)
rank([Omega,b])
sol = Omega\b;
JabejayiU = vpa(sol(1:N));
KhizW = vpa(sol(1+N:N+N));
TaghirzaviyePHI = vpa(sol(1+(2*N):N+(2*N)));
epsilon = 0.0001;
if abs((JabejayiU(i+1)-JabejayiU(i))/JabejayiU(i)) < epsilon
abs((KhizW(i+1)-KhizW(i))/KhizW(i)) < epsilon
abs((TaghirzaviyePHI(i+1)-TaghirzaviyePHI(i))/TaghirzaviyePHI(i)) < epsilon
end
end
ans = 21
ans = 22
Warning: Solution does not exist because the system is inconsistent.
Error using mupadengine/feval_internal
Unable to convert to matrix form because the system does not seem to be linear.
Unable to convert to matrix form because the system does not seem to be linear.
Error in sym/equationsToMatrix (line 61)
T = eng.feval_internal('symobj::equationsToMatrix',eqns,vars);
mahdi
2023-5-31
编辑:mahdi
2023-5-31
Man you are one lovely person, cause you are still answering my questions and helping me out.
But my teacher says something completely different from what I coded so I changed my code to this now:
Lets just forget about the i=1:100 loop for now.
{}; % for example only
C{1,1} = rand(7);
C{1,2} = rand(7);
Nu = 0.285; A_11 = 64; A_55 = 37; D_11 = 2; w0c = 0.005; R = 0.5;
A11 = 6.4*10^4; A55 = 3.7*10^4; D11 = 2*10^3;
N = 7; a = 0; b = 0.5;
X1 = a+(b-a)*(1-cos(((1:N)-1)*pi/(N-1)))/2;
syms r
w0 = w0c*cos(pi*r/(2*R));
w0onX1 = double(subs(w0,r,X1));
A = C{1,1}*w0onX1.';
B = C{1,2}*w0onX1.';
B1 = C{1,2}*w0onX1.';
C_1 = cell2mat(C(1));
al = zeros(N);
al(2:N, 1:N) = C_1(2:N, 1:N);
bl = zeros(N);
bl(1:N-1, 1:N) = C_1(1:N-1, 1:N);
syms U [N,1]
syms W [N,1]
syms PHI [N,1]
z1 = X1;
z2 = B.';
z3 = A.';
z4 = B1.';
Ug1 = 0; Ug2 = 0; Ug3 = 0;
atm=((2*A11*z1.^2).*(C{1,2}))+(((2*A11*z1).*(C{1,1}))-(2*A11));
atm2=(((A11.*z2).*(z1.^2)).*(C{1,2}))+(((((A11-(Nu*A11)).*z3).*z1).*(al)));
atm3=((D11*z1.^2).*C{1,2})+((D11*z1).*C{1,1})-(D11+(A55*z1.^2));
atm4=((A55*z1.^2).*(al));
atm5=((((2*A55*z1)+((2*A11*Ug2).*z1)+((A11*z2.^2).*z1)+(2*Nu*A11*Ug3)).*C{1,2})+...
((2*A55)+(2*A11*Ug2)+((2*A11*Ug1).*z1)+(A11*z3.^2)-(((A11*z4).*z3).*z1).*(al)));
atm6=(((2*A55*z1).*C{1,1})+(2*A55));
atm7=((((2*Nu*A11)-((2*A11*z4).*z1)).*C{1,1})-(2*Nu*A11*z4));
eq1 = (atm*U)+(atm2*W) == 0;
eq2 = (atm3*PHI)-(atm4*W) == 0;
eq3 = (atm5*W)+(atm6*PHI)-(atm7*U) == 0;
vars = [U;W;PHI];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
rank(Omega)
ans = 21
rank([Omega,b])
ans = 21
sol = Omega\b;
So this gives the right number of variables and equations now but the solutions are all zero cause the right hand side is zero.
Well he says that I must have made something like this: A matrix with only coefficients multiplicated by a 21*1 matrix of variables = a 21*1 matrix of zeroes except for 1 element which is like (coefficient*U7=P) and P=100,500,1000,5000,10000. I kinda dont know how to explain it well so I uploaded a picture too.
Then solve that like [A]*[U]=[B] and [U]=[A]\[B]
Do you have any suggestions?
Torsten
2023-5-31
{}; % for example only
C{1,1} = rand(7);
C{1,2} = rand(7);
Nu = 0.285; A_11 = 64; A_55 = 37; D_11 = 2; w0c = 0.005; R = 0.5;
A11 = 6.4*10^4; A55 = 3.7*10^4; D11 = 2*10^3;
N = 7; a = 0; b = 0.5;
X1 = a+(b-a)*(1-cos(((1:N)-1)*pi/(N-1)))/2;
syms r
w0 = w0c*cos(pi*r/(2*R));
w0onX1 = double(subs(w0,r,X1));
A = C{1,1}*w0onX1.';
B = C{1,2}*w0onX1.';
B1 = C{1,2}*w0onX1.';
C_1 = cell2mat(C(1));
al = zeros(N);
al(2:N, 1:N) = C_1(2:N, 1:N);
bl = zeros(N);
bl(1:N-1, 1:N) = C_1(1:N-1, 1:N);
syms U [N,1]
syms W [N,1]
syms PHI [N,1]
z1 = X1;
z2 = B.';
z3 = A.';
z4 = B1.';
Ug1 = 0; Ug2 = 0; Ug3 = 0;
atm=((2*A11*z1.^2).*(C{1,2}))+(((2*A11*z1).*(C{1,1}))-(2*A11));
atm2=(((A11.*z2).*(z1.^2)).*(C{1,2}))+(((((A11-(Nu*A11)).*z3).*z1).*(al)));
atm3=((D11*z1.^2).*C{1,2})+((D11*z1).*C{1,1})-(D11+(A55*z1.^2));
atm4=((A55*z1.^2).*(al));
atm5=((((2*A55*z1)+((2*A11*Ug2).*z1)+((A11*z2.^2).*z1)+(2*Nu*A11*Ug3)).*C{1,2})+...
((2*A55)+(2*A11*Ug2)+((2*A11*Ug1).*z1)+(A11*z3.^2)-(((A11*z4).*z3).*z1).*(al)));
atm6=(((2*A55*z1).*C{1,1})+(2*A55));
atm7=((((2*Nu*A11)-((2*A11*z4).*z1)).*C{1,1})-(2*Nu*A11*z4));
eq1 = (atm*U)+(atm2*W) == 0;
eq2 = (atm3*PHI)-(atm4*W) == 0;
eq3 = (atm5*W)+(atm6*PHI)-(atm7*U) == 0;
vars = [U;W;PHI];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
P = [100,500,1000,5000,10000];
for i = 1:numel(P)
b(7) = P(i);
sol(i,:) = double(Omega\b);
end
sol
sol = 5×21
0.0872 -0.0225 -0.0425 -0.0297 -0.0235 0.0155 0.0027 3.6938 -1.9387 -0.4951 0.0065 0.0098 0.0347 0.0037 0.0035 0.1348 -1.1442 -0.4329 -0.4464 0.5689 0.1580
0.4358 -0.1127 -0.2127 -0.1483 -0.1173 0.0775 0.0137 18.4692 -9.6934 -2.4757 0.0323 0.0490 0.1736 0.0184 0.0177 0.6740 -5.7209 -2.1647 -2.2318 2.8443 0.7900
0.8715 -0.2255 -0.4254 -0.2967 -0.2347 0.1549 0.0273 36.9384 -19.3869 -4.9513 0.0646 0.0980 0.3472 0.0368 0.0354 1.3480 -11.4419 -4.3294 -4.4635 5.6887 1.5799
4.3576 -1.1275 -2.1271 -1.4835 -1.1735 0.7747 0.1367 184.6918 -96.9343 -24.7566 0.3231 0.4901 1.7360 0.1840 0.1771 6.7398 -57.2093 -21.6470 -22.3175 28.4433 7.8997
8.7151 -2.2550 -4.2542 -2.9670 -2.3470 1.5494 0.2733 369.3836 -193.8685 -49.5133 0.6462 0.9802 3.4720 0.3681 0.3541 13.4797 -114.4186 -43.2940 -44.6350 56.8866 15.7994
mahdi
2023-6-8
编辑:mahdi
2023-6-8
I made a loop like this but I get an error:
I want the equations to be solved 50 times and some of their solutions be used in the next iteration and also be compared with that error function.
Any suggestions
{}; % for example only
C{1,1} = rand(7);
C{1,2} = rand(7);
Nu = 0.285; A_11 = 64; A_55 = 37; D_11 = 2; w0c = 0.005; R = 0.5;
A11 = 6.4*10^4; A55 = 3.7*10^4; D11 = 2*10^3;
N = 7; a = 0; b = 0.5;
X1 = a+(b-a)*(1-cos(((1:N)-1)*pi/(N-1)))/2;
syms r
w0 = w0c*cos(pi*r/(2*R));
w0onX1 = double(subs(w0,r,X1));
A = C{1,1}*w0onX1.';
B = C{1,2}*w0onX1.';
B1 = C{1,2}*w0onX1.';
C_1 = cell2mat(C(1));
al = zeros(N);
al(2:N,1:N) = C_1(2:N,1:N);
bl = zeros(N);
bl(1:N-1,1:N) = C_1(1:N-1,1:N);
syms U [N,1]
syms W [N,1]
syms PHI [N,1]
z1 = X1;
z2 = B.';
z3 = A.';
z4 = B1.';
Ug1 = 0; Ug2 = 0; Ug3 = 0;
atm=((2*A11*z1.^2).*(C{1,2}))+(((2*A11*z1).*(C{1,1}))-(2*A11));
atm2=(((A11.*z2).*(z1.^2)).*(C{1,2}))+(((((A11-(Nu*A11)).*z3).*z1).*(al)));
atm3=((D11*z1.^2).*C{1,2})+((D11*z1).*C{1,1})-(D11+(A55*z1.^2));
atm4=((A55*z1.^2).*(al));
atm5=((((2*A55*z1)+((2*A11*Ug2).*z1)+((A11*z2.^2).*z1)+(2*Nu*A11*Ug3)).*C{1,2})+...
((2*A55)+(2*A11*Ug2)+((2*A11*Ug1).*z1)+(A11*z3.^2)+(((A11*z4).*z3).*z1).*(al)));
atm6=(((2*A55*z1).*C{1,1})+(2*A55));
atm7=((((2*Nu*A11)+((2*A11*z4).*z1)).*C{1,1})+(2*Nu*A11*z4));
eq1 = (atm*U)+(atm2*W) == 0;
eq2 = (atm3*PHI)-(atm4*W) == 0;
eq3 = (atm5*W)+(atm6*PHI)+(atm7*U) == 0;
vars = [U;W;PHI];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
b(7) = -100;
sol = Omega\b;
JabejayiU = double(sol(1:N));
KhizW = double(sol(1+N:N+N));
TaghirzaviyePHI = double(sol(1+(2*N):N+(2*N)));
for i=1:50
Ug1 = (JabejayiU.')*C{1,2}; Ug2 = (JabejayiU.')*C{1,1}; Ug3 = (JabejayiU.');
z2 = (KhizW.')*C{1,2}; z3 = (KhizW.')*C{1,1};
eq1 = ((atm*U)+(atm2*W)) == 0;
eq2 = ((atm3*PHI)+(atm4*W)) == 0;
eq3 = ((atm5*W)+(atm6*PHI)+(atm7*U)) == 0;
vars = [U;W;PHI];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
b(7) = -100;
sol = Omega\b;
epsilon = 0.0001;
JabejayiU = double(sol(1:N));
KhizW = double(sol(1+N:N+N));
TaghirzaviyePHI = double(sol(1+(2*N):N+(2*N)));
if abs((JabejayiU(i+1)-JabejayiU(i))/JabejayiU(i)) < epsilon
abs((KhizW(i+1)-KhizW(i))/KhizW(i)) < epsilon;
abs((TaghirzaviyePHI(i+1)-TaghirzaviyePHI(i))/TaghirzaviyePHI(i)) < epsilon;
end
end
Index exceeds the number of array elements. Index must not exceed 7.
Torsten
2023-6-8
You must compare the complete solutions for JabejayiU, KhizW and TaghirzaviyePHI for iteration step i and iteration step i+1. What you do is to compare element i+1 with element i of the arrays where i is the iteration step. This doesn't make sense, and MATLAB errors as soon as the iteration number exceeds the number of array elements (7).
Torsten
2023-6-9
You do not change the equations you solve within the loop. So you get the same solution 50 times. What should the iteration loop be good for ?
mahdi
2023-6-9
编辑:mahdi
2023-6-9
Well yes the equations do not change but Ug1, Ug2, Ug3, z2 and z3 must be new in every iteration loop and they must be like I wrote:
Ug1 = (JabejayiU.')*C{1,2}; Ug2 = (JabejayiU.')*C{1,1}; Ug3 = (JabejayiU.');
z2 = (KhizW.')*C{1,2}; z3 = (KhizW.')*C{1,1};
As you can see in every iteration loop JabejayiU and KhizW which are the solutions must be used again in Ug1,...,z3 like the formula above and then equations be solved again and so on.
mahdi
2023-6-9
Is that so? Cause first Ug1, Ug2 and Ug3 are all just 0, z2 and z3 are 1*7 matrces and then equations are gonna be solved and giving solutions JabejayiU, KhizW and TaghirzaviyePHI then 2 of these solutions which are gonna be JabejayiU and KhizW (7*1 matrces and because of it I used Jabejayi.' and KhizW.' to make them 1*7) are gonna be used in the next iteration loop (just like the formula of my last comment) giving new solutions for JabejayiU and KhizW and then this process must be repeated 50 times untill at some point like iteration 5 or 7 the resulting solutions wont be much different and then we can tell the solutions are convergence.
Torsten
2023-6-9
I don't see where the new Ug1, Ug2, Ug3, z2 and z3 are used in eq1, eq2 and eq3 within the loop. If you want to use the new values for Ug1, Ug2, Ug3, z2 and z3 in eq1, eq2 and eq3, you have to recalculate atm, atm2,...,atm7.
mahdi
2023-6-9
So I just need to write atm,...,atm7 in the for loop before the eq1 then. Consider this done already, but this does not help me with the other problem.
Torsten
2023-6-9
编辑:Torsten
2023-6-9
Your system now is nonlinear in U,W and PHI. The solution method used is no longer adequate. Don't iterate because there is no guarantee for convergence.
Set up your system as usual and use matlabFunction and fsolve instead of equationsToMatrix to solve it.
rng(1)
{}; % for example only
C{1,1} = rand(7);
C{1,2} = rand(7);
Nu = 0.285; A_11 = 64; A_55 = 37; D_11 = 2; w0c = 0.005; R = 0.5;
A11 = 6.4*10^4; A55 = 3.7*10^4; D11 = 2*10^3;
N = 7; a = 0; b = 0.5;
X1 = a+(b-a)*(1-cos(((1:N)-1)*pi/(N-1)))/2;
syms r
w0 = w0c*cos(pi*r/(2*R));
w0onX1 = double(subs(w0,r,X1));
A = C{1,1}*w0onX1.';
B = C{1,2}*w0onX1.';
B1 = C{1,2}*w0onX1.';
C_1 = cell2mat(C(1));
al = zeros(N);
al(2:N,1:N) = C_1(2:N,1:N);
bl = zeros(N);
bl(1:N-1,1:N) = C_1(1:N-1,1:N);
syms U [N,1]
syms W [N,1]
syms PHI [N,1]
z1 = X1;
z4 = B1.';
Ug1 = U.'*C{1,2}; Ug2 = U.'*C{1,1}; Ug3 = U.';
z2 = W.'*C{1,2}; z3 = W.'*C{1,1};
atm=((2*A11*z1.^2).*(C{1,2}))+(((2*A11*z1).*(C{1,1}))-(2*A11));
atm2=(((A11.*z2).*(z1.^2)).*(C{1,2}))+(((((A11-(Nu*A11)).*z3).*z1).*(al)));
atm3=((D11*z1.^2).*C{1,2})+((D11*z1).*C{1,1})-(D11+(A55*z1.^2));
atm4=((A55*z1.^2).*(al));
atm5=((((2*A55*z1)+((2*A11*Ug2).*z1)+((A11*z2.^2).*z1)+(2*Nu*A11*Ug3)).*C{1,2})+...
((2*A55)+(2*A11*Ug2)+((2*A11*Ug1).*z1)+(A11*z3.^2)+(((A11*z4).*z3).*z1).*(al)));
atm6=(((2*A55*z1).*C{1,1})+(2*A55));
atm7=((((2*Nu*A11)+((2*A11*z4).*z1)).*C{1,1})+(2*Nu*A11*z4));
eq1 = (atm*U)+(atm2*W) == 0;
eq2 = (atm3*PHI)-(atm4*W) == 0;
eq3 = (atm5*W)+(atm6*PHI)+(atm7*U) == 0;
eq1(7) = subs(eq1(7),lhs(eq1(7)),lhs(eq1(7))+sym('100'));
f = matlabFunction([lhs(eq1);lhs(eq2);lhs(eq3)],"Vars",{[U;W;PHI]})
f = function_handle with value:
@(in1)[in1(1,:).*-1.28e+5-in1(2,:).*1.264360769370139e+5-in1(3,:).*1.267334399697794e+5-in1(4,:).*8.978841425467622e+4-in1(5,:).*1.135883543854816e+5-in1(6,:).*1.016404530840053e+5-in1(7,:).*7.966792822925179e+4+in1(11,:).*(in1(8,:).*3.264540082183265e+3+in1(9,:).*4.967796441173845e+2+in1(10,:).*5.032900754425425e+2+in1(11,:).*2.917595346427745e+3+in1(12,:).*1.437048064793434e+3+in1(13,:).*5.975251935053092e+2+in1(14,:).*3.351652111902892e+3)+in1(10,:).*(in1(8,:).*1.715444673058387e+2+in1(9,:).*2.06896871366325e+1+in1(10,:).*2.219571819441791e+2+in1(11,:).*2.749299812730016e+2+in1(12,:).*2.132563963196658e+2+in1(13,:).*3.912315660352179e+2+in1(14,:).*2.429389382124382e+2)+in1(12,:).*(in1(8,:).*1.088469838840782e+3+in1(9,:).*2.349961348990294e+3+in1(10,:).*2.272295821940149e+3+in1(11,:).*2.764653320039682e+3+in1(12,:).*1.952027111574643e+3+in1(13,:).*2.350369270700102e+3+in1(14,:).*1.092014386755356e+3)+in1(9,:).*(in1(8,:).*2.366500685597495e+1+in1(9,:).*6.048120973642112+in1(10,:).*2.429105019193534e+1+in1(11,:).*2.884389231232433e+1+in1(12,:).*4.218203616552341+in1(13,:).*1.706730066588052e+1+in1(14,:).*2.862302834161352e+1)+in1(14,:).*(in1(8,:).*1.442447763840166e+4+in1(9,:).*6.83498780961897e+3+in1(10,:).*8.786795619332904e+3+in1(11,:).*6.200344164878128e+3+in1(12,:).*3.600873148232448e+3+in1(13,:).*1.372398642095686e+4+in1(14,:).*8.71524016857475e+3)+in1(13,:).*(in1(8,:).*1.014823745306438e+3+in1(9,:).*3.368183263216046e+3+in1(10,:).*1.609456146962601e+3+in1(11,:).*3.62742280490641e+3+in1(12,:).*2.494281643906298e+3+in1(13,:).*2.337333778655289e+3+in1(14,:).*4.314001686522385e+2);in1(1,:).*-1.28e+5-in1(2,:).*1.262779143357423e+5-in1(3,:).*1.171726129192534e+5-in1(4,:).*1.168706286657343e+5-in1(5,:).*7.233454197893738e+4-in1(6,:).*6.205092789871104e+4-in1(7,:).*6.308893470774664e+4+in1(10,:).*(in1(8,:).*1.257231380658536e+2+in1(9,:).*2.573787970272283e+3+in1(10,:).*1.627164742771721e+3+in1(11,:).*2.175775762091544e+3+in1(12,:).*5.641148124792556e+2+in1(13,:).*8.069196768584876e+2+in1(14,:).*3.100215255798277e+3)+in1(12,:).*(in1(8,:).*4.909123265070745e+3+in1(9,:).*1.830613316655202e+4+in1(10,:).*6.387773766591589e+3+in1(11,:).*1.231441150570409e+4+in1(12,:).*1.86487116048835e+4+in1(13,:).*1.310857789722395e+4+in1(14,:).*1.278345846890527e+4)+in1(11,:).*(in1(8,:).*3.968551817850398e+3+in1(9,:).*1.199402321923257e+3+in1(10,:).*2.558960898349173e+3+in1(11,:).*3.586340443963706e+3+in1(12,:).*3.426359160815992e+3+in1(13,:).*3.958602712300056e+2+in1(14,:).*6.500696313183545e+2)+in1(14,:).*(in1(8,:).*1.18994448238035e+4+in1(9,:).*1.74921021640202e+4+in1(10,:).*6.027720710344429e+3+in1(11,:).*1.102639095855662e+4+in1(12,:).*1.811432466800663e+4+in1(13,:).*1.180535558027109e+4+in1(14,:).*9.326537471040024e+3)+in1(9,:).*(in1(8,:).*2.161882154389071e+2+in1(9,:).*2.42825333879213e+2+in1(10,:).*3.33869774335959e+2+in1(11,:).*2.622894850373104e+2+in1(12,:).*4.177691881801658e+2+in1(13,:).*1.286920741631465e+2+in1(14,:).*5.413101979355649e+2)+in1(13,:).*(in1(8,:).*7.992052932363333e+3+in1(9,:).*2.123958986209751e+4+in1(10,:).*1.757316612549308e+4+in1(11,:).*1.230737044129967e+4+in1(12,:).*1.927181624876517e+4+in1(13,:).*2.22492941905784e+4+in1(14,:).*1.239614947982919e+4);in1(1,:).*-1.28e+5-in1(2,:).*1.25605371392226e+5-in1(3,:).*1.202513303502949e+5-in1(4,:).*1.047314655205759e+5-in1(5,:).*1.102113882397055e+5-in1(6,:).*6.623707120760735e+4-in1(7,:).*1.028850679186427e+5+in1(14,:).*(in1(8,:).*9.449151333227075e+3+in1(9,:).*6.027720710344429e+3+in1(10,:).*5.596353057463494e+3+in1(11,:).*4.834833278070911e+3+in1(12,:).*4.339433980698522e+3+in1(13,:).*9.053546455030521e+3+in1(14,:).*5.988635082229147e+3)+in1(13,:).*(in1(8,:).*7.231007660972983e+3+in1(9,:).*1.757316612549308e+4+in1(10,:).*1.742306208105052e+4+in1(11,:).*6.078741218930922e+3+in1(12,:).*1.73211410786888e+4+in1(13,:).*2.132545412589933e+4+in1(14,:).*1.401426760882863e+4)+in1(10,:).*(in1(8,:).*2.873309141781141e+2+in1(9,:).*1.627164742771721e+3+in1(10,:).*1.283284622883137e+3+in1(11,:).*1.689308549739376e+3+in1(12,:).*6.110285263632219e+2+in1(13,:).*9.790699145189742e+2+in1(14,:).*2.225696786829062e+3)+in1(12,:).*(in1(8,:).*2.558907016284718e+3+in1(9,:).*6.387773766591589e+3+in1(10,:).*4.909630903876981e+3+in1(11,:).*6.482180547805791e+3+in1(12,:).*5.691627598651523e+3+in1(13,:).*5.806430867467453e+3+in1(14,:).*3.447329168577428e+3)+in1(11,:).*(in1(8,:).*8.17208705990996e+3+in1(9,:).*2.558960898349173e+3+in1(10,:).*5.560904915286832e+3+in1(11,:).*7.390953744169242e+3+in1(12,:).*7.306985414080078e+3+in1(13,:).*7.656843873815043e+2+in1(14,:).*8.260405792494059e+2)+in1(9,:).*(in1(8,:).*3.096647431905299e+2+in1(9,:).*3.33869774335959e+2+in1(10,:).*4.699035524481027e+2+in1(11,:).*3.757894696283357e+2+in1(12,:).*5.702031815397364e+2+in1(13,:).*1.863611797993887e+2+in1(14,:).*7.54554120360723e+2);in1(1,:).*-1.28e+5-in1(2,:).*1.261023539146309e+5-in1(3,:).*1.177333734544284e+5-in1(4,:).*9.349641681685059e+4-in1(5,:).*9.188732435805292e+4-in1(6,:).*1.000309789935144e+5-in1(7,:).*8.627443663638107e+4+in1(12,:).*(in1(8,:).*3.991883060401115e+3+in1(9,:).*1.231441150570409e+4+in1(10,:).*6.482180547805791e+3+in1(11,:).*1.006507765918129e+4+in1(12,:).*1.187994901428235e+4+in1(13,:).*9.82257086692599e+3+in1(14,:).*7.773303412913462e+3)+in1(9,:).*(in1(8,:).*2.508620278130593e+2+in1(9,:).*2.622894850373104e+2+in1(10,:).*3.757894696283357e+2+in1(11,:).*3.044829897495721e+2+in1(12,:).*4.453856664592989e+2+in1(13,:).*1.521601643777014e+2+in1(14,:).*5.990656179124037e+2)+in1(10,:).*(in1(8,:).*3.624526722160685e+2+in1(9,:).*2.175775762091544e+3+in1(10,:).*1.689308549739376e+3+in1(11,:).*2.226031666640478e+3+in1(12,:).*7.904159920059069e+2+in1(13,:).*1.260091022292998e+3+in1(14,:).*2.9482961724285e+3)+in1(14,:).*(in1(8,:).*9.074279194571261e+3+in1(9,:).*1.102639095855662e+4+in1(10,:).*4.834833278070911e+3+in1(11,:).*7.255135411172964e+3+in1(12,:).*1.085893946310561e+4+in1(13,:).*8.908139075841853e+3+in1(14,:).*6.695300615020857e+3)+in1(11,:).*(in1(8,:).*1.262528123770883e+4+in1(9,:).*3.586340443963706e+3+in1(10,:).*7.390953744169242e+3+in1(11,:).*1.139410567757531e+4+in1(12,:).*1.025355295895174e+4+in1(13,:).*1.386666223321988e+3+in1(14,:).*3.387008427019715e+3)+in1(13,:).*(in1(8,:).*3.750601969715906e+3+in1(9,:).*1.230737044129967e+4+in1(10,:).*6.078741218930922e+3+in1(11,:).*1.297313149915331e+4+in1(12,:).*9.20851154915234e+3+in1(13,:).*8.740707844468265e+3+in1(14,:).*1.8341006984106e+3);in1(1,:).*-1.28e+5-in1(2,:).*1.250476412488472e+5-in1(3,:).*1.247240307583596e+5-in1(4,:).*9.619117207199262e+4-in1(5,:).*7.079520282577496e+4-in1(6,:).*6.472567683183544e+4-in1(7,:).*6.2265024434265e+4+in1(9,:).*(in1(8,:).*3.671312443642008e+2+in1(9,:).*4.177691881801658e+2+in1(10,:).*5.702031815397364e+2+in1(11,:).*4.453856664592989e+2+in1(12,:).*7.203796413600615e+2+in1(13,:).*2.177609945481449e+2+in1(14,:).*9.273140105787016e+2)+in1(12,:).*(in1(8,:).*4.743594990136405e+3+in1(9,:).*1.86487116048835e+4+in1(10,:).*5.691627598651523e+3+in1(11,:).*1.187994901428235e+4+in1(12,:).*1.924590155652118e+4+in1(13,:).*1.297892003126547e+4+in1(14,:).*1.333103382705727e+4)+in1(10,:).*(in1(8,:).*2.352489976048832e+2+in1(9,:).*5.641148124792556e+2+in1(10,:).*6.110285263632219e+2+in1(11,:).*7.904159920059069e+2+in1(12,:).*3.778433749771653e+2+in1(13,:).*6.454396834262847e+2+in1(14,:).*9.4501933353004e+2)+in1(13,:).*(in1(8,:).*7.546816033147669e+3+in1(9,:).*1.927181624876517e+4+in1(10,:).*1.73211410786888e+4+in1(11,:).*9.20851154915234e+3+in1(12,:).*1.814305804063192e+4+in1(13,:).*2.158001402592649e+4+in1(14,:).*1.304112526190363e+4)+in1(14,:).*(in1(8,:).*9.430929381586814e+3+in1(9,:).*1.811432466800663e+4+in1(10,:).*4.339433980698522e+3+in1(11,:).*1.085893946310561e+4+in1(12,:).*1.978739626476403e+4+in1(13,:).*9.529851318680608e+3+in1(14,:).*8.158105402662894e+3)+in1(11,:).*(in1(8,:).*1.134652773234728e+4+in1(9,:).*3.426359160815992e+3+in1(10,:).*7.306985414080078e+3+in1(11,:).*1.025355295895174e+4+in1(12,:).*9.788260649189233e+3+in1(13,:).*1.133398128390026e+3+in1(14,:).*1.87509309115714e+3);in1(1,:).*-1.28e+5-in1(2,:).*1.270640191918784e+5-in1(3,:).*1.229411866626603e+5-in1(4,:).*1.239557516592317e+5-in1(5,:).*8.889110250880313e+4-in1(6,:).*5.163411411623892e+4-in1(7,:).*8.03005498460924e+4+in1(10,:).*(in1(8,:).*4.222655575663895e+2+in1(9,:).*8.069196768584876e+2+in1(10,:).*9.790699145189742e+2+in1(11,:).*1.260091022292998e+3+in1(12,:).*6.454396834262847e+2+in1(13,:).*1.116736076917481e+3+in1(14,:).*1.461412930344928e+3)+in1(11,:).*(in1(8,:).*1.539552475206125e+3+in1(9,:).*3.958602712300056e+2+in1(10,:).*7.656843873815043e+2+in1(11,:).*1.386666223321988e+3+in1(12,:).*1.133398128390026e+3+in1(13,:).*1.921080466498272e+2+in1(14,:).*6.514663707903619e+2)+in1(13,:).*(in1(8,:).*8.997725063829892e+3+in1(9,:).*2.22492941905784e+4+in1(10,:).*2.132545412589933e+4+in1(11,:).*8.740707844468265e+3+in1(12,:).*2.158001402592649e+4+in1(13,:).*2.625773536897018e+4+in1(14,:).*1.678703360262024e+4)+in1(12,:).*(in1(8,:).*3.904167528899096e+3+in1(9,:).*1.310857789722395e+4+in1(10,:).*5.806430867467453e+3+in1(11,:).*9.82257086692599e+3+in1(12,:).*1.297892003126547e+4+in1(13,:).*9.953222047088267e+3+in1(14,:).*8.688087305393721e+3)+in1(9,:).*(in1(8,:).*1.253514163321741e+2+in1(9,:).*1.286920741631465e+2+in1(10,:).*1.863611797993887e+2+in1(11,:).*1.521601643777014e+2+in1(12,:).*2.177609945481449e+2+in1(13,:).*7.637569506222814e+1+in1(14,:).*2.958079090488983e+2)+in1(14,:).*(in1(8,:).*1.56079787917089e+4+in1(9,:).*1.180535558027109e+4+in1(10,:).*9.053546455030521e+3+in1(11,:).*8.908139075841853e+3+in1(12,:).*9.529851318680608e+3+in1(13,:).*1.50299809854452e+4+in1(14,:).*1.022524835683261e+4);in1(1,:).*-1.28e+5-in1(2,:).*1.241356348965077e+5-in1(3,:).*1.140149768201874e+5-in1(4,:).*1.193301782933795e+5-in1(5,:).*8.850972837477611e+4-in1(6,:).*8.012851378372216e+4-in1(7,:).*9.122463475696053e+4+in1(10,:).*(in1(8,:).*3.683812066057913e+2+in1(9,:).*3.100215255798277e+3+in1(10,:).*2.225696786829062e+3+in1(11,:).*2.9482961724285e+3+in1(12,:).*9.4501933353004e+2+in1(13,:).*1.461412930344928e+3+in1(14,:).*4.011664279124778e+3)+in1(13,:).*(in1(8,:).*5.470607454893864e+3+in1(9,:).*1.239614947982919e+4+in1(10,:).*1.401426760882863e+4+in1(11,:).*1.8341006984106e+3+in1(12,:).*1.304112526190363e+4+in1(13,:).*1.678703360262024e+4+in1(14,:).*1.213259534395224e+4)+in1(14,:).*(in1(8,:).*1.05617674745316e+4+in1(9,:).*9.326537471040024e+3+in1(10,:).*5.988635082229147e+3+in1(11,:).*6.695300615020857e+3+in1(12,:).*8.158105402662894e+3+in1(13,:).*1.022524835683261e+4+in1(14,:).*7.160529543862636e+3)+in1(9,:).*(in1(8,:).*4.937006695786826e+2+in1(9,:).*5.413101979355649e+2+in1(10,:).*7.54554120360723e+2+in1(11,:).*5.990656179124037e+2+in1(12,:).*9.273140105787016e+2+in1(13,:).*2.958079090488983e+2+in1(14,:).*1.216446186255984e+3)+in1(11,:).*(in1(8,:).*3.784258525598922e+3+in1(9,:).*6.500696313183545e+2+in1(10,:).*8.260405792494059e+2+in1(11,:).*3.387008427019715e+3+in1(12,:).*1.87509309115714e+3+in1(13,:).*6.514663707903619e+2+in1(14,:).*3.458537825341188e+3)+in1(12,:).*(in1(8,:).*3.108344940981892e+3+in1(9,:).*1.278345846890527e+4+in1(10,:).*3.447329168577428e+3+in1(11,:).*7.773303412913462e+3+in1(12,:).*1.333103382705727e+4+in1(13,:).*8.688087305393721e+3+in1(14,:).*9.309958501586755e+3)+1.0e+2;in1(15,:).*-2.0e+3-in1(16,:).*2.017071209637813e+3-in1(17,:).*2.558334999527803e+3-in1(18,:).*3.715443972729316e+3-in1(19,:).*6.977943037273149e+3-in1(20,:).*9.640374571940614e+3-in1(21,:).*1.049481137858206e+4;in1(15,:).*-2.0e+3-in1(16,:).*2.014599918992944e+3-in1(17,:).*2.408947076863334e+3-in1(18,:).*4.138603572902099e+3-in1(19,:).*6.333352218420896e+3-in1(20,:).*9.021788240920392e+3-in1(21,:).*1.023576460480854e+4-in1(9,:).*1.646882891118387e+1-in1(10,:).*3.876140293218886e+2-in1(11,:).*7.247934119932489e+2-in1(12,:).*4.569085213156164e+3-in1(13,:).*5.527871941020398e+3-in1(14,:).*7.300833788176269e+3;in1(15,:).*-2.0e+3-in1(16,:).*2.004091435500502e+3-in1(17,:).*2.457052036723358e+3-in1(18,:).*3.948929148758998e+3-in1(19,:).*6.925177941245397e+3-in1(20,:).*9.087196730121897e+3-in1(21,:).*1.085757918622879e+4-in1(9,:).*2.236493962613788e+1-in1(10,:).*2.412543388684952e+2-in1(11,:).*1.600996048735288e+3-in1(12,:).*5.117108697875887e+2-in1(13,:).*6.720608300585919e+3-in1(14,:).*9.548405608431888e+2;in1(15,:).*-2.0e+3-in1(16,:).*2.011856787413079e+3-in1(17,:).*2.417708960225444e+3-in1(18,:).*3.773381512763291e+3-in1(19,:).*6.638864443094576e+3-in1(20,:).*9.615226539276695e+3-in1(21,:).*1.059803807244345e+4-in1(9,:).*1.739971944928397e+1-in1(10,:).*3.229925570702001e+2-in1(11,:).*2.026649914684588e+3-in1(12,:).*2.191075611354412e+3-in1(13,:).*1.472616439455817e+2-in1(14,:).*4.143015117127123e+3;in1(15,:).*-2.0e+3-in1(16,:).*1.995376902010209e+3-in1(17,:).*2.526937980599369e+3-in1(18,:).*3.815487063624885e+3-in1(19,:).*6.309300044152733e+3-in1(20,:).*9.06358119300046e+3-in1(21,:).*1.022289100678539e+4-in1(9,:).*2.844175354978936e+1-in1(10,:).*8.116119887536949e+1-in1(11,:).*2.068777909352646e+3-in1(12,:).*4.984018961564329e+3-in1(13,:).*6.040343928309443e+3-in1(14,:).*8.404508403611135e+3;in1(15,:).*-2.0e+3-in1(16,:).*2.026882807370071e+3-in1(17,:).*2.499081041604067e+3-in1(18,:).*4.249308619675495e+3-in1(19,:).*6.592048476700048e+3-in1(20,:).*8.859025525569265e+3-in1(21,:).*1.050469609134519e+4-in1(9,:).*8.48630328850359-in1(10,:).*1.145274233771955e+2-in1(11,:).*1.966647387926114e+2-in1(12,:).*2.774125623375229e+3-in1(13,:).*7.962549279275694e+3-in1(14,:).*2.715930872456535e+3;in1(15,:).*-2.0e+3-in1(16,:).*1.981126802754904e+3-in1(17,:).*2.359609012815429e+3-in1(18,:).*4.177034035834055e+3-in1(19,:).*6.586089505855876e+3-in1(20,:).*9.304250520373691e+3-in1(21,:).*1.067538491807751e+4-in1(9,:).*3.644846607421831e+1-in1(10,:).*4.629304537655445e+2-in1(11,:).*9.031418622604045e+1-in1(12,:).*3.599923108523556e+3-in1(13,:).*6.024411273628679e+3-in1(14,:).*2.661921881923726e+3;in1(15,:).*7.4e+4+in1(16,:).*7.485648263916063e+4+in1(17,:).*7.425333523708081e+4+in1(18,:).*9.191283915080885e+4+in1(19,:).*7.871279414291678e+4+in1(20,:).*8.489206338084523e+4+in1(21,:).*8.437642770638299e+4+in1(1,:).*1.555657856848332e+4+in1(2,:).*1.273703420060655e+4+in1(3,:).*1.401740048428899e+3+in1(4,:).*3.617596327209107e+4+in1(5,:).*6.51510288978945e+3+in1(6,:).*1.213115331247703e+4+in1(7,:).*1.079862406821265e+4+in1(13,:).*(in1(1,:).*6.158968390823141e+4+in1(2,:).*1.524331323587066e+5+in1(3,:).*1.458472400010591e+5+in1(4,:).*6.024898599941494e+4+in1(5,:).*1.477253717241136e+5+in1(6,:).*1.894830041283782e+5+in1(7,:).*1.146760575897481e+5+(in1(8,:).*2.699278917650261e-1+in1(9,:).*8.958862181960668e-1+in1(10,:).*4.280911898712949e-1+in1(11,:).*9.648400471483856e-1+in1(12,:).*6.634414978184481e-1+in1(13,:).*6.216957202091218e-1+in1(14,:).*1.147459729533752e-1).^2.*8.059076851576241e+3+(in1(8,:).*3.155156310060629e-1+in1(9,:).*6.865009276815837e-1+in1(10,:).*8.346256718973729e-1+in1(11,:).*1.828827734419181e-2+in1(12,:).*7.501443149449675e-1+in1(13,:).*9.888610889064947e-1+in1(14,:).*7.481656543798394e-1).^2.*6.4e+4+8.331830760963503e+4)+in1(12,:).*(in1(1,:).*4.126599362184813e+4+in1(2,:).*1.63099844558342e+5+in1(3,:).*4.907797824606964e+4+in1(4,:).*1.033299174307785e+5+in1(5,:).*1.81222445609162e+5+in1(6,:).*1.131903941281831e+5+in1(7,:).*1.168567304904159e+5+(in1(8,:).*1.698304195645689e-1+in1(9,:).*8.781425034294131e-1+in1(10,:).*9.83468338330501e-2+in1(11,:).*4.211076250050522e-1+in1(12,:).*9.578895301505019e-1+in1(13,:).*5.331652849730171e-1+in1(14,:).*6.918771139504734e-1).^2.*6.4e+4+(in1(8,:).*3.477658597455066e-1+in1(9,:).*7.508121031361555e-1+in1(10,:).*7.259979853504515e-1+in1(11,:).*8.833060912058098e-1+in1(12,:).*6.236722070556089e-1+in1(13,:).*7.509424340273372e-1+in1(14,:).*3.488983419778425e-1).^2.*8.346380633892157e+3+8.365050260793781e+4)+in1(9,:).*(in1(1,:).*4.754366839192668e+4+in1(2,:).*7.333568325638418e+4+in1(3,:).*7.28212213777185e+4+in1(4,:).*5.768867708503307e+4+in1(5,:).*8.983338762131439e+4+in1(6,:).*2.844825190556813e+4+in1(7,:).*1.175374092251168e+5+(in1(8,:).*3.455607270430477e-1+in1(9,:).*3.967674742306699e-1+in1(10,:).*5.388167340033569e-1+in1(11,:).*4.191945144032948e-1+in1(12,:).*6.852195003967595e-1+in1(13,:).*2.044522497315174e-1+in1(14,:).*8.781174363909454e-1).^2.*6.4e+4+(in1(8,:).*5.741176054920131e-1+in1(9,:).*1.467285749058101e-1+in1(10,:).*5.893055369032842e-1+in1(11,:).*6.997583600209312e-1+in1(12,:).*1.023344288278258e-1+in1(13,:).*4.140559878195683e-1+in1(14,:).*6.944001577277451e-1).^2.*1.230674790016598e+3+7.542296772595669e+4)+in1(8,:).*(in1(1,:).*5.812225891280424e+4+in1(2,:).*9.220153516059624e+4+in1(3,:).*1.463997662014549e+1+in1(4,:).*3.869856929687549e+4+in1(5,:).*1.878475402459047e+4+in1(6,:).*1.181934013040612e+4+in1(7,:).*2.384130705634188e+4+(in1(8,:).*4.17022004702574e-1+in1(9,:).*7.203244934421581e-1+in1(10,:).*1.143748173448866e-4+in1(11,:).*3.023325726318398e-1+in1(12,:).*1.46755890817113e-1+in1(13,:).*9.23385947687978e-2+in1(14,:).*1.862602113776709e-1).^2.*6.4e+4+7.4e+4)+in1(14,:).*(in1(1,:).*1.137059712620778e+5+in1(2,:).*1.777845142314967e+5+in1(3,:).*5.650263116845519e+4+in1(4,:).*1.106684126247391e+5+in1(5,:).*1.866828580560423e+5+in1(6,:).*1.132410430332701e+5+in1(7,:).*1.256754325898587e+5+(in1(8,:).*9.494892587070712e-1+in1(9,:).*4.499121334799405e-1+in1(10,:).*5.783896143871318e-1+in1(11,:).*4.081368027612812e-1+in1(12,:).*2.370269802430277e-1+in1(13,:).*9.033795205622538e-1+in1(14,:).*5.736794866722859e-1).^2.*3.038365627862628e+4+(in1(8,:).*2.804439920644052e-1+in1(9,:).*7.892793284514885e-1+in1(10,:).*1.03226006577642e-1+in1(11,:).*4.478935261759052e-1+in1(12,:).*9.085955030930956e-1+in1(13,:).*2.936141483736795e-1+in1(14,:).*2.877753385863487e-1).^2.*6.4e+4+1.091311025721616e+5)+in1(10,:).*(in1(1,:).*1.031397421528206e+4+in1(2,:).*9.106219654368827e+4+in1(3,:).*7.986404092141361e+4+in1(4,:).*8.583537628517718e+4+in1(5,:).*2.713807968797873e+4+in1(6,:).*4.178329917937937e+4+in1(7,:).*1.171866142470073e+5+(in1(8,:).*4.141792695269026e-1+in1(9,:).*4.995345894608716e-2+in1(10,:).*5.358964059155116e-1+in1(11,:).*6.637946452197888e-1+in1(12,:).*5.148891120583086e-1+in1(13,:).*9.445947559908133e-1+in1(14,:).*5.865550405019929e-1).^2.*3.313434156215221e+3+(in1(8,:).*2.738759319792616e-2+in1(9,:).*6.704675101784022e-1+in1(10,:).*4.17304802367127e-1+in1(11,:).*5.586898284457517e-1+in1(12,:).*1.403869385952338e-1+in1(13,:).*1.981014890848788e-1+in1(14,:).*8.007445686755367e-1).^2.*6.4e+4+7.783115824312385e+4)+in1(11,:).*(in1(1,:).*1.808376825654413e+5+in1(2,:).*5.357822142815359e+4+in1(3,:).*1.130883963814653e+5+in1(4,:).*1.963058470413962e+5+in1(5,:).*1.530973716556673e+5+in1(6,:).*1.863552467301278e+4+in1(7,:).*3.580831613764168e+4+(in1(8,:).*9.034019152878835e-1+in1(9,:).*1.374747041462375e-1+in1(10,:).*1.392763472507585e-1+in1(11,:).*8.073912887095238e-1+in1(12,:).*3.976768369855336e-1+in1(13,:).*1.653541971169328e-1+in1(14,:).*9.275085803960339e-1).^2.*1.445443064460614e+4+(in1(8,:).*9.682615757193975e-1+in1(9,:).*3.134241781592428e-1+in1(10,:).*6.923226156693141e-1+in1(11,:).*8.763891522960383e-1+in1(12,:).*8.946066635038473e-1+in1(13,:).*8.504421136977791e-2+in1(14,:).*3.905478323288236e-2).^2.*6.4e+4+9.071293543282584e+4);in1(15,:).*7.4e+4+in1(16,:).*7.498340009980316e+4+in1(17,:).*8.020182446915022e+4+in1(18,:).*7.979834729594599e+4+in1(19,:).*9.836845447016621e+4+in1(20,:).*9.769902115934526e+4+in1(21,:).*1.032033351527051e+5+in1(1,:).*2.662105335770335e+4+in1(2,:).*1.460581379454181e+4+in1(3,:).*2.497351328722489e+4+in1(4,:).*1.202230603638233e+4+in1(5,:).*3.259785567664053e+4+in1(6,:).*2.591342864698932e+4+in1(7,:).*2.970079765161742e+4+in1(13,:).*(in1(1,:).*7.338294274542873e+4+in1(2,:).*1.780929803551737e+5+in1(3,:).*1.770436550549501e+5+in1(4,:).*6.093256033682189e+4+in1(5,:).*1.757640652214512e+5+in1(6,:).*2.492793490806106e+5+in1(7,:).*1.426407932418876e+5+in1(8,:).*7.261167808147375e+1+in1(9,:).*1.579889535250644e+2+in1(10,:).*1.920778707955081e+2+in1(11,:).*4.208801012320202+in1(12,:).*1.72635622957092e+2+in1(13,:).*2.275730772070504e+2+in1(14,:).*1.721802608454607e+2+(in1(8,:).*2.699278917650261e-1+in1(9,:).*8.958862181960668e-1+in1(10,:).*4.280911898712949e-1+in1(11,:).*9.648400471483856e-1+in1(12,:).*6.634414978184481e-1+in1(13,:).*6.216957202091218e-1+in1(14,:).*1.147459729533752e-1).^2.*2.674794307286767e+4+(in1(8,:).*3.155156310060629e-1+in1(9,:).*6.865009276815837e-1+in1(10,:).*8.346256718973729e-1+in1(11,:).*1.828827734419181e-2+in1(12,:).*7.501443149449675e-1+in1(13,:).*9.888610889064947e-1+in1(14,:).*7.481656543798394e-1).^2.*6.4e+4+1.049273091780032e+5)+in1(9,:).*(in1(1,:).*4.69104986018903e+4+in1(2,:).*5.701753567515821e+4+in1(3,:).*7.183394923815476e+4+in1(4,:).*5.692058834722905e+4+in1(5,:).*8.857786217168482e+4+in1(6,:).*2.807363471702736e+4+in1(7,:).*1.1592843751822e+5+in1(8,:).*1.014053741302047+in1(9,:).*1.164320798585579+in1(10,:).*1.581166730570597+in1(11,:).*1.230133323602487+in1(12,:).*2.010788100650957+in1(13,:).*5.999685512064525e-1+in1(14,:).*2.576850324672089+(in1(8,:).*3.455607270430477e-1+in1(9,:).*3.967674742306699e-1+in1(10,:).*5.388167340033569e-1+in1(11,:).*4.191945144032948e-1+in1(12,:).*6.852195003967595e-1+in1(13,:).*2.044522497315174e-1+in1(14,:).*8.781174363909454e-1).^2.*6.4e+4+(in1(8,:).*5.741176054920131e-1+in1(9,:).*1.467285749058101e-1+in1(10,:).*5.893055369032842e-1+in1(11,:).*6.997583600209312e-1+in1(12,:).*1.023344288278258e-1+in1(13,:).*4.140559878195683e-1+in1(14,:).*6.944001577277451e-1).^2.*3.145264252206505e+2+7.436367117916138e+4)+in1(8,:).*(in1(1,:).*5.408532322503791e+4+in1(2,:).*9.220153516059624e+4+in1(3,:).*1.463997662014549e+1+in1(4,:).*3.869856929687549e+4+in1(5,:).*1.878475402459047e+4+in1(6,:).*1.181934013040612e+4+in1(7,:).*2.384130705634188e+4+(in1(8,:).*4.17022004702574e-1+in1(9,:).*7.203244934421581e-1+in1(10,:).*1.143748173448866e-4+in1(11,:).*3.023325726318398e-1+in1(12,:).*1.46755890817113e-1+in1(13,:).*9.23385947687978e-2+in1(14,:).*1.862602113776709e-1).^2.*6.4e+4+7.4e+4)+in1(10,:).*(in1(1,:).*1.015436992196411e+4+in1(2,:).*8.715497138588299e+4+in1(3,:).*6.414519247300252e+4+in1(4,:).*8.257954819511461e+4+in1(5,:).*2.63199589439018e+4+in1(6,:).*4.062884037235206e+4+in1(7,:).*1.125201848135028e+5+in1(8,:).*1.602133399035144+in1(9,:).*3.922135045828314e+1+in1(10,:).*2.44117092224358e+1+in1(11,:).*3.2682522607424e+1+in1(12,:).*8.212426754197266+in1(13,:).*1.15886419725817e+1+in1(14,:).*4.684236428881276e+1+(in1(8,:).*4.141792695269026e-1+in1(9,:).*4.995345894608716e-2+in1(10,:).*5.358964059155116e-1+in1(11,:).*6.637946452197888e-1+in1(12,:).*5.148891120583086e-1+in1(13,:).*9.445947559908133e-1+in1(14,:).*5.865550405019929e-1).^2.*3.996276715686973e+2+(in1(8,:).*2.738759319792616e-2+in1(9,:).*6.704675101784022e-1+in1(10,:).*4.17304802367127e-1+in1(11,:).*5.586898284457517e-1+in1(12,:).*1.403869385952338e-1+in1(13,:).*1.981014890848788e-1+in1(14,:).*8.007445686755367e-1).^2.*6.4e+4+7.446206949525131e+4)+in1(11,:).*(in1(1,:).*1.571059101383574e+5+in1(2,:).*4.589629801433371e+4+in1(3,:).*9.61197970141089e+4+in1(4,:).*1.468848128015762e+5+in1(5,:).*1.311708568765877e+5+in1(6,:).*1.655111905259055e+4+in1(7,:).*3.485109625913598e+4+in1(8,:).*6.144933465695334e+1+in1(9,:).*1.989101674201863e+1+in1(10,:).*4.39372636151885e+1+in1(11,:).*5.561878283680822e+1+in1(12,:).*5.67749311038648e+1+in1(13,:).*5.397209118016917+in1(14,:).*2.47855590370721+(in1(8,:).*9.034019152878835e-1+in1(9,:).*1.374747041462375e-1+in1(10,:).*1.392763472507585e-1+in1(11,:).*8.073912887095238e-1+in1(12,:).*3.976768369855336e-1+in1(13,:).*1.653541971169328e-1+in1(14,:).*9.275085803960339e-1).^2.*2.1995952663398e+3+(in1(8,:).*9.682615757193975e-1+in1(9,:).*3.134241781592428e-1+in1(10,:).*6.923226156693141e-1+in1(11,:).*8.763891522960383e-1+in1(12,:).*8.946066635038473e-1+in1(13,:).*8.504421136977791e-2+in1(14,:).*3.905478323288236e-2).^2.*6.4e+4+7.654328202670539e+4)+in1(12,:).*(in1(1,:).*4.455157022755808e+4+in1(2,:).*1.800885823424485e+5+in1(3,:).*5.098061769850685e+4+in1(4,:).*1.114767580542079e+5+in1(5,:).*2.144571138501517e+5+in1(6,:).*1.235051268584877e+5+in1(7,:).*1.302419371303696e+5+in1(8,:).*2.563744124830702e+1+in1(9,:).*1.325635707491928e+2+in1(10,:).*1.484634602456024e+1+in1(11,:).*6.357001309283396e+1+in1(12,:).*1.446021072936531e+2+in1(13,:).*8.048613260320949e+1+in1(14,:).*1.044451218187657e+2+(in1(8,:).*1.698304195645689e-1+in1(9,:).*8.781425034294131e-1+in1(10,:).*9.83468338330501e-2+in1(11,:).*4.211076250050522e-1+in1(12,:).*9.578895301505019e-1+in1(13,:).*5.331652849730171e-1+in1(14,:).*6.918771139504734e-1).^2.*6.4e+4+(in1(8,:).*3.477658597455066e-1+in1(9,:).*7.508121031361555e-1+in1(10,:).*7.259979853504515e-1+in1(11,:).*8.833060912058098e-1+in1(12,:).*6.236722070556089e-1+in1(13,:).*7.509424340273372e-1+in1(14,:).*3.488983419778425e-1).^2.*1.801949047526773e+4+9.483503586202832e+4)+in1(14,:).*(in1(1,:).*1.047393534481414e+5+in1(2,:).*1.525489367653211e+5+in1(3,:).*5.320219266511214e+4+in1(4,:).*9.634794157099178e+4+in1(5,:).*1.576323921725513e+5+in1(6,:).*1.038533366543559e+5+in1(7,:).*9.824983657439013e+4+in1(8,:).*7.392142279423409e+1+in1(9,:).*2.080438611350695e+2+in1(10,:).*2.720904526930947e+1+in1(11,:).*1.180589624028952e+2+in1(12,:).*2.394940673846202e+2+in1(13,:).*7.739290630021854e+1+in1(14,:).*7.585387127319872e+1+(in1(8,:).*9.494892587070712e-1+in1(9,:).*4.499121334799405e-1+in1(10,:).*5.783896143871318e-1+in1(11,:).*4.081368027612812e-1+in1(12,:).*2.370269802430277e-1+in1(13,:).*9.033795205622538e-1+in1(14,:).*5.736794866722859e-1).^2.*1.43971882713581e+4+(in1(8,:).*2.804439920644052e-1+in1(9,:).*7.892793284514885e-1+in1(10,:).*1.03226006577642e-1+in1(11,:).*4.478935261759052e-1+in1(12,:).*9.085955030930956e-1+in1(13,:).*2.936141483736795e-1+in1(14,:).*2.877753385863487e-1).^2.*6.4e+4+9.06467489387578e+4);in1(15,:).*7.4e+4+in1(16,:).*7.533547345588732e+4+in1(17,:).*7.786006942189593e+4+in1(18,:).*8.680796838988231e+4+in1(19,:).*7.672912463886714e+4+in1(20,:).*1.028125050686063e+5+in1(21,:).*7.781936224337276e+4+in1(1,:).*3.477882302701652e+2+in1(2,:).*1.978987200124965e+4+in1(3,:).*1.569396064769172e+4+in1(4,:).*2.599796323137807e+4+in1(5,:).*3.882804588764487e+3+in1(6,:).*3.141633139523973e+4+in1(7,:).*4.215344794495254e+3+in1(12,:).*(in1(1,:).*4.434928900989388e+4+in1(2,:).*1.790426464298374e+5+in1(3,:).*5.08634789824148e+4+in1(4,:).*1.109751861722558e+5+in1(5,:).*2.124109740192406e+5+in1(6,:).*1.228700856018204e+5+in1(7,:).*1.294178577608272e+5+in1(8,:).*2.871243749738971+in1(9,:).*1.484634602456024e+1+in1(10,:).*1.662704082541604+in1(11,:).*7.119470347910682+in1(12,:).*1.619459183717886e+1+in1(13,:).*9.013976977631795+in1(14,:).*1.169724389841984e+1+(in1(8,:).*1.698304195645689e-1+in1(9,:).*8.781425034294131e-1+in1(10,:).*9.83468338330501e-2+in1(11,:).*4.211076250050522e-1+in1(12,:).*9.578895301505019e-1+in1(13,:).*5.331652849730171e-1+in1(14,:).*6.918771139504734e-1).^2.*6.4e+4+(in1(8,:).*3.477658597455066e-1+in1(9,:).*7.508121031361555e-1+in1(10,:).*7.259979853504515e-1+in1(11,:).*8.833060912058098e-1+in1(12,:).*6.236722070556089e-1+in1(13,:).*7.509424340273372e-1+in1(14,:).*3.488983419778425e-1).^2.*1.742395164841083e+4+9.414644409347503e+4)+in1(11,:).*(in1(1,:).*1.571617329156779e+5+in1(2,:).*4.591436772663357e+4+in1(3,:).*9.615971119864089e+4+in1(4,:).*1.470010628371687e+5+in1(5,:).*1.312224332582377e+5+in1(6,:).*1.655602207073435e+4+in1(7,:).*3.48533478681251e+4+in1(8,:).*1.357354252332954e+2+in1(9,:).*4.39372636151885e+1+in1(10,:).*9.705301438475652e+1+in1(11,:).*1.228563202751979e+2+in1(12,:).*1.254101359924508e+2+in1(13,:).*1.192189433452472e+1+in1(14,:).*5.474881728701097+(in1(8,:).*9.034019152878835e-1+in1(9,:).*1.374747041462375e-1+in1(10,:).*1.392763472507585e-1+in1(11,:).*8.073912887095238e-1+in1(12,:).*3.976768369855336e-1+in1(13,:).*1.653541971169328e-1+in1(14,:).*9.275085803960339e-1).^2.*2.228421556012137e+3+(in1(8,:).*9.682615757193975e-1+in1(9,:).*3.134241781592428e-1+in1(10,:).*6.923226156693141e-1+in1(11,:).*8.763891522960383e-1+in1(12,:).*8.946066635038473e-1+in1(13,:).*8.504421136977791e-2+in1(14,:).*3.905478323288236e-2).^2.*6.4e+4+7.657661242413903e+4)+in1(8,:).*(in1(1,:).*7.814273684357669e+4+in1(2,:).*9.220153516059624e+4+in1(3,:).*1.463997662014549e+1+in1(4,:).*3.869856929687549e+4+in1(5,:).*1.878475402459047e+4+in1(6,:).*1.181934013040612e+4+in1(7,:).*2.384130705634188e+4+(in1(8,:).*4.17022004702574e-1+in1(9,:).*7.203244934421581e-1+in1(10,:).*1.143748173448866e-4+in1(11,:).*3.023325726318398e-1+in1(12,:).*1.46755890817113e-1+in1(13,:).*9.23385947687978e-2+in1(14,:).*1.862602113776709e-1).^2.*6.4e+4+7.4e+4)+in1(13,:).*(in1(1,:).*6.456953198173973e+4+in1(2,:).*1.589167048533268e+5+in1(3,:).*1.53729758190645e+5+in1(4,:).*6.042170735663615e+4+in1(5,:).*1.548100162617372e+5+in1(6,:).*2.045919766158416e+5+in1(7,:).*1.217420149141445e+5+in1(8,:).*8.82789347583449e+1+in1(9,:).*1.920778707955081e+2+in1(10,:).*2.335220762347971e+2+in1(11,:).*5.116924436873309+in1(12,:).*2.098848187876232e+2+in1(13,:).*2.766760026255716e+2+in1(14,:).*2.093312042285572e+2+(in1(8,:).*2.699278917650261e-1+in1(9,:).*8.958862181960668e-1+in1(10,:).*4.280911898712949e-1+in1(11,:).*9.648400471483856e-1+in1(12,:).*6.634414978184481e-1+in1(13,:).*6.216957202091218e-1+in1(14,:).*1.147459729533752e-1).^2.*1.27812645669783e+4+(in1(8,:).*3.155156310060629e-1+in1(9,:).*6.865009276815837e-1+in1(10,:).*8.346256718973729e-1+in1(11,:).*1.828827734419181e-2+in1(12,:).*7.501443149449675e-1+in1(13,:).*9.888610889064947e-1+in1(14,:).*7.481656543798394e-1).^2.*6.4e+4+8.877833715556866e+4)+in1(10,:).*(in1(1,:).*1.036731084594811e+4+in1(2,:).*9.236791470977652e+4+in1(3,:).*8.511697238559552e+4+in1(4,:).*8.692341030194295e+4+in1(5,:).*2.741147962641354e+4+in1(6,:).*4.216909671483106e+4+in1(7,:).*1.18746043621054e+5+in1(8,:).*9.971817445296262e-1+in1(9,:).*2.44117092224358e+1+in1(10,:).*1.519405987294103e+1+in1(11,:).*2.034188596837629e+1+in1(12,:).*5.111485749254122+in1(13,:).*7.212871428750974+in1(14,:).*2.915509443067665e+1+(in1(8,:).*4.141792695269026e-1+in1(9,:).*4.995345894608716e-2+in1(10,:).*5.358964059155116e-1+in1(11,:).*6.637946452197888e-1+in1(12,:).*5.148891120583086e-1+in1(13,:).*9.445947559908133e-1+in1(14,:).*5.865550405019929e-1).^2.*4.287171247324093e+3+(in1(8,:).*2.738759319792616e-2+in1(9,:).*6.704675101784022e-1+in1(10,:).*4.17304802367127e-1+in1(11,:).*5.586898284457517e-1+in1(12,:).*1.403869385952338e-1+in1(13,:).*1.981014890848788e-1+in1(14,:).*8.007445686755367e-1).^2.*6.4e+4+7.895704175471848e+4)+in1(14,:).*(in1(1,:).*1.070453206568448e+5+in1(2,:).*1.590388324358213e+5+in1(3,:).*5.405097457162123e+4+in1(4,:).*1.000307724162858e+5+in1(5,:).*1.651033721022154e+5+in1(6,:).*1.062675962474286e+5+in1(7,:).*1.053029447142869e+5+in1(8,:).*9.667823545508282+in1(9,:).*2.720904526930947e+1+in1(10,:).*3.55853876402862+in1(11,:).*1.544035779254528e+1+in1(12,:).*3.132226486110302e+1+in1(13,:).*1.012184199791851e+1+in1(14,:).*9.920559088186083+(in1(8,:).*9.494892587070712e-1+in1(9,:).*4.499121334799405e-1+in1(10,:).*5.783896143871318e-1+in1(11,:).*4.081368027612812e-1+in1(12,:).*2.370269802430277e-1+in1(13,:).*9.033795205622538e-1+in1(14,:).*5.736794866722859e-1).^2.*1.850846766038822e+4+(in1(8,:).*2.804439920644052e-1+in1(9,:).*7.892793284514885e-1+in1(10,:).*1.03226006577642e-1+in1(11,:).*4.478935261759052e-1+in1(12,:).*9.085955030930956e-1+in1(13,:).*2.936141483736795e-1+in1(14,:).*2.877753385863487e-1).^2.*6.4e+4+9.540041573232388e+4)+in1(9,:).*(in1(1,:).*4.75661690614679e+4+in1(2,:).*7.391557391451064e+4+in1(3,:).*7.285630562290698e+4+in1(4,:).*5.771597230843068e+4+in1(5,:).*8.987800466351579e+4+in1(6,:).*2.846156450780604e+4+in1(7,:).*1.175945865277104e+5+in1(8,:).*1.37710160353164+in1(9,:).*1.581166730570597+in1(10,:).*2.147250339339834+in1(11,:).*1.670541218373265+in1(12,:).*2.730683202463505+in1(13,:).*8.147671275036125e-1+in1(14,:).*3.499404981841071+(in1(8,:).*3.455607270430477e-1+in1(9,:).*3.967674742306699e-1+in1(10,:).*5.388167340033569e-1+in1(11,:).*4.191945144032948e-1+in1(12,:).*6.852195003967595e-1+in1(13,:).*2.044522497315174e-1+in1(14,:).*8.781174363909454e-1).^2.*6.4e+4+(in1(8,:).*5.741176054920131e-1+in1(9,:).*1.467285749058101e-1+in1(10,:).*5.893055369032842e-1+in1(11,:).*6.997583600209312e-1+in1(12,:).*1.023344288278258e-1+in1(13,:).*4.140559878195683e-1+in1(14,:).*6.944001577277451e-1).^2.*1.263231541667393e+3+7.546061147005292e+4);in1(15,:).*7.4e+4+in1(16,:).*7.503898619235476e+4+in1(17,:).*7.91678809131232e+4+in1(18,:).*9.021319931747671e+4+in1(19,:).*8.568573659389019e+4+in1(20,:).*7.463133821714077e+4+in1(21,:).*9.057206046850848e+4+in1(1,:).*1.137270808654294e+4+in1(2,:).*1.542428396358817e+4+in1(3,:).*2.087635817883369e+4+in1(4,:).*3.278725176439519e+4+in1(5,:).*1.576808809189636e+4+in1(6,:).*1.089019670454789e+3+in1(7,:).*1.701902570680486e+4+in1(10,:).*(in1(1,:).*1.042335604512885e+4+in1(2,:).*9.373994053477716e+4+in1(3,:).*9.063666094703115e+4+in1(4,:).*8.806669742794967e+4+in1(5,:).*2.769876350269591e+4+in1(6,:).*4.257448602135124e+4+in1(7,:).*1.20384664748511e+5+in1(8,:).*1.335032078759463+in1(9,:).*3.2682522607424e+1+in1(10,:).*2.034188596837629e+1+in1(11,:).*2.723382217858333e+1+in1(12,:).*6.843283566723337+in1(13,:).*9.656629586507696+in1(14,:).*3.903299126537337e+1+(in1(8,:).*4.141792695269026e-1+in1(9,:).*4.995345894608716e-2+in1(10,:).*5.358964059155116e-1+in1(11,:).*6.637946452197888e-1+in1(12,:).*5.148891120583086e-1+in1(13,:).*9.445947559908133e-1+in1(14,:).*5.865550405019929e-1).^2.*5.310357161758311e+3+(in1(8,:).*2.738759319792616e-2+in1(9,:).*6.704675101784022e-1+in1(10,:).*4.17304802367127e-1+in1(11,:).*5.586898284457517e-1+in1(12,:).*1.403869385952338e-1+in1(13,:).*1.981014890848788e-1+in1(14,:).*8.007445686755367e-1).^2.*6.4e+4+8.014010046828305e+4)+in1(14,:).*(in1(1,:).*1.039895524550821e+5+in1(2,:).*1.504387028468915e+5+in1(3,:).*5.29262054290451e+4+in1(4,:).*9.51504439592961e+4+in1(5,:).*1.552031520041544e+5+in1(6,:).*1.030683236526287e+5+in1(7,:).*9.595647027347783e+4+in1(8,:).*4.194830083819231e+1+in1(9,:).*1.180589624028952e+2+in1(10,:).*1.544035779254528e+1+in1(11,:).*6.699509674356226e+1+in1(12,:).*1.359060581879922e+2+in1(13,:).*4.3918268798214e+1+in1(14,:).*4.304491027948459e+1+(in1(8,:).*9.494892587070712e-1+in1(9,:).*4.499121334799405e-1+in1(10,:).*5.783896143871318e-1+in1(11,:).*4.081368027612812e-1+in1(12,:).*2.370269802430277e-1+in1(13,:).*9.033795205622538e-1+in1(14,:).*5.736794866722859e-1).^2.*1.3060377688361e+4+(in1(8,:).*2.804439920644052e-1+in1(9,:).*7.892793284514885e-1+in1(10,:).*1.03226006577642e-1+in1(11,:).*4.478935261759052e-1+in1(12,:).*9.085955030930956e-1+in1(13,:).*2.936141483736795e-1+in1(14,:).*2.877753385863487e-1).^2.*6.4e+4+8.91010617021674e+4)+in1(13,:).*(in1(1,:).*7.468205535135412e+4+in1(2,:).*1.809195978710377e+5+in1(3,:).*1.80480165405114e+5+in1(4,:).*6.100786098535946e+4+in1(5,:).*1.788527297440325e+5+in1(6,:).*2.558663481949428e+5+in1(7,:).*1.457213107738593e+5+in1(8,:).*1.934363747571071+in1(9,:).*4.208801012320202+in1(10,:).*5.116924436873309+in1(11,:).*1.121218007086634e-1+in1(12,:).*4.598985995239639+in1(13,:).*6.062511184200635+in1(14,:).*4.586855219804749+(in1(8,:).*2.699278917650261e-1+in1(9,:).*8.958862181960668e-1+in1(10,:).*4.280911898712949e-1+in1(11,:).*9.648400471483856e-1+in1(12,:).*6.634414978184481e-1+in1(13,:).*6.216957202091218e-1+in1(14,:).*1.147459729533752e-1).^2.*2.880665661707941e+4+(in1(8,:).*3.155156310060629e-1+in1(9,:).*6.865009276815837e-1+in1(10,:).*8.346256718973729e-1+in1(11,:).*1.828827734419181e-2+in1(12,:).*7.501443149449675e-1+in1(13,:).*9.888610889064947e-1+in1(14,:).*7.481656543798394e-1).^2.*6.4e+4+1.073076967134981e+5)+in1(12,:).*(in1(1,:).*4.563164268757339e+4+in1(2,:).*1.856733152563445e+5+in1(3,:).*5.160607518147794e+4+in1(4,:).*1.141548810290939e+5+in1(5,:).*2.253823955259635e+5+in1(6,:).*1.268959042139934e+5+in1(7,:).*1.346420759182778e+5+in1(8,:).*1.229427109281124e+1+in1(9,:).*6.357001309283396e+1+in1(10,:).*7.119470347910682+in1(11,:).*3.048459348057874e+1+in1(12,:).*6.934301635025082e+1+in1(13,:).*3.859661047496907e+1+in1(14,:).*5.008599062304605e+1+(in1(8,:).*1.698304195645689e-1+in1(9,:).*8.781425034294131e-1+in1(10,:).*9.83468338330501e-2+in1(11,:).*4.211076250050522e-1+in1(12,:).*9.578895301505019e-1+in1(13,:).*5.331652849730171e-1+in1(14,:).*6.918771139504734e-1).^2.*6.4e+4+(in1(8,:).*3.477658597455066e-1+in1(9,:).*7.508121031361555e-1+in1(10,:).*7.259979853504515e-1+in1(11,:).*8.833060912058098e-1+in1(12,:).*6.236722070556089e-1+in1(13,:).*7.509424340273372e-1+in1(14,:).*3.488983419778425e-1).^2.*2.119934618893944e+4+9.851174403096122e+4)+in1(9,:).*(in1(1,:).*4.772980309467273e+4+in1(2,:).*7.813277496416012e+4+in1(3,:).*7.311145254338513e+4+in1(4,:).*5.791447429005375e+4+in1(5,:).*9.020247796648623e+4+in1(6,:).*2.855837917337191e+4+in1(7,:).*1.180104031599235e+5+in1(8,:).*1.071372512296264+in1(9,:).*1.230133323602487+in1(10,:).*1.670541218373265+in1(11,:).*1.299665861569746+in1(12,:).*2.124446675107888+in1(13,:).*6.338814086711649e-1+in1(14,:).*2.722505222070922+(in1(8,:).*3.455607270430477e-1+in1(9,:).*3.967674742306699e-1+in1(10,:).*5.388167340033569e-1+in1(11,:).*4.191945144032948e-1+in1(12,:).*6.852195003967595e-1+in1(13,:).*2.044522497315174e-1+in1(14,:).*8.781174363909454e-1).^2.*6.4e+4+(in1(8,:).*5.741176054920131e-1+in1(9,:).*1.467285749058101e-1+in1(10,:).*5.893055369032842e-1+in1(11,:).*6.997583600209312e-1+in1(12,:).*1.023344288278258e-1+in1(13,:).*4.140559878195683e-1+in1(14,:).*6.944001577277451e-1).^2.*1.499997499716282e+3+7.573437210904695e+4)+in1(8,:).*(in1(1,:).*6.109901027361163e+4+in1(2,:).*9.220153516059624e+4+in1(3,:).*1.463997662014549e+1+in1(4,:).*3.869856929687549e+4+in1(5,:).*1.878475402459047e+4+in1(6,:).*1.181934013040612e+4+in1(7,:).*2.384130705634188e+4+(in1(8,:).*4.17022004702574e-1+in1(9,:).*7.203244934421581e-1+in1(10,:).*1.143748173448866e-4+in1(11,:).*3.023325726318398e-1+in1(12,:).*1.46755890817113e-1+in1(13,:).*9.23385947687978e-2+in1(14,:).*1.862602113776709e-1).^2.*6.4e+4+7.4e+4)+in1(11,:).*(in1(1,:).*1.778628537469911e+5+in1(2,:).*5.261527577279755e+4+in1(3,:).*1.109613458814723e+5+in1(4,:).*1.901108138913882e+5+in1(5,:).*1.503488357739335e+5+in1(6,:).*1.837423993635759e+4+in1(7,:).*3.56883265629514e+4+in1(8,:).*1.71823152334474e+2+in1(9,:).*5.561878283680822e+1+in1(10,:).*1.228563202751979e+2+in1(11,:).*1.555199035006239e+2+in1(12,:).*1.58752697491327e+2+in1(13,:).*1.509154638765639e+1+in1(14,:).*6.930478433813341+(in1(8,:).*9.034019152878835e-1+in1(9,:).*1.374747041462375e-1+in1(10,:).*1.392763472507585e-1+in1(11,:).*8.073912887095238e-1+in1(12,:).*3.976768369855336e-1+in1(13,:).*1.653541971169328e-1+in1(14,:).*9.275085803960339e-1).^2.*1.291826061935238e+4+(in1(8,:).*9.682615757193975e-1+in1(9,:).*3.134241781592428e-1+in1(10,:).*6.923226156693141e-1+in1(11,:).*8.763891522960383e-1+in1(12,:).*8.946066635038473e-1+in1(13,:).*8.504421136977791e-2+in1(14,:).*3.905478323288236e-2).^2.*6.4e+4+8.893673884112619e+4);in1(15,:).*7.4e+4+in1(16,:).*7.569833710886666e+4+in1(17,:).*7.529857918200591e+4+in1(18,:).*9.055022327482117e+4+in1(19,:).*1.005814344616764e+5+in1(20,:).*9.989608444154801e+4+in1(21,:).*1.076180336144445e+5+in1(1,:).*5.697270733941707e+3+in1(2,:).*2.513281052708899e+4+in1(3,:).*5.543674628417637e+3+in1(4,:).*3.345920410407609e+4+in1(5,:).*3.553444539582633e+4+in1(6,:).*2.827780991585708e+4+in1(7,:).*3.413314517763083e+4+in1(10,:).*(in1(1,:).*1.03581054184545e+4+in1(2,:).*9.214255938164812e+4+in1(3,:).*8.421036319190318e+4+in1(4,:).*8.67356253195239e+4+in1(5,:).*2.736429323164901e+4+in1(6,:).*4.210251149571497e+4+in1(7,:).*1.184768999977143e+5+in1(8,:).*3.354653278454494e-1+in1(9,:).*8.212426754197266+in1(10,:).*5.111485749254122+in1(11,:).*6.843283566723337+in1(12,:).*1.719572437078377+in1(13,:).*2.426506794600803+in1(14,:).*9.808165226856897+(in1(8,:).*4.141792695269026e-1+in1(9,:).*4.995345894608716e-2+in1(10,:).*5.358964059155116e-1+in1(11,:).*6.637946452197888e-1+in1(12,:).*5.148891120583086e-1+in1(13,:).*9.445947559908133e-1+in1(14,:).*5.865550405019929e-1).^2.*4.119112896466468e+3+(in1(8,:).*2.738759319792616e-2+in1(9,:).*6.704675101784022e-1+in1(10,:).*4.17304802367127e-1+in1(11,:).*5.586898284457517e-1+in1(12,:).*1.403869385952338e-1+in1(13,:).*1.981014890848788e-1+in1(14,:).*8.007445686755367e-1).^2.*6.4e+4+7.876272428653935e+4)+in1(8,:).*(in1(1,:).*6.306595873582828e+4+in1(2,:).*9.220153516059624e+4+in1(3,:).*1.463997662014549e+1+in1(4,:).*3.869856929687549e+4+in1(5,:).*1.878475402459047e+4+in1(6,:).*1.181934013040612e+4+in1(7,:).*2.384130705634188e+4+(in1(8,:).*4.17022004702574e-1+in1(9,:).*7.203244934421581e-1+in1(10,:).*1.143748173448866e-4+in1(11,:).*3.023325726318398e-1+in1(12,:).*1.46755890817113e-1+in1(13,:).*9.23385947687978e-2+in1(14,:).*1.862602113776709e-1).^2.*6.4e+4+7.4e+4)+in1(12,:).*(in1(1,:).*4.351514357660984e+4+in1(2,:).*1.74729528904212e+5+in1(3,:).*5.038043499930438e+4+in1(4,:).*1.089068582294782e+5+in1(5,:).*2.039733236257744e+5+in1(6,:).*1.20251371070077e+5+in1(7,:).*1.260196082806974e+5+in1(8,:).*2.796566212852254e+1+in1(9,:).*1.446021072936531e+2+in1(10,:).*1.619459183717886e+1+in1(11,:).*6.934301635025082e+1+in1(12,:).*1.577339031801275e+2+in1(13,:).*8.779534465286993e+1+in1(14,:).*1.139301289651462e+2+(in1(8,:).*1.698304195645689e-1+in1(9,:).*8.781425034294131e-1+in1(10,:).*9.83468338330501e-2+in1(11,:).*4.211076250050522e-1+in1(12,:).*9.578895301505019e-1+in1(13,:).*5.331652849730171e-1+in1(14,:).*6.918771139504734e-1).^2.*6.4e+4+(in1(8,:).*3.477658597455066e-1+in1(9,:).*7.508121031361555e-1+in1(10,:).*7.259979853504515e-1+in1(11,:).*8.833060912058098e-1+in1(12,:).*6.236722070556089e-1+in1(13,:).*7.509424340273372e-1+in1(14,:).*3.488983419778425e-1).^2.*1.496813296933461e+4+9.130690374579315e+4)+in1(9,:).*(in1(1,:).*4.684472940769155e+4+in1(2,:).*5.532252205618574e+4+in1(3,:).*7.173139839935012e+4+in1(4,:).*5.684080472812175e+4+in1(5,:).*8.844744707858093e+4+in1(6,:).*2.803472213621407e+4+in1(7,:).*1.157613089307283e+5+in1(8,:).*1.751276107845692+in1(9,:).*2.010788100650957+in1(10,:).*2.730683202463505+in1(11,:).*2.124446675107888+in1(12,:).*3.472641552595522+in1(13,:).*1.036148821696117+in1(14,:).*4.450233970142676+(in1(8,:).*3.455607270430477e-1+in1(9,:).*3.967674742306699e-1+in1(10,:).*5.388167340033569e-1+in1(11,:).*4.191945144032948e-1+in1(12,:).*6.852195003967595e-1+in1(13,:).*2.044522497315174e-1+in1(14,:).*8.781174363909454e-1).^2.*6.4e+4+(in1(8,:).*5.741176054920131e-1+in1(9,:).*1.467285749058101e-1+in1(10,:).*5.893055369032842e-1+in1(11,:).*6.997583600209312e-1+in1(12,:).*1.023344288278258e-1+in1(13,:).*4.140559878195683e-1+in1(14,:).*6.944001577277451e-1).^2.*2.19363420498529e+2+7.425363895495142e+4)+in1(13,:).*(in1(1,:).*6.900360847187912e+4+in1(2,:).*1.685643968367179e+5+in1(3,:).*1.654591107365917e+5+in1(4,:).*6.067872036216763e+4+in1(5,:).*1.65352116063586e+5+in1(6,:).*2.270744448639454e+5+in1(7,:).*1.322563077504055e+5+in1(8,:).*7.934328318445599e+1+in1(9,:).*1.72635622957092e+2+in1(10,:).*2.098848187876232e+2+in1(11,:).*4.598985995239639+in1(12,:).*1.886401400149477e+2+in1(13,:).*2.486706764955485e+2+in1(14,:).*1.881425626840108e+2+(in1(8,:).*2.699278917650261e-1+in1(9,:).*8.958862181960668e-1+in1(10,:).*4.280911898712949e-1+in1(11,:).*9.648400471483856e-1+in1(12,:).*6.634414978184481e-1+in1(13,:).*6.216957202091218e-1+in1(14,:).*1.147459729533752e-1).^2.*1.980797902166436e+4+(in1(8,:).*3.155156310060629e-1+in1(9,:).*6.865009276815837e-1+in1(10,:).*8.346256718973729e-1+in1(11,:).*1.828827734419181e-2+in1(12,:).*7.501443149449675e-1+in1(13,:).*9.888610889064947e-1+in1(14,:).*7.481656543798394e-1).^2.*6.4e+4+9.690297574379942e+4)+in1(14,:).*(in1(1,:).*1.009184022657413e+5+in1(2,:).*1.417952823151372e+5+in1(3,:).*5.179577447449992e+4+in1(4,:).*9.024554912594994e+4+in1(5,:).*1.452530966265453e+5+in1(6,:).*9.985294670433628e+4+in1(7,:).*8.656294797293096e+4+in1(8,:).*8.50962009417591e+1+in1(9,:).*2.394940673846202e+2+in1(10,:).*3.132226486110302e+1+in1(11,:).*1.359060581879922e+2+in1(12,:).*2.756986339298449e+2+in1(13,:).*8.909247220960998e+1+in1(14,:).*8.732077966142624e+1+(in1(8,:).*9.494892587070712e-1+in1(9,:).*4.499121334799405e-1+in1(10,:).*5.783896143871318e-1+in1(11,:).*4.081368027612812e-1+in1(12,:).*2.370269802430277e-1+in1(13,:).*9.033795205622538e-1+in1(14,:).*5.736794866722859e-1).^2.*7.584863367776887e+3+(in1(8,:).*2.804439920644052e-1+in1(9,:).*7.892793284514885e-1+in1(10,:).*1.03226006577642e-1+in1(11,:).*4.478935261759052e-1+in1(12,:).*9.085955030930956e-1+in1(13,:).*2.936141483736795e-1+in1(14,:).*2.877753385863487e-1).^2.*6.4e+4+8.276999826899203e+4)+in1(11,:).*(in1(1,:).*1.651681094071102e+5+in1(2,:).*4.850601448282774e+4+in1(3,:).*1.0188439929283e+5+in1(4,:).*1.636742130595426e+5+in1(5,:).*1.386197708572581e+5+in1(6,:).*1.725923697848022e+4+in1(7,:).*3.517628467176869e+4+in1(8,:).*1.753948421428355e+2+in1(9,:).*5.67749311038648e+1+in1(10,:).*1.254101359924508e+2+in1(11,:).*1.58752697491327e+2+in1(12,:).*1.6205269160723e+2+in1(13,:).*1.540525453287934e+1+in1(14,:).*7.074542367298479+(in1(8,:).*9.034019152878835e-1+in1(9,:).*1.374747041462375e-1+in1(10,:).*1.392763472507585e-1+in1(11,:).*8.073912887095238e-1+in1(12,:).*3.976768369855336e-1+in1(13,:).*1.653541971169328e-1+in1(14,:).*9.275085803960339e-1).^2.*6.362829391768537e+3+(in1(8,:).*9.682615757193975e-1+in1(9,:).*3.134241781592428e-1+in1(10,:).*6.923226156693141e-1+in1(11,:).*8.763891522960383e-1+in1(12,:).*8.946066635038473e-1+in1(13,:).*8.504421136977791e-2+in1(14,:).*3.905478323288236e-2).^2.*6.4e+4+8.135702148423237e+4);in1(15,:).*7.4e+4+in1(16,:).*7.450674104065815e+4+in1(17,:).*7.583243877403512e+4+in1(18,:).*7.557331791034089e+4+in1(19,:).*8.879533665800122e+4+in1(20,:).*1.08136938385199e+5+in1(21,:).*8.486372348982614e+4+in1(1,:).*3.712127774099168e+3+in1(2,:).*7.587309654252092e+3+in1(3,:).*7.659172657198806e+3+in1(4,:).*3.598517996957556e+3+in1(5,:).*1.989447857625686e+4+in1(6,:).*3.714624846365734e+4+in1(7,:).*1.128786796083629e+4+in1(12,:).*(in1(1,:).*4.45526326667548e+4+in1(2,:).*1.800940758990112e+5+in1(3,:).*5.098123294477077e+4+in1(4,:).*1.114793924541463e+5+in1(5,:).*2.144678607656757e+5+in1(6,:).*1.235084622780112e+5+in1(7,:).*1.302462654324904e+5+in1(8,:).*1.556580351793818e+1+in1(9,:).*8.048613260320949e+1+in1(10,:).*9.013976977631795+in1(11,:).*3.859661047496907e+1+in1(12,:).*8.779534465286993e+1+in1(13,:).*4.88672529323888e+1+in1(14,:).*6.341398226491821e+1+(in1(8,:).*1.698304195645689e-1+in1(9,:).*8.781425034294131e-1+in1(10,:).*9.83468338330501e-2+in1(11,:).*4.211076250050522e-1+in1(12,:).*9.578895301505019e-1+in1(13,:).*5.331652849730171e-1+in1(14,:).*6.918771139504734e-1).^2.*6.4e+4+(in1(8,:).*3.477658597455066e-1+in1(9,:).*7.508121031361555e-1+in1(10,:).*7.259979853504515e-1+in1(11,:).*8.833060912058098e-1+in1(12,:).*6.236722070556089e-1+in1(13,:).*7.509424340273372e-1+in1(14,:).*3.488983419778425e-1).^2.*1.802261841665609e+4+9.483865254425861e+4)+in1(14,:).*(in1(1,:).*1.128783745236217e+5+in1(2,:).*1.754553325893837e+5+in1(3,:).*5.619800882394126e+4+in1(4,:).*1.093466684193804e+5+in1(5,:).*1.840015716095027e+5+in1(6,:).*1.123745808574542e+5+in1(7,:).*1.231441176299036e+5+in1(8,:).*2.749897890127469e+1+in1(9,:).*7.739290630021854e+1+in1(10,:).*1.012184199791851e+1+in1(11,:).*4.3918268798214e+1+in1(12,:).*8.909247220960998e+1+in1(13,:).*2.879038060972011e+1+in1(14,:).*2.821785521536796e+1+(in1(8,:).*9.494892587070712e-1+in1(9,:).*4.499121334799405e-1+in1(10,:).*5.783896143871318e-1+in1(11,:).*4.081368027612812e-1+in1(12,:).*2.370269802430277e-1+in1(13,:).*9.033795205622538e-1+in1(14,:).*5.736794866722859e-1).^2.*2.890814465799212e+4+(in1(8,:).*2.804439920644052e-1+in1(9,:).*7.892793284514885e-1+in1(10,:).*1.03226006577642e-1+in1(11,:).*4.478935261759052e-1+in1(12,:).*9.085955030930956e-1+in1(13,:).*2.936141483736795e-1+in1(14,:).*2.877753385863487e-1).^2.*6.4e+4+1.074250422608034e+5)+in1(9,:).*(in1(1,:).*4.730653974893007e+4+in1(2,:).*6.722436800880693e+4+in1(3,:).*7.245147773800748e+4+in1(4,:).*5.740101989897878e+4+in1(5,:).*8.93631803391551e+4+in1(6,:).*2.830795388920703e+4+in1(7,:).*1.169348326595885e+5+in1(8,:).*5.225367053079866e-1+in1(9,:).*5.999685512064525e-1+in1(10,:).*8.147671275036125e-1+in1(11,:).*6.338814086711649e-1+in1(12,:).*1.036148821696117+in1(13,:).*3.091607251833457e-1+in1(14,:).*1.327837789935138+(in1(8,:).*3.455607270430477e-1+in1(9,:).*3.967674742306699e-1+in1(10,:).*5.388167340033569e-1+in1(11,:).*4.191945144032948e-1+in1(12,:).*6.852195003967595e-1+in1(13,:).*2.044522497315174e-1+in1(14,:).*8.781174363909454e-1).^2.*6.4e+4+(in1(8,:).*5.741176054920131e-1+in1(9,:).*1.467285749058101e-1+in1(10,:).*5.893055369032842e-1+in1(11,:).*6.997583600209312e-1+in1(12,:).*1.023344288278258e-1+in1(13,:).*4.140559878195683e-1+in1(14,:).*6.944001577277451e-1).^2.*8.875677404601925e+2+7.50262501999071e+4)+in1(10,:).*(in1(1,:).*1.054640307238852e+4+in1(2,:).*9.675221815280066e+4+in1(3,:).*1.027551167436376e+5+in1(4,:).*9.057678007937706e+4+in1(5,:).*2.832949418922945e+4+in1(6,:).*4.346451674261426e+4+in1(7,:).*1.239822513658452e+5+in1(8,:).*4.733786607751187e-1+in1(9,:).*1.15886419725817e+1+in1(10,:).*7.212871428750974+in1(11,:).*9.656629586507696+in1(12,:).*2.426506794600803+in1(13,:).*3.424069319375522+in1(14,:).*1.384040535446825e+1+(in1(8,:).*4.141792695269026e-1+in1(9,:).*4.995345894608716e-2+in1(10,:).*5.358964059155116e-1+in1(11,:).*6.637946452197888e-1+in1(12,:).*5.148891120583086e-1+in1(13,:).*9.445947559908133e-1+in1(14,:).*5.865550405019929e-1).^2.*7.556758047926507e+3+(in1(8,:).*2.738759319792616e-2+in1(9,:).*6.704675101784022e-1+in1(10,:).*4.17304802367127e-1+in1(11,:).*5.586898284457517e-1+in1(12,:).*1.403869385952338e-1+in1(13,:).*1.981014890848788e-1+in1(14,:).*8.007445686755367e-1).^2.*6.4e+4+8.273750149291502e+4)+in1(11,:).*(in1(1,:).*1.579697386757674e+5+in1(2,:).*4.617591744383248e+4+in1(3,:).*9.673744832601309e+4+in1(4,:).*1.486837218319002e+5+in1(5,:).*1.319689746425619e+5+in1(6,:).*1.662699071636918e+4+in1(7,:).*3.48859387408836e+4+in1(8,:).*1.667360264224098e+1+in1(9,:).*5.397209118016917+in1(10,:).*1.192189433452472e+1+in1(11,:).*1.509154638765639e+1+in1(12,:).*1.540525453287934e+1+in1(13,:).*1.464473467666928+in1(14,:).*6.725289459308766e-1+(in1(8,:).*9.034019152878835e-1+in1(9,:).*1.374747041462375e-1+in1(10,:).*1.392763472507585e-1+in1(11,:).*8.073912887095238e-1+in1(12,:).*3.976768369855336e-1+in1(13,:).*1.653541971169328e-1+in1(14,:).*9.275085803960339e-1).^2.*2.645667153870924e+3+(in1(8,:).*9.682615757193975e-1+in1(9,:).*3.134241781592428e-1+in1(10,:).*6.923226156693141e-1+in1(11,:).*8.763891522960383e-1+in1(12,:).*8.946066635038473e-1+in1(13,:).*8.504421136977791e-2+in1(14,:).*3.905478323288236e-2).^2.*6.4e+4+7.705905264666326e+4)+in1(13,:).*(in1(1,:).*6.821710441894342e+4+in1(2,:).*1.668531162886068e+5+in1(3,:).*1.633785910756851e+5+in1(4,:).*6.063313211601148e+4+in1(5,:).*1.634821879649063e+5+in1(6,:).*2.230865676957336e+5+in1(7,:).*1.303913119733001e+5+in1(8,:).*1.045925215242796e+2+in1(9,:).*2.275730772070504e+2+in1(10,:).*2.766760026255716e+2+in1(11,:).*6.062511184200635+in1(12,:).*2.486706764955485e+2+in1(13,:).*3.278045984478902e+2+in1(14,:).*2.480147562259647e+2+(in1(8,:).*2.699278917650261e-1+in1(9,:).*8.958862181960668e-1+in1(10,:).*4.280911898712949e-1+in1(11,:).*9.648400471483856e-1+in1(12,:).*6.634414978184481e-1+in1(13,:).*6.216957202091218e-1+in1(14,:).*1.147459729533752e-1).^2.*1.856160011734854e+4+(in1(8,:).*3.155156310060629e-1+in1(9,:).*6.865009276815837e-1+in1(10,:).*8.346256718973729e-1+in1(11,:).*1.828827734419181e-2+in1(12,:).*7.501443149449675e-1+in1(13,:).*9.888610889064947e-1+in1(14,:).*7.481656543798394e-1).^2.*6.4e+4+9.546185013568425e+4)+in1(8,:).*(in1(1,:).*7.131140545247621e+4+in1(2,:).*9.220153516059624e+4+in1(3,:).*1.463997662014549e+1+in1(4,:).*3.869856929687549e+4+in1(5,:).*1.878475402459047e+4+in1(6,:).*1.181934013040612e+4+in1(7,:).*2.384130705634188e+4+(in1(8,:).*4.17022004702574e-1+in1(9,:).*7.203244934421581e-1+in1(10,:).*1.143748173448866e-4+in1(11,:).*3.023325726318398e-1+in1(12,:).*1.46755890817113e-1+in1(13,:).*9.23385947687978e-2+in1(14,:).*1.862602113776709e-1).^2.*6.4e+4+7.4e+4);in1(15,:).*7.4e+4+in1(16,:).*7.617644043595094e+4+in1(17,:).*8.140688726024871e+4+in1(18,:).*7.472251348980832e+4+in1(19,:).*9.319958991212564e+4+in1(20,:).*9.982777817026809e+4+in1(21,:).*8.46476875276949e+4+in1(1,:).*7.138388347990858e+3+in1(2,:).*3.217258060225795e+4+in1(3,:).*2.974875381336473e+4+in1(4,:).*1.902199356904943e+3+in1(5,:).*2.573885352603942e+4+in1(6,:).*2.820430176082637e+4+in1(7,:).*1.107096833097911e+4+in1(11,:).*(in1(1,:).*1.81584612406816e+5+in1(2,:).*5.382000100254342e+4+in1(3,:).*1.136224632436329e+5+in1(4,:).*1.978613164181036e+5+in1(5,:).*1.537874831193667e+5+in1(6,:).*1.870112890734488e+4+in1(7,:).*3.58384435162102e+4+in1(8,:).*7.657004826260836+in1(9,:).*2.47855590370721+in1(10,:).*5.474881728701097+in1(11,:).*6.930478433813341+in1(12,:).*7.074542367298479+in1(13,:).*6.725289459308766e-1+in1(14,:).*3.088449146405182e-1+(in1(8,:).*9.034019152878835e-1+in1(9,:).*1.374747041462375e-1+in1(10,:).*1.392763472507585e-1+in1(11,:).*8.073912887095238e-1+in1(12,:).*3.976768369855336e-1+in1(13,:).*1.653541971169328e-1+in1(14,:).*9.275085803960339e-1).^2.*1.484013728633654e+4+(in1(8,:).*9.682615757193975e-1+in1(9,:).*3.134241781592428e-1+in1(10,:).*6.923226156693141e-1+in1(11,:).*8.763891522960383e-1+in1(12,:).*8.946066635038473e-1+in1(13,:).*8.504421136977791e-2+in1(14,:).*3.905478323288236e-2).^2.*6.4e+4+9.115890873732663e+4)+in1(13,:).*(in1(1,:).*5.866600599211001e+4+in1(2,:).*1.460717752517399e+5+in1(3,:).*1.381133073354301e+5+in1(4,:).*6.007952044100623e+4+in1(5,:).*1.407742727864306e+5+in1(6,:).*1.746588358798022e+5+in1(7,:).*1.0774329361119e+5+in1(8,:).*7.913399888753188e+1+in1(9,:).*1.721802608454607e+2+in1(10,:).*2.093312042285572e+2+in1(11,:).*4.586855219804749+in1(12,:).*1.881425626840108e+2+in1(13,:).*2.480147562259647e+2+in1(14,:).*1.876462978160536e+2+(in1(8,:).*2.699278917650261e-1+in1(9,:).*8.958862181960668e-1+in1(10,:).*4.280911898712949e-1+in1(11,:).*9.648400471483856e-1+in1(12,:).*6.634414978184481e-1+in1(13,:).*6.216957202091218e-1+in1(14,:).*1.147459729533752e-1).^2.*3.425902408207364e+3+(in1(8,:).*3.155156310060629e-1+in1(9,:).*6.865009276815837e-1+in1(10,:).*8.346256718973729e-1+in1(11,:).*1.828827734419181e-2+in1(12,:).*7.501443149449675e-1+in1(13,:).*9.888610889064947e-1+in1(14,:).*7.481656543798394e-1).^2.*6.4e+4+7.796119965948976e+4)+in1(14,:).*(in1(1,:).*1.069607813275978e+5+in1(2,:).*1.588009056236809e+5+in1(3,:).*5.401985728045488e+4+in1(4,:).*9.989575573078559e+4+in1(5,:).*1.648294776471304e+5+in1(6,:).*1.061790867986141e+5+in1(7,:).*1.050443699049821e+5+in1(8,:).*2.695213431616288e+1+in1(9,:).*7.585387127319872e+1+in1(10,:).*9.920559088186083+in1(11,:).*4.304491027948459e+1+in1(12,:).*8.732077966142624e+1+in1(13,:).*2.821785521536796e+1+in1(14,:).*2.765671505873189e+1+(in1(8,:).*9.494892587070712e-1+in1(9,:).*4.499121334799405e-1+in1(10,:).*5.783896143871318e-1+in1(11,:).*4.081368027612812e-1+in1(12,:).*2.370269802430277e-1+in1(13,:).*9.033795205622538e-1+in1(14,:).*5.736794866722859e-1).^2.*1.835774357351315e+4+(in1(8,:).*2.804439920644052e-1+in1(9,:).*7.892793284514885e-1+in1(10,:).*1.03226006577642e-1+in1(11,:).*4.478935261759052e-1+in1(12,:).*9.085955030930956e-1+in1(13,:).*2.936141483736795e-1+in1(14,:).*2.877753385863487e-1).^2.*6.4e+4+9.522614100687458e+4)+in1(8,:).*(in1(1,:).*5.532548224780056e+4+in1(2,:).*9.220153516059624e+4+in1(3,:).*1.463997662014549e+1+in1(4,:).*3.869856929687549e+4+in1(5,:).*1.878475402459047e+4+in1(6,:).*1.181934013040612e+4+in1(7,:).*2.384130705634188e+4+(in1(8,:).*4.17022004702574e-1+in1(9,:).*7.203244934421581e-1+in1(10,:).*1.143748173448866e-4+in1(11,:).*3.023325726318398e-1+in1(12,:).*1.46755890817113e-1+in1(13,:).*9.23385947687978e-2+in1(14,:).*1.862602113776709e-1).^2.*6.4e+4+7.4e+4)+in1(12,:).*(in1(1,:).*4.127522545861615e+4+in1(2,:).*1.631475796359065e+5+in1(3,:).*4.908332429608192e+4+in1(4,:).*1.033528084821331e+5+in1(5,:).*1.813158286189225e+5+in1(6,:).*1.132193765383654e+5+in1(7,:).*1.168943403402643e+5+in1(8,:).*2.019940817200098e+1+in1(9,:).*1.044451218187657e+2+in1(10,:).*1.169724389841984e+1+in1(11,:).*5.008599062304605e+1+in1(12,:).*1.139301289651462e+2+in1(13,:).*6.341398226491821e+1+in1(14,:).*8.229095980204068e+1+(in1(8,:).*1.698304195645689e-1+in1(9,:).*8.781425034294131e-1+in1(10,:).*9.83468338330501e-2+in1(11,:).*4.211076250050522e-1+in1(12,:).*9.578895301505019e-1+in1(13,:).*5.331652849730171e-1+in1(14,:).*6.918771139504734e-1).^2.*6.4e+4+(in1(8,:).*3.477658597455066e-1+in1(9,:).*7.508121031361555e-1+in1(10,:).*7.259979853504515e-1+in1(11,:).*8.833060912058098e-1+in1(12,:).*6.236722070556089e-1+in1(13,:).*7.509424340273372e-1+in1(14,:).*3.488983419778425e-1).^2.*8.373560207468221e+3+8.368192898988513e+4)+in1(10,:).*(in1(1,:).*1.038950953516438e+4+in1(2,:).*9.291135420738029e+4+in1(3,:).*8.730324083921916e+4+in1(4,:).*8.737624972380562e+4+in1(5,:).*2.752526859638171e+4+in1(6,:).*4.232966552997648e+4+in1(7,:).*1.193950776450803e+5+in1(8,:).*1.913440395090443+in1(9,:).*4.684236428881276e+1+in1(10,:).*2.915509443067665e+1+in1(11,:).*3.903299126537337e+1+in1(12,:).*9.808165226856897+in1(13,:).*1.384040535446825e+1+in1(14,:).*5.594420045530195e+1+(in1(8,:).*4.141792695269026e-1+in1(9,:).*4.995345894608716e-2+in1(10,:).*5.358964059155116e-1+in1(11,:).*6.637946452197888e-1+in1(12,:).*5.148891120583086e-1+in1(13,:).*9.445947559908133e-1+in1(14,:).*5.865550405019929e-1).^2.*4.692440324015943e+3+(in1(8,:).*2.738759319792616e-2+in1(9,:).*6.704675101784022e-1+in1(10,:).*4.17304802367127e-1+in1(11,:).*5.586898284457517e-1+in1(12,:).*1.403869385952338e-1+in1(13,:).*1.981014890848788e-1+in1(14,:).*8.007445686755367e-1).^2.*6.4e+4+7.942563412464343e+4)+in1(9,:).*(in1(1,:).*4.772186500647169e+4+in1(2,:).*7.792819335458922e+4+in1(3,:).*7.309907505247257e+4+in1(4,:).*5.790484471479154e+4+in1(5,:).*9.018673736749582e+4+in1(6,:).*2.855368257487486e+4+in1(7,:).*1.179902313836902e+5+in1(8,:).*2.244282431168015+in1(9,:).*2.576850324672089+in1(10,:).*3.499404981841071+in1(11,:).*2.722505222070922+in1(12,:).*4.450233970142676+in1(13,:).*1.327837789935138+in1(14,:).*5.703030989250677+(in1(8,:).*3.455607270430477e-1+in1(9,:).*3.967674742306699e-1+in1(10,:).*5.388167340033569e-1+in1(11,:).*4.191945144032948e-1+in1(12,:).*6.852195003967595e-1+in1(13,:).*2.044522497315174e-1+in1(14,:).*8.781174363909454e-1).^2.*6.4e+4+(in1(8,:).*5.741176054920131e-1+in1(9,:).*1.467285749058101e-1+in1(10,:).*5.893055369032842e-1+in1(11,:).*6.997583600209312e-1+in1(12,:).*1.023344288278258e-1+in1(13,:).*4.140559878195683e-1+in1(14,:).*6.944001577277451e-1).^2.*1.488511691897548e+3+7.572109164375654e+4)]
[sol,fval] = fsolve(f,zeros(21,1))
Equation solved, inaccuracy possible.
The vector of function values is near zero, as measured by the value
of the function tolerance. However, the last step was ineffective.
sol = 21×1
0.0262
-0.0316
0.0022
0.0029
-0.0010
0.0016
0.0007
-0.0667
0.1455
0.0164
fval = 21×1
1.0e-11 *
-0.0317
-0.0511
-0.0306
-0.0144
-0.0160
-0.0477
-0.0540
0
0.0004
0.0165
U = sol(1:7)
U = 7×1
0.0262
-0.0316
0.0022
0.0029
-0.0010
0.0016
0.0007
W = sol(8:14)
W = 7×1
-0.0667
0.1455
0.0164
0.0220
0.0082
0.0003
-0.0031
PHI = sol(15:21)
PHI = 7×1
-0.3024
-0.1237
0.0826
0.1193
0.1232
0.0419
-0.1014
mahdi
2023-6-10
编辑:mahdi
2023-6-10
Ok thanks, can you explain this line please:
eq1(7) = subs(eq1(7),lhs(eq1(7)),lhs(eq1(7))+sym('100'));
Is this like the line b(7) = -100; ?
Also I must write this whole process after the end of this right?
{}; % for example only
C{1,1} = rand(7);
C{1,2} = rand(7);
Nu = 0.285; A_11 = 64; A_55 = 37; D_11 = 2; w0c = 0.005; R = 0.5;
A11 = 6.4*10^4; A55 = 3.7*10^4; D11 = 2*10^3;
N = 7; a = 0; b = 0.5;
X1 = a+(b-a)*(1-cos(((1:N)-1)*pi/(N-1)))/2;
syms r
w0 = w0c*cos(pi*r/(2*R));
w0onX1 = double(subs(w0,r,X1));
A = C{1,1}*w0onX1.';
B = C{1,2}*w0onX1.';
B1 = C{1,2}*w0onX1.';
C_1 = cell2mat(C(1));
al = zeros(N);
al(2:N,1:N) = C_1(2:N,1:N);
bl = zeros(N);
bl(1:N-1,1:N) = C_1(1:N-1,1:N);
syms U [N,1]
syms W [N,1]
syms PHI [N,1]
z1 = X1;
z2 = B.';
z3 = A.';
z4 = B1.';
Ug1 = 0; Ug2 = 0; Ug3 = 0;
atm=((2*A11*z1.^2).*(C{1,2}))+(((2*A11*z1).*(C{1,1}))-(2*A11));
atm2=(((A11.*z2).*(z1.^2)).*(C{1,2}))+(((((A11-(Nu*A11)).*z3).*z1).*(al)));
atm3=((D11*z1.^2).*C{1,2})+((D11*z1).*C{1,1})-(D11+(A55*z1.^2));
atm4=((A55*z1.^2).*(al));
atm5=((((2*A55*z1)+((2*A11*Ug2).*z1)+((A11*z2.^2).*z1)+(2*Nu*A11*Ug3)).*C{1,2})+...
((2*A55)+(2*A11*Ug2)+((2*A11*Ug1).*z1)+(A11*z3.^2)+(((A11*z4).*z3).*z1).*(al)));
atm6=(((2*A55*z1).*C{1,1})+(2*A55));
atm7=((((2*Nu*A11)+((2*A11*z4).*z1)).*C{1,1})+(2*Nu*A11*z4));
eq1 = (atm*U)+(atm2*W) == 0;
eq2 = (atm3*PHI)-(atm4*W) == 0;
eq3 = (atm5*W)+(atm6*PHI)+(atm7*U) == 0;
vars = [U;W;PHI];
[Omega,b] = equationsToMatrix([eq1,eq2,eq3],vars);
b(7) = -100;
sol = Omega\b;
JabejayiU = double(sol(1:N));
KhizW = double(sol(1+N:N+N));
TaghirzaviyePHI = double(sol(1+(2*N):N+(2*N)));
Then again for example after I wrote all these and the nonlinear system, I must write the iteration procces? Or I can just combine the loop with the nonlinear system?
Torsten
2023-6-10
编辑:Torsten
2023-6-10
Is this like the line b(7) = -100; ?
Yes.
And the code I gave you is complete - no modifications are necessary.
Of course I don't know whether "fsolve" returns a valid and reasonable solution. Nonlinear equations can have none or multiple solutions in contrast to linear systems.
Instead of
f = matlabFunction([lhs(eq1);lhs(eq2);lhs(eq3)],"Vars",{[U;W;PHI]})
[sol,fval] = fsolve(f,zeros(21,1))
you can also try
sol = vpasolve([eq1;eq2;eq3],[U;W;PHI])
to see if you get a different solution.
回答(0 个)
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改,无法完成操作。请重新加载页面以查看其更新后的状态。
您也可以从以下列表中选择网站:
如何获得最佳网站性能
选择中国网站(中文或英文)以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
- América Latina (Español)
- Canada (English)
- United States (English)
欧洲
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom(English)
亚太
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)