octave-nkf: liboctave/Quad.cc annotate

annotate liboctave/Quad.cc @ 8710:739141cde75a ss-3-1-52

fix typo in Array-f.cc

author	Jaroslav Hajek <highegg@gmail.com>
date	Mon, 09 Feb 2009 21:51:31 +0100
parents	4976f66d469b
children	eb63fbe60fab

rev	line source
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	1 /*
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	2
7017 a1dbe9d80eee [project @ 2007-10-12 21:27:11 by jwe] jwe parents: 7016 diff changeset	3 Copyright (C) 1993, 1994, 1995, 1996, 1997, 1998, 2000, 2002, 2003,
a1dbe9d80eee [project @ 2007-10-12 21:27:11 by jwe] jwe parents: 7016 diff changeset	4 2004, 2005, 2006, 2007 John W. Eaton
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	5
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	6 This file is part of Octave.
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	7
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	8 Octave is free software; you can redistribute it and/or modify it
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	9 under the terms of the GNU General Public License as published by the
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5775 diff changeset	10 Free Software Foundation; either version 3 of the License, or (at your
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5775 diff changeset	11 option) any later version.
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	12
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	13 Octave is distributed in the hope that it will be useful, but WITHOUT
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	16 for more details.
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	17
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	18 You should have received a copy of the GNU General Public License
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5775 diff changeset	19 along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5775 diff changeset	20 <http://www.gnu.org/licenses/>.
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	21
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	22 */
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	23
238 780cbbc57b7c [project @ 1993-11-30 20:23:04 by jwe] jwe parents: 3 diff changeset	24 #ifdef HAVE_CONFIG_H
1192 b6360f2d4fa6 [project @ 1995-03-30 21:38:35 by jwe] jwe parents: 1011 diff changeset	25 #include <config.h>
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	26 #endif
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	27
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	28 #include "Quad.h"
1847 2ffe49eb95a5 [project @ 1996-02-03 12:47:55 by jwe] jwe parents: 1546 diff changeset	29 #include "f77-fcn.h"
2292 a265f5e924e3 [project @ 1996-06-06 14:59:05 by jwe] jwe parents: 1993 diff changeset	30 #include "lo-error.h"
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4066 diff changeset	31 #include "quit.h"
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	32 #include "sun-utils.h"
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	33
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	34 static integrand_fcn user_fcn;
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	35 static float_integrand_fcn float_user_fcn;
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	36
5775 ace8d8d26933 [project @ 2006-04-24 19:13:06 by jwe] jwe parents: 5307 diff changeset	37 // FIXME -- would be nice to not have to have this global
260 74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	38 // variable.
74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	39 // Nonzero means an error occurred in the calculation of the integrand
74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	40 // function, and the user wants us to quit.
74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	41 int quad_integration_error = 0;
74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	42
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	43 typedef octave_idx_type (quad_fcn_ptr) (double, int&, double*);
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	44 typedef octave_idx_type (quad_float_fcn_ptr) (float, int&, float*);
3507 00fdd363c098 [project @ 2000-02-01 21:31:44 by jwe] jwe parents: 3136 diff changeset	45
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	46 extern "C"
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	47 {
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	48 F77_RET_T
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	49 F77_FUNC (dqagp, DQAGP) (quad_fcn_ptr, const double&, const double&,
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	50 const octave_idx_type&, const double*, const double&,
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	51 const double&, double&, double&, octave_idx_type&,
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	52 octave_idx_type&, const octave_idx_type&, const octave_idx_type&, octave_idx_type&, octave_idx_type*,
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	53 double*);
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	54
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	55 F77_RET_T
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	56 F77_FUNC (dqagi, DQAGI) (quad_fcn_ptr, const double&, const octave_idx_type&,
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	57 const double&, const double&, double&,
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	58 double&, octave_idx_type&, octave_idx_type&, const octave_idx_type&,
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	59 const octave_idx_type&, octave_idx_type&, octave_idx_type, double);
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	60
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	61 F77_RET_T
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	62 F77_FUNC (qagp, QAGP) (quad_float_fcn_ptr, const float&, const float&,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	63 const octave_idx_type&, const float*, const float&,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	64 const float&, float&, float&, octave_idx_type&,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	65 octave_idx_type&, const octave_idx_type&, const octave_idx_type&, octave_idx_type&, octave_idx_type*,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	66 float*);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	67
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	68 F77_RET_T
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	69 F77_FUNC (qagi, QAGI) (quad_float_fcn_ptr, const float&, const octave_idx_type&,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	70 const float&, const float&, float&,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	71 float&, octave_idx_type&, octave_idx_type&, const octave_idx_type&,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	72 const octave_idx_type&, octave_idx_type&, octave_idx_type, float);
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4192 diff changeset	73 }
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	74
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	75 static octave_idx_type
3136 af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	76 user_function (double x, int& ierr, double result)
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	77 {
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4066 diff changeset	78 BEGIN_INTERRUPT_WITH_EXCEPTIONS;
84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4066 diff changeset	79
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	80 #if defined (sun) && defined (__GNUC__)
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	81 double xx = access_double (x);
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	82 #else
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	83 double xx = *x;
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	84 #endif
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	85
260 74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	86 quad_integration_error = 0;
74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	87
3136 af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	88 double xresult = (*user_fcn) (xx);
af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	89
af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	90 #if defined (sun) && defined (__GNUC__)
af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	91 assign_double (result, xresult);
af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	92 #else
af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	93 *result = xresult;
af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	94 #endif
260 74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	95
74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	96 if (quad_integration_error)
1251 97eac19837dc [project @ 1995-04-11 15:58:32 by jwe] jwe parents: 1192 diff changeset	97 ierr = -1;
260 74d73a4b3fc7 [project @ 1993-12-14 23:57:42 by jwe] jwe parents: 238 diff changeset	98
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4066 diff changeset	99 END_INTERRUPT_WITH_EXCEPTIONS;
84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4066 diff changeset	100
3136 af7ec9d3a5e6 [project @ 1998-02-01 20:11:06 by jwe] jwe parents: 2847 diff changeset	101 return 0;
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	102 }
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	103
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	104 static octave_idx_type
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	105 float_user_function (float x, int& ierr, float result)
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	106 {
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	107 BEGIN_INTERRUPT_WITH_EXCEPTIONS;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	108
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	109 #if defined (sun) && defined (__GNUC__)
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	110 float xx = access_float (x);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	111 #else
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	112 float xx = *x;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	113 #endif
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	114
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	115 quad_integration_error = 0;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	116
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	117 float xresult = (*float_user_fcn) (xx);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	118
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	119 #if defined (sun) && defined (__GNUC__)
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	120 assign_float (result, xresult);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	121 #else
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	122 *result = xresult;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	123 #endif
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	124
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	125 if (quad_integration_error)
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	126 ierr = -1;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	127
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	128 END_INTERRUPT_WITH_EXCEPTIONS;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	129
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	130 return 0;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	131 }
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	132
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	133 double
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	134 DefQuad::do_integrate (octave_idx_type& ier, octave_idx_type& neval, double& abserr)
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	135 {
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	136 octave_idx_type npts = singularities.capacity () + 2;
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	137 double *points = singularities.fortran_vec ();
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	138 double result = 0.0;
1935 23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	139
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	140 octave_idx_type leniw = 183*npts - 122;
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	141 Array<octave_idx_type> iwork (leniw);
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	142 octave_idx_type *piwork = iwork.fortran_vec ();
1935 23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	143
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	144 octave_idx_type lenw = 2*leniw - npts;
1935 23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	145 Array<double> work (lenw);
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	146 double *pwork = work.fortran_vec ();
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	147
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	148 user_fcn = f;
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	149 octave_idx_type last;
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	150
289 c23f50e61c58 [project @ 1994-01-13 06:25:58 by jwe] jwe parents: 260 diff changeset	151 double abs_tol = absolute_tolerance ();
c23f50e61c58 [project @ 1994-01-13 06:25:58 by jwe] jwe parents: 260 diff changeset	152 double rel_tol = relative_tolerance ();
c23f50e61c58 [project @ 1994-01-13 06:25:58 by jwe] jwe parents: 260 diff changeset	153
1935 23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	154 F77_XFCN (dqagp, DQAGP, (user_function, lower_limit, upper_limit,
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	155 npts, points, abs_tol, rel_tol, result,
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	156 abserr, neval, ier, leniw, lenw, last,
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	157 piwork, pwork));
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	158
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	159 return result;
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	160 }
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	161
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	162 float
7924 4976f66d469b miscellaneous cleanup John W. Eaton <jwe@octave.org> parents: 7805 diff changeset	163 DefQuad::do_integrate (octave_idx_type&, octave_idx_type&, float&)
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	164 {
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	165 (*current_liboctave_error_handler) ("incorrect integration function called");
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	166 return 0.0;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	167 }
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	168
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	169 double
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	170 IndefQuad::do_integrate (octave_idx_type& ier, octave_idx_type& neval, double& abserr)
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	171 {
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	172 double result = 0.0;
1935 23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	173
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	174 octave_idx_type leniw = 128;
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	175 Array<octave_idx_type> iwork (leniw);
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	176 octave_idx_type *piwork = iwork.fortran_vec ();
1935 23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	177
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	178 octave_idx_type lenw = 8*leniw;
1935 23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	179 Array<double> work (lenw);
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	180 double *pwork = work.fortran_vec ();
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	181
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	182 user_fcn = f;
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	183 octave_idx_type last;
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	184
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	185 octave_idx_type inf;
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	186 switch (type)
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	187 {
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	188 case bound_to_inf:
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	189 inf = 1;
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	190 break;
1360 7eb93d12654c [project @ 1995-09-05 21:51:54 by jwe] jwe parents: 1315 diff changeset	191
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	192 case neg_inf_to_bound:
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	193 inf = -1;
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	194 break;
1360 7eb93d12654c [project @ 1995-09-05 21:51:54 by jwe] jwe parents: 1315 diff changeset	195
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	196 case doubly_infinite:
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	197 inf = 2;
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	198 break;
1360 7eb93d12654c [project @ 1995-09-05 21:51:54 by jwe] jwe parents: 1315 diff changeset	199
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	200 default:
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	201 assert (0);
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	202 break;
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	203 }
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	204
289 c23f50e61c58 [project @ 1994-01-13 06:25:58 by jwe] jwe parents: 260 diff changeset	205 double abs_tol = absolute_tolerance ();
c23f50e61c58 [project @ 1994-01-13 06:25:58 by jwe] jwe parents: 260 diff changeset	206 double rel_tol = relative_tolerance ();
c23f50e61c58 [project @ 1994-01-13 06:25:58 by jwe] jwe parents: 260 diff changeset	207
1935 23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	208 F77_XFCN (dqagi, DQAGI, (user_function, bound, inf, abs_tol, rel_tol,
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	209 result, abserr, neval, ier, leniw, lenw,
23c350d0cf9d [project @ 1996-02-12 04:03:05 by jwe] jwe parents: 1882 diff changeset	210 last, piwork, pwork));
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	211
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	212 return result;
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	213 }
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	214
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	215 float
7924 4976f66d469b miscellaneous cleanup John W. Eaton <jwe@octave.org> parents: 7805 diff changeset	216 IndefQuad::do_integrate (octave_idx_type&, octave_idx_type&, float&)
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	217 {
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	218 (*current_liboctave_error_handler) ("incorrect integration function called");
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	219 return 0.0;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	220 }
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	221
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	222 double
7924 4976f66d469b miscellaneous cleanup John W. Eaton <jwe@octave.org> parents: 7805 diff changeset	223 FloatDefQuad::do_integrate (octave_idx_type&, octave_idx_type&, double&)
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	224 {
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	225 (*current_liboctave_error_handler) ("incorrect integration function called");
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	226 return 0.0;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	227 }
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	228
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	229 float
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	230 FloatDefQuad::do_integrate (octave_idx_type& ier, octave_idx_type& neval, float& abserr)
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	231 {
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	232 octave_idx_type npts = singularities.capacity () + 2;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	233 float *points = singularities.fortran_vec ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	234 float result = 0.0;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	235
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	236 octave_idx_type leniw = 183*npts - 122;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	237 Array<octave_idx_type> iwork (leniw);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	238 octave_idx_type *piwork = iwork.fortran_vec ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	239
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	240 octave_idx_type lenw = 2*leniw - npts;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	241 Array<float> work (lenw);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	242 float *pwork = work.fortran_vec ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	243
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	244 float_user_fcn = ff;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	245 octave_idx_type last;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	246
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	247 float abs_tol = single_precision_absolute_tolerance ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	248 float rel_tol = single_precision_relative_tolerance ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	249
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	250 F77_XFCN (qagp, QAGP, (float_user_function, lower_limit, upper_limit,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	251 npts, points, abs_tol, rel_tol, result,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	252 abserr, neval, ier, leniw, lenw, last,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	253 piwork, pwork));
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	254
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	255 return result;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	256 }
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	257
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	258 double
7924 4976f66d469b miscellaneous cleanup John W. Eaton <jwe@octave.org> parents: 7805 diff changeset	259 FloatIndefQuad::do_integrate (octave_idx_type&, octave_idx_type&, double&)
7805 62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	260 {
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	261 (*current_liboctave_error_handler) ("incorrect integration function called");
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	262 return 0.0;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	263 }
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	264
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	265 float
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	266 FloatIndefQuad::do_integrate (octave_idx_type& ier, octave_idx_type& neval, float& abserr)
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	267 {
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	268 float result = 0.0;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	269
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	270 octave_idx_type leniw = 128;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	271 Array<octave_idx_type> iwork (leniw);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	272 octave_idx_type *piwork = iwork.fortran_vec ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	273
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	274 octave_idx_type lenw = 8*leniw;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	275 Array<float> work (lenw);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	276 float *pwork = work.fortran_vec ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	277
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	278 float_user_fcn = ff;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	279 octave_idx_type last;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	280
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	281 octave_idx_type inf;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	282 switch (type)
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	283 {
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	284 case bound_to_inf:
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	285 inf = 1;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	286 break;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	287
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	288 case neg_inf_to_bound:
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	289 inf = -1;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	290 break;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	291
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	292 case doubly_infinite:
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	293 inf = 2;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	294 break;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	295
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	296 default:
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	297 assert (0);
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	298 break;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	299 }
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	300
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	301 float abs_tol = single_precision_absolute_tolerance ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	302 float rel_tol = single_precision_relative_tolerance ();
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	303
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	304 F77_XFCN (qagi, QAGI, (float_user_function, bound, inf, abs_tol, rel_tol,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	305 result, abserr, neval, ier, leniw, lenw,
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	306 last, piwork, pwork));
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	307
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	308 return result;
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	309 }
62affb34e648 Make quad work with single precision David Bateman <dbateman@free.fr> parents: 7482 diff changeset	310
3 9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	311 /*
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	312 ;;; Local Variables: ***
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	313 ;;; mode: C++ ***
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	314 ;;; End: ***
9a4c07481e61 [project @ 1993-08-08 01:20:23 by jwe] jwe parents: diff changeset	315 */

Mercurial > octave-nkf

annotate liboctave/Quad.cc @ 8710:739141cde75a ss-3-1-52