From help-octave-request at bevo dot che dot wisc dot edu  Mon Nov  3 23:13:05 1997
Subject: Ack!  How to compute double integral?
From: "John W. Eaton" <jwe at bevo dot che dot wisc dot edu>
To: jmr at everest dot radiology dot uiowa dot edu (Joe Reinhardt)
cc: help-octave at bevo dot che dot wisc dot edu
Date: Mon, 3 Nov 1997 23:13:08 -0600

On  3-Nov-1997, Joe Reinhardt <jmr at everest dot radiology dot uiowa dot edu> wrote:

| I am trying to compute an M-dimensional integral using quad().   The
| documentation says quad only does integration over a single variable,
| so I tried nesting subroutine calls over each of my variables.
| 
| This method doesn't seem to work. 

It won't work because Octave's quad function is ultimately implemented
with some Fortran subroutines that can't be called recursively.

| I attached some sample code for 2-D showing that my technique fails.
| The outer quad() call only evaluates one value of "y".
| 
| 1) Is there a way to compute M-dimensional integrals in Octave?

Yes, you can use an M-file specifically designed for multidimensional
integration.  Here is some information that was posted to this list
some time ago:

  Date: Wed, 29 Nov 95 07:49:41 CST
  From: vdp at us0 dot mayo dot EDU (Vinayak Dutt)
  Subject: Re: How to quad("sin(a*sin(x))*sin(x)",0,pi/2)?
  To: help-octave at bevo dot che dot wisc dot edu

  Francesco Potorti wrote:
  # 
  # I found the following in the matlab archives, form the
  # contribution/integration directory or some like that, if I am not
  # wrong.  They should be public domain.  I also have the sources.
  # 
  # ===================================================================
  # Index for integration
  # ===================================================================
  # 
  # nit: <directory>
  # 
  # ===================================================================
  # 
  #  Numerical Integration Toolbox
  #  
  #  MATLAB Toolbox for 1-D, 2-D, and n-D Numerical Integration
  #  
  #  The original 1-D routines were obtained from NETLIB and were 
  #  written by
  #           Howard Wilson
  #           Department of Engineering Mechanics
  #           University of Alabama
  #           Box 870278
  #           Tuscaloosa, Alabama 35487-0278
  #           Phone 205 348-1617
  #           Email address: HWILSON  at  UA1VM.UA.EDU
  #  
  #  The rest of the routines were written by
  #           Bryce Gardner
  #           Ray W. Herrick Laboratories
  #           Purdue University
  #           West Lafayette, IN 47906
  #           Phone: 317-494-0231
  #           Fax:  317-494-0787
  #           Email:  gardner at ecn dot purdue dot edu
  #  
  #  These are the general purpose integration routines: 
  #  
  #         quadg.m      -- High accuracy replacement for QUAD and QUAD8 (1-D)
  #         quad2dg.m    -- 2-D integration over a rectangular region
  #         quad2dggen.m -- 2-D integration over a general region
  #         quadndg.m    -- n-D integration over a n-D hyper-rectangular region

  I had fixed these routines sometime back (to work well in Octave)
  and have used them myself and found them quite useful. Gaussian
  Quadrature works a lot faster too :-)

  I had checked with Gardner some time back and he did not have any
  problems if Octave ftp sites kept his software.

  If anyone is interested, I can mail the fixed version of that MATLAB
  package.

I don't remember whether Vinayak ever made his versions available.
Perhaps something like these functions will eventually be distributed
as part of Octave.  Would someone like to volunteer to make that
happen?

| 2) I used global variables to pass parameters from integrand to
|    integrand.   Is there a better way to do this?

Nope, not with quad.  I don't know about the set of M-files for
quadrature mentioned in the message above.  They may allow for passing
additional parameters.

Thanks,

jwe