None.
Open in the MATLAB editor: d01rc_example
function d01rc_example fprintf('d01rc example results\n\n'); % Setup phase. % set problem parameters ni = int64(2); nx = int64(0); % lower (a) and upper (b) bounds a = 0; b = pi; iopts = zeros(100, 1, 'int64'); opts = zeros(100, 1); % initialize option arrays [iopts, opts, ifail] = d01zk('Initialize = d01ra', iopts, opts); % set any non-default options required [iopts, opts, ifail] = d01zk('Quadrature Rule = gk41', iopts, opts); [iopts, opts, ifail] = d01zk('Absolute Tolerance = 1.0e-7', iopts, opts); [iopts, opts, ifail] = d01zk('Relative Tolerance = 1.0e-7', iopts, opts); % determine maximum required array lengths [lenxrq, ldfmrq, sdfmrq, licmin, licmax, lcmin, lcmax, ifail] = ... d01rc(ni, iopts, opts); % allocate remaining arrays needi = zeros(ni, 1, 'int64'); comm = zeros(lcmax, 1); icomm = zeros(licmax, 1, 'int64'); fm = zeros(ldfmrq, sdfmrq); dinest = zeros(ni, 1); errest = zeros(ni, 1); x = zeros(1, lenxrq); % Solve phase. % Use d01ra to evaluate the definate integrals of: % f_1 = (x*sin(2*x))*cos(15*x) % f_2 = (x*sin(2*x))*(x*cos(50*x)) % set initial irevcm irevcm = int64(1); while irevcm ~= 0 [irevcm, sid, needi, x, nx, dinest, errest, icomm, comm, ifail] = ... d01ra(irevcm, a, b, needi, x, nx, fm, dinest, errest, ... iopts, opts, icomm, comm); switch irevcm case 11 % Initial returns. % These will occur during the non-adaptive phase. % All values must be supplied. % dinest and errest do not contain approximations % over the complete interval at this stage. % Calculate x*sin(2*x), storing the result in fm(2,1:nx) for re-use. fm(2, :) = x.*sin(2*x); % Calculate f_1 fm(1, :) = fm(2, :).*cos(15*x); % Calculate f_2 fm(2, :) = fm(2, :).*x.*cos(50*x); case 12 % Intermediate returns. % These will occur during the adaptive phase. % All requested values must be supplied. % dinest and errest do not contain approximations % over the complete interval at this stage. % Calculate x*sin(2*x). fm(2, :) = x.*sin(2*x); % Calculate f_1 if required if needi(1) == 1 fm(1, :) = fm(2, :).*cos(15*x); end % Complete f_2 calculation if required. if needi(2) == 1 fm(2, :) = fm(2, :).*x.*cos(50*x); end case 0 % Final return end end % query some currently set options and statistics. [ivalue, rvalue, cvalue, optype, ifail] = ... d01zl('Quadrature rule', iopts, opts); display_option('Quadrature rule',optype,ivalue,rvalue,cvalue); [ivalue, rvalue, cvalue, optype, ifail] = ... d01zl('Maximum Subdivisions', iopts, opts); display_option('Maximum Subdivisions',optype,ivalue,rvalue,cvalue); [ivalue, rvalue, cvalue, optype, ifail] = ... d01zl('Extrapolation', iopts, opts); display_option('Extrapolation',optype,ivalue,rvalue,cvalue); [ivalue, rvalue, cvalue, optype, ifail] = ... d01zl('Extrapolation Safeguard', iopts, opts); display_option('Extrapolation Safeguard',optype,ivalue,rvalue,cvalue); % print solution fprintf('\nIntegral | needi | dinest | errest \n'); for j=1:ni fprintf('%9d %9d %12.4e %12.4e\n', j, needi(j), dinest(j), errest(j)); end function [dinest, errest, user] = monit(ni, ns, dinest, errest, fcount, ... sinfoi, evals, ldi, sinfor, fs, ... es, ldr, user) % Display information on individual segments fprintf('\nInformation on splitting and evaluations over subregions.\n'); for k=1:ns sid = sinfoi(1,k); parent = sinfoi(2,k); child1 = sinfoi(3,k); child2 = sinfoi(4,k); level = sinfoi(5,k); lbnd = sinfor(1,k); ubnd = sinfor(2,k); fprintf('\nSegment %3d Sid = %3d', k, sid); fprintf(' Parent = %3d Level = %3d.\n', parent, level); if (child1>0) fprintf('Children = (%3d, %3d)\n', child1, child2); end fprintf('Bounds (%11.4e, %11.4e)\n', lbnd, ubnd); for j = 1:ni if (evals(j,k) ~= 0) fprintf('Integral %2d approximation %11.4e\n', j, fs(j,k)); fprintf('Integral %2d error estimate %11.4e\n', j, es(j,k)); end if (evals(j,k) ~= 1) fprintf('Integral %2d evaluation', j); fprintf(' has been superseded by descendants.\n'); end end end function display_option(optstr,optype,ivalue,rvalue,cvalue) % Query optype and print the appropriate option values switch optype case 1 fprintf('%30s: %13d\n', optstr, ivalue); case 2 fprintf('%30s: %13.4e\n', optstr, rvalue); case 3 fprintf('%30s: %16s\n', optstr, cvalue); case 4 fprintf('%30s: %3d %16s\n', optstr, ivalue, cvalue); case 5 fprintf('%30s: %14.4e %16s\n', optstr, rvalue, cvalue); end
d01rc example results Quadrature rule: GK41 Maximum Subdivisions: 50 Extrapolation: ON Extrapolation Safeguard: 1.0000e-12 Integral | needi | dinest | errest 1 0 -2.8431e-02 1.1234e-14 2 0 7.9083e-03 2.6600e-09