From: <tr...@us...> - 2008-03-11 16:50:05
|
Revision: 4715 http://octave.svn.sourceforge.net/octave/?rev=4715&view=rev Author: treichl Date: 2008-03-11 09:50:08 -0700 (Tue, 11 Mar 2008) Log Message: ----------- Changes in help section, tests section and demos section. Modified Paths: -------------- trunk/octave-forge/main/odepkg/inst/odepkg.m trunk/octave-forge/main/odepkg/inst/odepkg_event_handle.m trunk/octave-forge/main/odepkg/inst/odepkg_structure_check.m Modified: trunk/octave-forge/main/odepkg/inst/odepkg.m =================================================================== --- trunk/octave-forge/main/odepkg/inst/odepkg.m 2008-03-11 16:47:41 UTC (rev 4714) +++ trunk/octave-forge/main/odepkg/inst/odepkg.m 2008-03-11 16:50:08 UTC (rev 4715) @@ -17,7 +17,7 @@ %# -*- texinfo -*- %# @deftypefn {Function File} {[@var{}] =} odepkg () %# -%# OdePkg is part of the GNU Octave Repository (resp. the Octave--Forge project). The package includes commands for setting up various options, output functions etc. before solving a set of differential equations with the solver functions that are also included. At this time OdePkg is under development with the main target to make a package that is mostly compatible to proprietary solver products. +%# OdePkg is part of the GNU Octave Repository (the Octave--Forge project). The package includes commands for setting up various options, output functions etc. before solving a set of differential equations with the solver functions that are also included. At this time OdePkg is under development with the main target to make a package that is mostly compatible to proprietary solver products. %# %# If this function is called without any input argument then open the OdePkg tutorial in the Octave window. The tutorial can also be opened with the following command %# @@ -32,7 +32,7 @@ if (nargin == 0) doc ('odepkg'); elseif (nargin > 1) - error ('Number of input arguments must at most be one'); + error ('Number of input arguments must be zero or one'); elseif (ischar (vstr)) feval (vstr); else @@ -45,10 +45,11 @@ %# octave --quiet --eval "odepkg ('odepkg_validate_mfiles')" vfun = {'ode23', 'ode45', 'ode54', 'ode78', ... - 'odeget', 'odeset', ... + 'ode23d', 'ode45d', 'ode54d', 'ode78d', ... + 'odeget', 'odeset', 'odeplot', 'odephas2', 'odephas3', 'odeprint', ... 'odepkg_structure_check', 'odepkg_event_handle', ... - 'odeplot', 'odephas2', 'odephas3', 'odeprint', ... 'odepkg_testsuite_calcscd', 'odepkg_testsuite_calcmescd'}; + %# vfun = sort (vfun); for vcnt=1:length(vfun) printf ('Testing function %s ... ', vfun{vcnt}); @@ -60,11 +61,10 @@ %# From command line in the 'src' directory do something like this: %# octave --quiet --eval "odepkg ('odepkg_validate_ccfiles')" - vfile = {'odepkg_octsolver_mebdfdae.cc', 'odepkg_octsolver_mebdfi.cc', ... - 'odepkg_octsolver_ddaskr.cc', ... + vfile = {'odepkg_octsolver_mebdfdae.cc', 'odepkg_octsolver_mebdfi.cc', ... %#'odepkg_octsolver_ddaskr.cc', ... 'odepkg_octsolver_radau.cc', 'odepkg_octsolver_radau5.cc', ... 'odepkg_octsolver_rodas.cc', 'odepkg_octsolver_seulex.cc'}; - vsolv = {'odebda', 'odebdi', 'odekdi', ... + vsolv = {'odebda', 'odebdi', ...#'odekdi', ... 'ode2r', 'ode5r', 'oders', 'odesx'}; for vcnt=1:length(vfile) @@ -75,16 +75,18 @@ function [] = odepkg_internal_mhelpextract () - vfun = {'odepkg', 'ode23d', 'ode45d', 'ode54d', 'ode78d', ... + vfun = {'odepkg', 'odeget', 'odeset', ... 'ode23', 'ode45', 'ode54', 'ode78', ... - 'odeget', 'odeset', ... + 'ode23d', 'ode45d', 'ode54d', 'ode78d', ... 'odeplot', 'odephas2', 'odephas3', 'odeprint', ... 'odepkg_structure_check', 'odepkg_event_handle', ... 'odepkg_testsuite_calcscd', 'odepkg_testsuite_calcmescd', ... - 'odepkg_testsuite_chemakzo', 'odepkg_testsuite_hires', ... + 'odepkg_testsuite_oregonator', 'odepkg_testsuite_pollution', ... + 'odepkg_testsuite_hires', ... + 'odepkg_testsuite_robertson', 'odepkg_testsuite_chemakzo', ... + 'odepkg_testsuite_transistor', ... 'odepkg_testsuite_implrober', 'odepkg_testsuite_implakzo', ... - 'odepkg_testsuite_oregonator', 'odepkg_testsuite_pollution', ... - 'odepkg_testsuite_robertson', 'odepkg_testsuite_transistor', ... + 'odepkg_testsuite_imptrans', ... }; vfun = sort (vfun); @@ -206,8 +208,8 @@ % test ('ode45', 'verbose'); clear ('all'); % test ('ode23', 'verbose'); clear ('all'); - demo ('odeget.m'); - demo ('odeset.m'); + demo ('odeget.m'); clear ('all'); + demo ('odeset.m'); clear ('all'); demo ('odepkg_structure_check.m'); clear ('all'); demo ('odepkg_event_handle.m', 'verbose'); clear ('all'); Modified: trunk/octave-forge/main/odepkg/inst/odepkg_event_handle.m =================================================================== --- trunk/octave-forge/main/odepkg/inst/odepkg_event_handle.m 2008-03-11 16:47:41 UTC (rev 4714) +++ trunk/octave-forge/main/odepkg/inst/odepkg_event_handle.m 2008-03-11 16:50:08 UTC (rev 4715) @@ -17,7 +17,7 @@ %# -*- texinfo -*- %# @deftypefn {Function File} {[@var{sol}] =} odepkg_event_handle (@var{@@fun}, @var{time}, @var{y}, @var{flag}, [@var{par1}, @var{par2}, @dots{}]) %# -%# Return the solution of the event function that is specified as the first input argument @var{@@fun} in form of a function handle. The second input argument @var{time} is of type double scalar and specifies the time of the event evaluation, the third input argument @var{y} is of type double column vector and specifies the solutions, the third input argument @var{flag} is of type string and can be of the form +%# Return the solution of the event function that is specified as the first input argument @var{@@fun} in form of a function handle. The second input argument @var{time} is of type double scalar and specifies the time of the event evaluation, the third input argument @var{y} either is of type double column vector (for ODEs and DAEs) and specifies the solutions or is of type cell array (for IDEs and DDEs) and specifies the derivatives or the history values, the third input argument @var{flag} is of type string and can be of the form %# @table @option %# @item @code{"init"} %# then initialize internal persistent variables of the function @command{odepkg_event_handle} and return an empty cell array of size 4, @@ -58,89 +58,106 @@ %# a value for veveold if (strcmp (vflag, 'init')) - if (isempty (varargin)) - [veveold, vterm, vdir] = feval (vevefun, vt, vy); + if (~iscell (vy)) + vinpargs = {vevefun, vt, vy}; else - [veveold, vterm, vdir] = feval (vevefun, vt, vy, varargin); + vinpargs = {vevefun, vt, vy{1}, vy{2}}; + vy = vy{1}; %# Delete cell element 2 end + if (nargin > 4) + vinpargs = {vinpargs{:}, varargin{:}}; + end + [veveold, vterm, vdir] = feval (vinpargs{:}); - %# My first implementation was to return row vectors from the event - %# functions but this is faulty if we use LabMat events. So this has - %# changed 20060928 and by now these return vectors must be column - %# vectors. The not so nice implementation then is: + %# We assume that all return values must be column vectors veveold = veveold(:)'; vterm = vterm(:)'; vdir = vdir(:)'; + vtold = vt; vyold = vy; vevecnt = 1; vretcell = cell (1,4); - vtold = vt; %# veveold has been set in the code lines before - for vcnt = 1:4, vretcell{vcnt} = []; end - vyold = vy; vevecnt = 1; - %# Process the event, find the zero crossings either for a rising %# or for a falling edge elseif (isempty (vflag)) - %# Call the event function if an event function has been defined to - %# to get the value(s) veve - if (isempty (varargin)) - [veve, vterm, vdir] = feval (vevefun, vt, vy); + if (~iscell (vy)) + vinpargs = {vevefun, vt, vy}; else - [veve, vterm, vdir] = feval (vevefun, vt, vy, varargin); + vinpargs = {vevefun, vt, vy{1}, vy{2}}; + vy = vy{1}; %# Delete cell element 2 end + if (nargin > 4) + vinpargs = {vinpargs{:}, varargin{:}}; + end + [veve, vterm, vdir] = feval (vinpargs{:}); + + %# We assume that all return values must be column vectors veve = veve(:)'; vterm = vterm(:)'; vdir = vdir(:)'; - %# Check if one or more signs of the event function results changed + %# Check if one or more signs of the event has changed vsignum = (sign (veveold) ~= sign (veve)); if (any (vsignum)) %# One or more values have changed vindex = find (vsignum); %# Get the index of the changed values - if (any (vdir(vindex) == 0)) %# Rising or falling (both are possible) + if (any (vdir(vindex) == 0)) + %# Rising or falling (both are possible) %# Don't change anything, keep the index elseif (any (vdir(vindex) == sign (veve(vindex)))) %# Detected rising or falling, need a new index vindex = find (vdir == sign (veve)); else - %# Found a zero crossing that will not be notified + %# Found a zero crossing but must not be notified vindex = []; end - %# We create a new output values if we found a valid index + %# Create new output values if a valid index has been found if (~isempty (vindex)) %# Change the persistent result cell array vretcell{1} = any (vterm(vindex)); %# Stop integration or not - vretcell{2}(vevecnt,1) = vindex(1,1); %# Take first event that occurs + vretcell{2}(vevecnt,1) = vindex(1,1); %# Take first event found %# Calculate the time stamp when the event function returned 0 and %# calculate new values for the integration results, we do both by %# a linear interpolation vtnew = vt - veve(1,vindex) * (vt - vtold) / (veve(1,vindex) - veveold(1,vindex)); - vynew = (vy - (vt - vtnew) * (vy - vyold) / (vt - vtold)); + vynew = (vy - (vt - vtnew) * (vy - vyold) / (vt - vtold))'; vretcell{3}(vevecnt,1) = vtnew; - vretcell{4}(vevecnt,:) = vynew'; + vretcell{4}(vevecnt,:) = vynew; vevecnt = vevecnt + 1; - end + end %# if (~isempty (vindex)) - end %# Check if one or more signs ... - veveold = veve; vtold = vt; - vretval = vretcell; vyold = vy; + end %# Check for one or more signs ... + veveold = veve; vtold = vt; vretval = vretcell; vyold = vy; - %# Reset this event handling function - elseif (strcmp (vflag, 'done')) + elseif (strcmp (vflag, 'done')) %# Clear this event handling function clear ('veveold', 'vtold', 'vretcell', 'vyold', 'vevecnt'); - for vcnt = 1:4, vretcell{vcnt} = []; end + vretcell = cell (1,4); - end %# if (strcmp (vflag, ...) == true) + end -%!function [veve, vterm, vdir] = feve (vt, vy, varargin) -%! veve = vy(1); %# Which event component should be tread +%!function [veve, vterm, vdir] = feveode (vt, vy, varargin) +%! veve = vy(1); %# Which event component should be treaded %! vterm = 1; %# Terminate if an event is found %! vdir = -1; %# In which direction, -1 for falling -%!test -%! odepkg_event_handle (@feve, 0.0, [0 1 2 3], 'init', 123, 456); -%!test -%! A = odepkg_event_handle (@feve, 2.0, [0 0 3 2], '', 123, 456); -%! if ~(iscell (A) && isempty (A{1})), error; end -%! A = odepkg_event_handle (@feve, 3.0, [-1 0 3 2], '', 123, 456); -%! if ~(iscell (A) && A{1} == 1), error; end -%!test -%! odepkg_event_handle (@feve, 4.0, [0 1 2 3], 'done', 123, 456); +%!function [veve, vterm, vdir] = feveide (vt, vy, vyd, varargin) +%! veve = vy(1); %# Which event component should be treaded +%! vterm = 1; %# Terminate if an event is found +%! vdir = -1; %# In which direction, -1 for falling +%! +%!test %# First call to initialize the odepkg_event_handle function +%! odepkg_event_handle (@feveode, 0.0, [0 1 2 3], 'init', 123, 456); +%!test %# Two calls to find the event that may occur with ODE/DAE syntax +%! A = odepkg_event_handle (@feveode, 2.0, [ 0 0 3 2], '', 123, 456); +%! B = odepkg_event_handle (@feveode, 3.0, [-1 0 3 2], '', 123, 456); +%! assert (A{:}, {[], [], [], []}); +%! assert (B{:}, {1, 1, 2, [0 0 3 2]}); +%!test %# Last call to cleanup the odepkg_event_handle function +%! odepkg_event_handle (@feveode, 4.0, [0 1 2 3], 'done', 123, 456); +%!test %# First call to initialize the odepkg_event_handle function +%! odepkg_event_handle (@feveide, 0.0, {[0 1 2 3], [0 1 2 3]}, 'init', 123, 456); +%!test %# Two calls to find the event that may occur with IDE/DDE syntax +%! A = odepkg_event_handle (@feveide, 2.0, {[0 0 3 2], [0 0 3 2]}, '', 123, 456); +%! B = odepkg_event_handle (@feveide, 3.0, {[-1 0 3 2], [0 0 3 2]}, '', 123, 456); +%! assert (A{:}, {[], [], [], []}); +%! assert (B{:}, {1, 1, 2, [0 0 3 2]}); +%!test %# Last call to cleanup the odepkg_event_handle function +%! odepkg_event_handle (@feveide, 4.0, {[0 1 2 3], [0 1 2 3]}, 'done', 123, 456); %# Local Variables: *** %# mode: octave *** Modified: trunk/octave-forge/main/odepkg/inst/odepkg_structure_check.m =================================================================== --- trunk/octave-forge/main/odepkg/inst/odepkg_structure_check.m 2008-03-11 16:47:41 UTC (rev 4714) +++ trunk/octave-forge/main/odepkg/inst/odepkg_structure_check.m 2008-03-11 16:50:08 UTC (rev 4715) @@ -64,16 +64,17 @@ vfld{vcntarg}); end - switch (vsol) - case {'ode23', 'ode45', 'ode54', 'ode78'} - if (~isscalar (vret.(vfld{vcntarg})) && ... - ~isempty (vret.(vfld{vcntarg}))) + switch (vsol) + case {'ode23', 'ode45', 'ode54', 'ode78', ... + 'ode23d', 'ode45d', 'ode54d', 'ode78d',} + if (~isscalar (vret.(vfld{vcntarg})) && ... + ~isempty (vret.(vfld{vcntarg}))) error ('OdePkg:InvalidParameter', ... - 'Value of option "RelTol" must be a scalar for solver "ode23"'); - end - otherwise + 'Value of option "RelTol" must be a scalar'); + end + otherwise - end + end case 'AbsTol' if (isnumeric (vret.(vfld{vcntarg})) && ... @@ -124,7 +125,7 @@ case 'Refine' if (isscalar (vret.(vfld{vcntarg})) && ... - mod (vret.(vfld{vcntarg}), 1) == 0 && ... + mod (vret.(vfld{vcntarg}), 1) == 0 && ... vret.(vfld{vcntarg}) >= 0 && ... vret.(vfld{vcntarg}) <= 5) else @@ -145,7 +146,7 @@ case 'InitialStep' if (isempty (vret.(vfld{vcntarg})) || ... (isscalar (vret.(vfld{vcntarg})) && ... - isreal (vret.(vfld{vcntarg})) && ... + isreal (vret.(vfld{vcntarg})) && ... vret.(vfld{vcntarg}) > 0) ) else error ('OdePkg:InvalidParameter', ... @@ -174,18 +175,19 @@ case 'Jacobian' if (isempty (vret.(vfld{vcntarg})) || ... - ismatrix (vret.(vfld{vcntarg})) || ... - isa (vret.(vfld{vcntarg}), 'function_handle') || ... - iscell (vret.(vfld{vcntarg}))) + ismatrix (vret.(vfld{vcntarg})) || ... + isa (vret.(vfld{vcntarg}), 'function_handle') || ... + iscell (vret.(vfld{vcntarg}))) else error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or not valid parameter value', ... - vfld{vcntarg}); + 'Unknown parameter name "%s" or not valid parameter value', ... + vfld{vcntarg}); end case 'JPattern' if (isempty (vret.(vfld{vcntarg})) || ... - issparse (vret.(vfld{vcntarg}))) + isvector (vret.(vfld{vcntarg}) || ... + ismatrix (vret.(vfld{vcntarg})))) else error ('OdePkg:InvalidParameter', ... 'Unknown parameter name "%s" or not valid parameter value', ... @@ -223,7 +225,8 @@ case 'MvPattern' if (isempty (vret.(vfld{vcntarg})) || ... - issparse (vret.(vfld{vcntarg}))) + isvector (vret.(vfld{vcntarg}) || ... + ismatrix (vret.(vfld{vcntarg})))) else error ('OdePkg:InvalidParameter', ... 'Unknown parameter name "%s" or not valid parameter value', ... @@ -333,10 +336,10 @@ %!test A = odeset ('Jacobian', [1, 2; 3, 4]); %!test A = odeset ('Jacobian', []); %!error A = odeset ('Jacobian', '12'); -%#!test A = odeset ('JPattern', ); -%#!test A = odeset ('JPattern', ); -%#!error A = odeset ('JPattern', ); -%#!error A = odeset ('JPattern', ); +%!test A = odeset ('JPattern', []); +%!test A = odeset ('JPattern', [1, 2, 4]); +%!test A = odeset ('JPattern', [1, 2; 3, 4]); +%!error A = odeset ('JPattern', 12); %!test A = odeset ('Vectorized', 'on'); %!test A = odeset ('Vectorized', 'off'); %!error A = odeset ('Vectorized', []); @@ -351,10 +354,10 @@ %!error A = odeset ('MStateDependence', [1, 2; 3, 4]); %!error A = odeset ('MStateDependence', []); %!error A = odeset ('MStateDependence', '12'); -%#!test A = odeset ('MvPattern', ); -%#!test A = odeset ('MvPattern', ); -%#!error A = odeset ('MvPattern', ); -%#!error A = odeset ('MvPattern', ); +%!test A = odeset ('MvPattern', []); +%!test A = odeset ('MvPattern', [1, 2, 3 ]); +%!test A = odeset ('MvPattern', [1, 2; 3, 4]); +%!error A = odeset ('MvPattern', 12); %!test A = odeset ('MassSingular', 'yes'); %!test A = odeset ('MassSingular', 'no'); %!test A = odeset ('MassSingular', 'maybe'); @@ -375,18 +378,18 @@ %!error A = odeset ('BDF', []); %!demo +%! # Return the checked OdePkg options structure that is created by +%! # the command odeset. %! %! odepkg_structure_check (odeset); %! -%! %---------------------------------------------------------------- -%! % Returns the checked OdePkg options structure created by odeset. %!demo +%! # Create the OdePkg options structure A with odeset and check it +%! # with odepkg_structure_check. This actually is unnecessary +%! # because odeset automtically calls odepkg_structure_check before +%! # returning. %! %! A = odeset (); odepkg_structure_check (A); -%! -%! %---------------------------------------------------------------- -%! % Create the OdePkg options structure A with odeset and check it -%! % with odepkg_structure_check. %# Local Variables: *** %# mode: octave *** This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |