[Help-gsl] Setting maximum step size in ode solvers

2006-06-12 Thread Bernard Siow
Hello,

Is it possible to set the maximum step size in the ode (ordinary
differential equation) solvers? It's pretty important that the solver
doesn't start running away... Also, I will have large and sudden changes
in some of the parameters and I think I will need to reset the step size
to something small whenever there is a large sudden change in parameters
(predictable). 

Any help would be apprecaited,

Thanks,

Bernard


___
Help-gsl mailing list
Help-gsl@gnu.org
http://lists.gnu.org/mailman/listinfo/help-gsl


Re: [Help-gsl] Re: free

2006-06-12 Thread Jochen Küpper
John D Lamb <[EMAIL PROTECTED]> writes:
> Brian Gough wrote:
>> Jochen Küpper writes:

>>> Shouldn't the GSL free-routines be changed like the following patch?

>> My thinking on that was that for most people calling the free()
>> functions on a null pointer is usually an error rather than by design,
>> so it's helpful to get a segmentation fault.  

> A check in the _free codes would also add a small overhead.

Two good arguments; I agree that the current code is what it should
be. 

Maybe a comment could be added in the code:
,
| Index: block/init_source.c
| ===
| RCS file: /cvs/gsl/gsl/block/init_source.c,v
| retrieving revision 1.7
| diff -u -r1.7 init_source.c
| --- block/init_source.c   26 Jun 2005 13:26:59 -  1.7
| +++ block/init_source.c   12 Jun 2006 21:20:39 -
| @@ -74,6 +74,11 @@
|  void
|  FUNCTION (gsl_block, free) (TYPE (gsl_block) * b)
|  {
| +  // no check for b == 0 performed:
| +  // - For most people calling this function on a null pointer is
| +  //   usually an error, rather than by design, so it is helpful
| +  //   to get a segmentation fault.
| +  // - The check imposes a small performance overhead.
|free (b->data);
|free (b);
|  }
`
This might be rephrased shorter/better by a native speaker... Thinking
of it, maybe it should even go into the documentation.

Greetings,
Jochen
-- 
Einigkeit und Recht und Freiheithttp://www.Jochen-Kuepper.de
Liberté, Égalité, FraternitéGnuPG key: CC1B0B4D
(Part 3 you find in my messages before fall 2003.)


___
Help-gsl mailing list
Help-gsl@gnu.org
http://lists.gnu.org/mailman/listinfo/help-gsl


Re: [Help-gsl] Re: free

2006-06-12 Thread John D Lamb
Brian Gough wrote:
> Jochen Küpper writes:
>  > > I checked the source code for gsl_matrix_free and it will fail if passed
>  > > a null pointer.
>  > 
>  > Shouldn't the GSL free-routines be changed like the following patch?
>
> My thinking on that was that for most people calling the free()
> functions on a null pointer is usually an error rather than by design,
> so it's helpful to get a segmentation fault.  
>
> In general, I've never been convinced that free'ing a null pointer is
> particularly useful as opposed to just doing "if (p) free(p)" in the
> application.
A check in the _free codes would also add a small overhead. I think
there's an arguable case for a delete() function that checks for a null
pointer. But GSL is a C library, not a C++ library and there's no new()
function to go along with it, which ought to can throw a bad_alloc
exception, not available in C.

In any case there's no difficulty to write your own inline new() and
delete() functions as needed.

-- 
JDL



___
Help-gsl mailing list
Help-gsl@gnu.org
http://lists.gnu.org/mailman/listinfo/help-gsl


Re: [Help-gsl] Comments/questions about BFGS algorithm implemented for the GSL

2006-06-12 Thread Brian Gough
James Bergstra writes:
 > I am guessing, that the todo list has:
 > - implement revised line search as an alternative implementation of 
 >   "directional_minimize.c" ?

Yes, I want to put this in the next release.

 > - add testing framework for minimization algos
 > - add rosenbrock to test cases

I think we have the rosenbrock function in the tests already.

-- 
Brian Gough


___
Help-gsl mailing list
Help-gsl@gnu.org
http://lists.gnu.org/mailman/listinfo/help-gsl


[Help-gsl] Re: free

2006-06-12 Thread Brian Gough
Jochen Küpper writes:
 > > I checked the source code for gsl_matrix_free and it will fail if passed
 > > a null pointer.
 > 
 > Shouldn't the GSL free-routines be changed like the following patch?

My thinking on that was that for most people calling the free()
functions on a null pointer is usually an error rather than by design,
so it's helpful to get a segmentation fault.  

In general, I've never been convinced that free'ing a null pointer is
particularly useful as opposed to just doing "if (p) free(p)" in the
application.

-- 
Brian Gough


___
Help-gsl mailing list
Help-gsl@gnu.org
http://lists.gnu.org/mailman/listinfo/help-gsl


Re: [Help-gsl] Numerical Integration With gsl_vector ?

2006-06-12 Thread Brian Gough
Daniel J Farrell writes:
 > Is it possible to use two gsl_vector (or even just C arrays) as  
 > arguments to gsl's more advanced integration function?
 > For example, I have made this function which simply takes in two  
 > columns of numbers an X axis and Y axis then computes the area via  
 > the trapezium rule.
 > double integrate_trapezium(gsl_vector *x, gsl_vector *y);
 > I would like to use a better algorithm i.e.:
 > double integrate_fancy_quadrature_algorithm(gsl_vector *x, gsl_vector  
 > *y, other, arguments, if, needed);
 >  What do you think?

Hello,

Take a look at the spline/interpolation chapter of the GSL
manual. Integration is provided there--you can choose the
interpolation type.

-- 
best regards,

Brian Gough

Network Theory Ltd,
Publishing the GSL Manual - http://www.network-theory.co.uk/gsl/manual/


___
Help-gsl mailing list
Help-gsl@gnu.org
http://lists.gnu.org/mailman/listinfo/help-gsl


Re: [Help-gsl] Bug (or not?) in multiroots/fdjacobian.c

2006-06-12 Thread Brian Gough
eknecronzontas writes:
 > The stepsize for finite-differencing in
 > gsl-1.8/multiroots/fdjacobian.c is specified by the
 > lines
 > 
 > > double xj = gsl_vector_get (x, j);
 > > double dx = epsrel * fabs (xj);
 > 
 > This is, of course mathematically correct.
 > Nevertheless, the behavior is less than ideal if one
 > of the elements of the vector 'x' is sufficiently
 > small. This occurs often if one component of the root
 > happens to be zero. If this occurs, then 'dx' can be
 > so small that one of the columns of the Jacobian is
 > identically zero. This immediately leads to NAN's
 > which are not easy to trace.

Hello, 

Thanks for your email.  I agree the behavior is not good for
small/zero x.  Probably the best thing would be to signal an error if
the step size goes to zero or the derivative overflows.  I'll add a
note to the documentation to make it clearer anyway.  

You can change the behavior with a new drop-in routine copying from
multiroots/hybrid.c into your own application, no need to modify the source.

-- 
Brian Gough

Network Theory Ltd,
Publishing the GSL Manual - http://www.network-theory.co.uk/gsl/manual/


___
Help-gsl mailing list
Help-gsl@gnu.org
http://lists.gnu.org/mailman/listinfo/help-gsl


Re: [Help-gsl] Freeing up memory used by GSL matrix

2006-06-12 Thread James Bergstra
On Mon, Jun 12, 2006 at 12:15:42AM -0500, [EMAIL PROTECTED] wrote:
> I am trying to figure out a way to efficiently free the memory after a gsl
> matrix is no longer being used. Of course, one could call
> gsl_matrix_free(). However in my case, the block pointer may or may not be
> allocated.
> 
> A code snippet follows:
> 
> 
> gsl_matrix *matQ_temp = gsl_matrix_alloc(1,1);
> gsl_block_free(matQ_temp->block);
> matQ_temp->blcok = 0;
> 
>   if (matQ_temp->block!=0){
> std::cout<<"true"< gsl_matrix_free(matQ_temp);
>   }
>   else{
> std::cout<<"false"< free(matQ_temp->data);
>   }
> 
> When running this I get the following o/p:
> 
> false
> *** glibc detected *** double free or corruption (!prev): 0x081b39a8 ***
> Aborted
> 
> 
> I thought even after gsl_block_free has been called, the data pointer
> needs to be freed. If yes, then why am I getting this error? Or does
> gsl_block_free() automatically free the data pointer as well?
> 
> Also, is there anything else I need to do to free the memory?

gsl_block_free not only deallocates the block itself, but also the data that it
encapsulates.  Your matrix buffer was free'd on line 2, so that's why you get a
double free in the else block.



-- 
James Bergstra
http://www-etud.iro.umontreal.ca/~bergstrj



___
Help-gsl mailing list
Help-gsl@gnu.org
http://lists.gnu.org/mailman/listinfo/help-gsl