I want help in writing a matlab code so that I wouldn't have to manually input the various steps as I have shown below. Thanks

2 views (last 30 days)
George Adomako Kumi
George Adomako Kumi on 1 Dec 2021
E = 2900; % elastic modulus
A = 15; % area
Ix = 19; % modulus of inertia about x-axis
Iy = 11; % modulus of inertia about y-axis
Ixy = 12; % modulus of inertia about x-y axis
K = [ E*A 1 1 1; 1 E*Ix 1 1; 1 1 E*Iy 1; 1 1 1 E*Ixy ]; % stiffness matrix
u = [ 4; 6; 7; 8 ]; % deformation vector
F = [ 15;29;34;45]; % force vector
Q = [ 100;78;85;79]; % maximum force vector
R = Q - F; % initialize R
dF = 2*ones(4,1); % incremental force vector
T = F + dF; % new force vector
R = Q - T; % residual force vector
m = pinv(K)*R; % current deformation vector
q = u + m; %% end of first step
T1 = T + dF; % new force vector (2)
R1 = Q - T1; % residual force vector (2)
K1 = T1\q; % updated stiffness matrix (2)
m1 = pinv(K1)*R1; % updated current deformation vector (2)
q1 = q + m1; %% end of second step
T2 = T1 + dF; % new force vector (3)
R2 = Q - T2; % residual force vector (3)
K2 = T2\q1; % updated stiffness matrix (3)
m2 = pinv(K2)*R2; % updated current deformation vector (3)
q2 = q1 + m2; %% end of third step
T3 = T2 + dF; % new force vector (4)
R3 = Q - T3; %residual force vector (4)
K3 = T3\q2; % updated stiffness matrix (4)
m3 = pinv(K3)*R3; % updated current deformation vector (4)
q3 = q2 + m3; %% end of fourth step
T4 = T3 + dF; % new force vector (5)
R4 = Q - T4; % residual force vector (5)
K4 = T4\q3; % updated stiffness matrix (5)
m4 = pinv(K4)*R4; % updated current deformation vector (5)
q4 = q3 + m4; %% end of fifth step

Answers (2)

Walter Roberson
Walter Roberson on 1 Dec 2021
Unfortunately, someone who was not familiar with the field would have little idea what is going on. What is K -- is that spring constants? Are Ix an Iy and Ixy momentums ? Is Q "quality factor" ?
You are not updating K during your loop, so it is not clear why you are calculating pinv(K) inside the loop instead of before the loop?
In order to "write clearly" then you need comments -- even if you were using variable names that explained more, you would at least have comments at the beginning explaining the purpose of the code.
  6 Comments

Sign in to comment.


John D'Errico
John D'Errico on 1 Dec 2021
Edited: John D'Errico on 1 Dec 2021
R is a matrix. So the statemen
while R
is meaningless. As well, this if statement
if R <= 0.0001;
break
end
is also useless.
The problem is, there is no indication of what you are trying to do. What problem you are trying to solve. (I recall you have asked essentially the same question before.) So it is difficult to know if your implementation is correct.
If your goal is to iterate until the MATRIX R is essentially zero, then a good idea is to test the matrix norm, applied to R.
Next, you can put that test into the while statement. So there is absolutely no need to have an if statement ainside the loop, since the loop will iterate until that test in the while loop is not satisfied.
What else? It is often a bad idea to use i as a loop counter, since i is defined in MATLAB as sqrt(-1).
What else? This is illegal syntax:
dF = 2*ones*(4,1);
What else? R was never defined initially in your code, you you use it in the call to the while loop before it is defined inside the loopitself. You NEED to initialize R.
Given all of that, I might try this code:
E = 2900;
A = 15;
Ix = 19;
Iy = 11;
Ixy = 12;
K = [ E*A 1 1 1; 1 E*Ix 1 1; 1 1 E*Iy 1; 1 1 1 E*Ixy ];
u = [ 4; 6; 7; 8 ];
F = [ 15;29;34;45];
Q = [ 100;78;85;79];
R = Q - F; % initialize R
dF = 2*ones(4,1);
tol = 1e-8;
iter = 0;
while norm(R) > tol
iter = iter + 1
T = F + dF;
R = Q - T;
m = pinv(K)*R;
q = u + m;
S = F/q;
end
However, this code seems never to converge.

Community Treasure Hunt

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

Start Hunting!

Translated by