help-octave
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: clobbered


From: Andy Jacobson
Subject: Re: clobbered
Date: 02 Nov 2000 20:52:09 +0100
User-agent: Gnus/5.0806 (Gnus v5.8.6) Emacs/20.7

>>>>> "John" == John W Eaton <address@hidden> writes:

    John> Can you please send a *complete* example that demonstrates
    John> the problem?

Sorry, since the code is copyrighted, I was reluctant to send it the
way it was.  Here's a fragment which does manifest the problem.  Note
that the problem goes if either the F77_XFCN call is removed *or* if
optimization in the mkoctfile script is turned off.

        Thanks,
                Andy

#include <octave/pager.h>
#include <iostream.h>
#include <octave/oct.h>
#include <octave/pt-fvc.h>
#include <math.h>
#include <octave/f77-fcn.h>



extern "C"
{

  int F77_FCN (dgenunf, DGENUNF) (const double&, const double&,
                                  double&);
}


double tt;
static tree_fvc *funk;
octave_value_list *fin,fout;

DEFUN_DLD(amebsa, args, nargout,"") {

  Matrix p=args(0).matrix_value();
  ColumnVector y=args(1).vector_value();        
  int ndim=(int)args(2).double_value();
  double ftol=args(4).double_value();
  funk = is_valid_function (args(5), "amebsa", 1);
  int iter=(int)args(6).double_value();

  double ranval=0;
  int ihi,ilo,n,mpts;
  double rtol,swap,yhi,ylo,ynhi;
  mpts=ndim+1;

  octave_value_list retval;


  ColumnVector psum(ndim),pb(ndim);

    for (;;) {
      ilo=0;
      ihi=1;
      F77_XFCN (dgenunf, DGENUNF, (0.0, 1.0, ranval));
      ynhi=ylo=y(0)+tt*log(ranval);
      yhi=y(1)+tt*log(ranval);
      if (ylo > yhi) {
        ihi=0;
        ilo=1;
        ynhi=yhi;
        yhi=ylo;
        ylo=ynhi;
      }

      if (rtol < ftol || iter < 0) {
        swap=y(0);
        y(0)=y(ilo);
        y(ilo)=swap;
        for (n=0;n<ndim;n++) {
          swap=p(0,n);
          p(0,n)=p(ilo,n);
          p(ilo,n)=swap;
        }
        break;
      }
      iter -= 2;
    }

  return (retval);
}

-- 
address@hidden



-------------------------------------------------------------
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
-------------------------------------------------------------



reply via email to

[Prev in Thread] Current Thread [Next in Thread]