Trying to code from Polymath

Im supposed to turn this Polymath code into a Matlab code but Im not sure how. Please help. I am supposed to graph as well<</matlabcentral/answers/uploaded_files/1577366/IMG_0357.jpeg>>

Answer 1

[Alan Stevens]

Like this for example (I'll leave you to extract the numerical values):

Wspan = [0 2];      % Integration range
% Initial values
Fe0 = 5.58E-04;     
Fo0 = 0.001116;
Fd0 = 0;
Fu10 = 1E-7;
Fu20 = 0;
F0 = [Fe0, Fo0, Fd0, Fu10, Fu20];

% Call integrating function 
[W, F] = ode45(@fn, Wspan, F0);

% Extract individual variables
Fe = F(:,1); Fo = F(:,2); Fd = F(:,3); Fu1 = F(:,4); Fu2 = F(:,5);
% Note that Fu1 and Fu2 are virtually indistinguishable!

% Plot F vs W 
plot(W,F), grid
xlabel('W'), ylabel('F')
legend('Fe','Fo','Fd','Fu1','Fu2')

% Calculate X and S
X = 1 - Fe/0.000558;
S = Fd./Fu1;

% Plot X vs W and S vs W
figure
plot(W,X,W,S), grid
xlabel('W'), ylabel('X and S')
legend('X','S')

function dbydW = fn(~,F)
        
        Fe = F(1); Fo = F(2); Fd = F(3); Fu1 = F(4); Fu2 = F(5);
        Finert = 0.007626;
        Ft = Fe + Fo + Fd + Fu1 + Fu2 + Finert;
        
        K1 = 6.5;
        K2 = 4.33;
        Pto = 2;
        Pe = Pto*Fe/Ft;
        Po = Pto*Fo/Ft;
        k1 = 0.15;
        k2 = 0.088;
        r1e = -k1*Pe*Po^0.58/(1+K1*Pe)^2;
        r2e = -k2*Pe*Po^0.3/(1+K2*Pe)^2;
        
        dbydW = [ r1e + r2e;
                  r1e/2 + 3*r2e;
                  -r1e;
                  -2*r2e;
                  -2*r2e];

end

[Houtan]

Thanks Alan , I tried to run the script on matlab, but it gives me unrecognized function fn

[abdurahman]

also can you plot X Vs S?

[Alan Stevens]

The function is listed below the second figure.

plot(X,S) (or plot(S,X) if you want S on the x-axis and X on the y-axis)) will plot X vs S.