Re: ODE system solving

[Geraint Paul Bevan]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Daniele wrote:
| there was an error in the 2nd equation: so i changed the m-file this way:
|
| function SysO=f(y,z)
| g=9810;
| rho=9e-07;
| Cp=1500;
| Ka=3e-5;
| R1=0.53;
| w0=57;
| Va=1000;
| Ta=25+273;
| Etap=0.434*exp(508.3/y(3));
| Nup=Etap/rho;
| dNup=-Etap*508.3/(y(3)^2*rho);
| a=0.41;
| b=0.61;
| Nua=10.2+0.131*(y(3)-273);
| Rz=R1*sqrt(w0/y(1));
| Re=2*y(1)*Rz/Nua;
| Cd=a*Re^b;
| tau=(1/2)*Cd*rho*y(1)^2;
| modified therm
| hc=0.21*Re^(1/3)*Ka/Rz*(1+(8*Va/Rz))^(1/6);
| ^      ^      ^
| SysO(1)=y(2);
| |       |       |
|
SysO(2)=y(1)*y(2)/(3*Nup)-(g/(3*Nup))+(2*tau/(3*Nup*rho*Rz))+((y(2)^2)/y(1))
| -(1/Nup)*dNup*y(2)*(-2*hc*(y(3)-Ta)/(Rz*rho*Cp*y(1)));
| SysO(3)=-2*hc*(y(3)-Ta)/(Rz*rho*Cp*y(1));
| endfunction
|
|
| i changed the initial conditions too (a very small change):
| Y0=[57 0.186 520];
|
| if you run this, you don't find tollerance problems and so on... it works!
|
| now other questions for you:
| how can i plot the 2 functions ( y(:,1) and y(:,3) ) in the same plot
with 2
| y-axes?
| and: it's possible to write an m-file including the commands Y0=[57 0.186
| 520]; z=linspace....  ecc?
| so i can just run an m-file that plot my functions (like matlab does)
|
| thanks all
|


You can put any commands into an m-file (say doplot.m) and then execute
all the commands just by typing the name (doplot) at the Octave prompt.
The only constraint is that the script must not start with the keyword
'function'. See the Octave manual section on functions and scripts for
more information.

There are various plotting interfaces within Octave. Many of the most
matlab-like ones are in the octave-forge package. If you are familiar
with gnuplot, you can do pretty much anything in Octave using the gset,
gplot, graw, etc. commands.

The script 'doplot.m' will plot random data against 2 y axes if you type
'doplot' at the Octave prompt. Have a look at the Gnuplot manual to see
the vast number of options that Gnuplot provides.

doplot.m:
## create random data
x = [1:100]';
y = rand(100,3);
y(:,3) *= 100;
data = [x, y];

## set axis labels and tics
gset ylabel 'left axis'
gset y2label 'right axis'
gset ytics
gset y2tics

## plot first data
gplot data using 1:2 title 'y(1)' axes x1y1 with points

## plot second data over first
hold on
gplot data using 1:4 title 'y(3)' axes x1y2 with lines
hold off

- --
Geraint Bevan
http://homepage.ntlworld.com/geraint.bevan

-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.2.4 (GNU/Linux)

iEYEARECAAYFAkHhbvUACgkQcXV3N50QmNMAEACfU/UfWmDaMfWyF9LAb1L+X1BX
vYkAmwVTt7NZPMcUy8I8/fSMJkaP+RqD
=ZaNv
-----END PGP SIGNATURE-----



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------