calling a recursive function in fortran 77 ?
in [Fortran]
|
Prev: Selecting the first element of a array returned by pack()
Next: size of a derived type containing pointers...
From: Ingo Thies on 29 Mar 2010 15:02 Am 2010-03-25 21:36, schrieb Richard Maine: >> function wrapfunc(x,y) >> wrapfunc = func(x,y) >> end >> >> and call wrapfunc inside func. But I also never felt quite comfortable >> with this solution. > > Recursion is near the top of the list of things that I recommend *NOT* Hmm a similar problem: I want to do two-dimensional quadrature integration using e.g. the SLATEC library (there exists a Fortran 90 transcription). The natural attempt would be to use an implicitly recursive code like the one above. Let f1(x) be the function to be integrated inside the recursion and f2(x) the outer one. Let DQAG be the quadrature routine (see slatec): function f2(x) .... call DQAG(f1...result1...) f2=result1 end and finally in the calling routine call DQAG(f2...result2...) Previously, I had used an own implementation of a Gauss quadrature from the Numerical Recipes and used two copies of it that only differ in their names, e.g. quad1, quad1. Then, you can safely call quad1 inside f2 and integrate f2 with quad2. However, the more efficient DQAG in SLATEC has too complicated dependencies for a convenient duplication. So, isn't there any efficient way to use a given quadrature routine for a two-dimensional quadrature? Or is duplication the only way to do this in a safe way? Since it is not really recursive routine which calls itself explicitely but calls a third routine/function that calls it again, I do not see how to make use of RECURSIVE here. Ingo
From: Gordon Sande on 29 Mar 2010 15:45 On 2010-03-29 16:02:43 -0300, Ingo Thies <ingo.thies(a)gmx.de> said: > Am 2010-03-25 21:36, schrieb Richard Maine: > >>> function wrapfunc(x,y) >>> wrapfunc = func(x,y) >>> end >>> >>> and call wrapfunc inside func. But I also never felt quite comfortable >>> with this solution. >> >> Recursion is near the top of the list of things that I recommend *NOT* > > Hmm a similar problem: > > I want to do two-dimensional quadrature integration using e.g. the > SLATEC library (there exists a Fortran 90 transcription). The natural > attempt would be to use an implicitly recursive code like the one above. > > Let f1(x) be the function to be integrated inside the recursion and > f2(x) the outer one. Let DQAG be the quadrature routine (see slatec): > > function f2(x) > ... > call DQAG(f1...result1...) > f2=result1 > end > > and finally in the calling routine > > call DQAG(f2...result2...) > > Previously, I had used an own implementation of a Gauss quadrature from > the Numerical Recipes and used two copies of it that only differ in > their names, e.g. quad1, quad1. Then, you can safely call quad1 inside > f2 and integrate f2 with quad2. > However, the more efficient DQAG in SLATEC has too complicated > dependencies for a convenient duplication. How so? Any routine called by either copy of DQAC will have completed before it can be called again by either copy. It is very common for a utility routine to be called from two different places. The dynamic call tree will be a tree but the static call dependencies (loosely but incorrectly the "static call tree") need only be a directed acyclic graph. It would be a lattice except for the requirement on the lowest level being a single item so instead it is a DAG. > So, isn't there any efficient way to use a given quadrature routine > for a two-dimensional quadrature? Or is duplication the only way to do > this in a safe way? > > Since it is not really recursive routine which calls itself explicitely > but calls a third routine/function that calls it again, I do not see > how to make use of RECURSIVE here. > > Ingo
From: glen herrmannsfeldt on 29 Mar 2010 16:56 Ingo Thies <ingo.thies(a)gmx.de> wrote: > Am 2010-03-25 21:36, schrieb Richard Maine: (snip) >> Recursion is near the top of the list of things that I recommend *NOT* > Hmm a similar problem: > I want to do two-dimensional quadrature integration using e.g. the > SLATEC library (there exists a Fortran 90 transcription). The natural > attempt would be to use an implicitly recursive code like the one above. > Let f1(x) be the function to be integrated inside the recursion and > f2(x) the outer one. Let DQAG be the quadrature routine (see slatec): (snip) > Previously, I had used an own implementation of a Gauss quadrature from > the Numerical Recipes and used two copies of it that only differ in > their names, e.g. quad1, quad1. Then, you can safely call quad1 inside > f2 and integrate f2 with quad2. However, the more efficient DQAG in > SLATEC has too complicated dependencies for a convenient duplication. > So, isn't there any efficient way to use a given quadrature routine for > a two-dimensional quadrature? Or is duplication the only way to do this > in a safe way? I remember seeing those in libraries in the Fortran 66 days! A better solution is a two dimensional integration routine. In rectangular coordinates, Gaussian quadrature is separable, so all you need is nested DO loops to loop over the two separately. A lot less overhead than nesting subroutines. Though I always liked the higher dimension Gaussian quadrature for other coordinate systems. There is, for example, one for spherical coordinates. (Though the ones I have seen, such as in Abramowitz and Stegun, are of low order.) Many of the fancier integration routines won't work nested, as the numerical methods fail. You can nest a simpler ons inside the more complicated one, though. > Since it is not really recursive routine which calls itself explicitely > but calls a third routine/function that calls it again, I do not see how > to make use of RECURSIVE here. Well, REENTRANT is sometimes a better name. Though in the case above it is, in fact, recursive, there are cases in multithread systems where a routine is called by different threads at the same time. The requirement that local variables not be in static storage is the same in both cases. IBM has, in the past, also used the term REFRESHABLE. In that case, the code and static storage can be reloaded from the original storage medium without change. There are recursive (and reentrant) programs that do things with static storage: DATA PI/0.0/ if(PI.EQ.0) PI=4.*ATAN(1.) You can put that at the beginning of a REENTRANT, and so RECURSIVE routine even though PI is in static storage. It will fail if it is in REFRESHABLE storage. It will also fail if loaded into write protected memory. -- glen
From: Ingo Thies on 30 Mar 2010 06:45 Gordon Sande wrote: >> However, the more efficient DQAG in SLATEC has too complicated >> dependencies for a convenient duplication. > > How so? Any routine called by either copy of DQAC will have completed > before it > can be called again by either copy. Hmm, but wouldn't the recursivity of DQAG transfer to all the routines called within it, so that I would have to duplicate all of them? > It is very common for a utility routine > to be called from two different places. The dynamic call tree will be a > tree > but the static call dependencies (loosely but incorrectly the "static > call tree") > need only be a directed acyclic graph. It would be a lattice except for the > requirement on the lowest level being a single item so instead it is a DAG. Hmm, I am not completely sure that I am understanding this completely. Just a hint would be nice: Given the case that I want do do 2D quadrature with a single 1D routine without modifying the routine itself. What would be the best way for it? Ingo
From: Ingo Thies on 30 Mar 2010 06:49
glen herrmannsfeldt wrote: > A better solution is a two dimensional integration routine. > In rectangular coordinates, Gaussian quadrature is separable, > so all you need is nested DO loops to loop over the two separately. If ever possible I would like not to modify the routines themselves but use them as a blackbox. But each routine does already contain loops over the integrated coordinate. Ingo |