How to extract variable from function file, and how to check is that result real?

10 visualizzazioni (ultimi 30 giorni)
I am using matlab function file to solve system of differential equations, and I have variable R in my function file.
function f=fun_z(z,p)
ri=2;
R=ri-z*(ri-1);
sig=1; beta=1;
f=zeros(4,1);
f(1)=-32*1*beta/(R.^4*p(1));
f(2)=(-(2-sig)*8*f(1)/(sig.*R)-f(1)*p(2))/p(1);
f(3)=(-p(2)*f(2)+(2-sig)*(-8*f(2)/R-8*f(1)/(R.*R*p(1)))/sig-f(1)*p(3))/p(1);
f(4)=(-f(2)*p(3)-f(3)*p(2)+(2-sig)*8*(-f(3)/R-f(2)/p(1)+p(2)*f(1)/(p(1).*p(1)))/(sig.*R.*R) -f(1)*p(4))/p(1);
I am not sure did I properly tell to Matlab what is necessary to do with variable R so I extracted it from function file into workspace with command:
[zv,pv,R]=ode45(@fun_z,[1 0],[1; 0; 0; 0])
where I previous in function file changed the first line in:
function f=fun_z(z,p,R)
and added line of code in function file:
global R;
As a result I need to get R where it is dependent on z, so it means R need to have the same number of elements as z, but that is not the case. I cannot see why, but I only get 6 variables? What would be correct way to tell Matlab exactly what R (R is radial coordinate, z is longitudinal coordinate, R is linearly dependent on z, where z is going from 1 to 0, with some step) is? Or how to extract it properly?

Risposta accettata

Jan
Jan il 12 Set 2018
[zv,pv,R]=ode45(@fun_z,[1 0],[1; 0; 0; 0])
This is a bold and wrong guess. Why do you assume, that ode45 replies an arbitrary global variable from the function to be integrated as 3rd output?
You can either modify the function to accept vectors, then:
function [f, R] = fun_z(z, p)
ri = 2;
R = ri - z * (ri - 1);
sig = 1;
beta = 1;
f = zeros(4, size(p, 2));
f(1, :) = -32 * beta / (R .^ 4 * p(1, :));
f(2, :) = ((-2+sig) * 8 * f(1, :) / (sig .* R) - f(1, :) *. p(2, :)) / p(1, :);
... etc. I cannot test this currently, please debug this by your own.
Then:
[zv, pv] = ode45(@fun_z, [1 0], [1; 0; 0; 0]);
[~, R] = fun_z(zv.', pv.');
Or in your case simply:
[zv, pv] = ode45(@fun_z, [1 0], [1; 0; 0; 0]);
ri = 2;
R = ri - zv * (ri - 1);

Più risposte (0)

Categorie

Scopri di più su Programming in Help Center e File Exchange

Prodotti


Release

R2015a

Community Treasure Hunt

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

Start Hunting!

Translated by