[FLASH-USERS] Fwd: Error in ut_hunt ?

Klaus Weide klaus at flash.uchicago.edu
Mon Jan 11 16:03:24 EST 2021


On Fri, 8 Jan 2021, RIGON Gabriel wrote:

> I was working on some simulations and I saw something that might be an error
> in the ut_hunt.F90 subroutine.
> Directory: ./source/flashUtilities/interpolation/oneDim/ut_hunt.F90
> 
> Take a monotonous array (xx) an a value (x) as argument return the highest
> index (i) such as xx <= x < xx(i)

More precisely (dealing here only with the strictly monotonically 
increasing case), ut_hunt attempts to implement:

 Take a monotonous array xx(1:n) and a value (x) as argument, return the 
 highest index (i) such that
                      x < xx(i+1)     if x < xx(1),
             xx(i) <= x < xx(i+1)     if x < xx(n),
             xx(n)                    otherwise.

Therefore I agree that your suggested change appears more appropriate.

Possibly the case 'x.eq.xx(1)' should be handled explicitly in addition, 
as in the three lines commented out with '!!$'.

Ususally this doesn't matter in practice; the subroutine is often used 
like this

   call ut_hunt(xx,n,x,i)

   i  = min(max(i,0),n-1)
   dx = x - xx(i)
   ! Apply linear or higher polynomial interpolation:
   Approximate slope_at_i based on xx(i) and neighboring point(s)
   x_interpolated = x(i) + slope_at_i * dx  [ + higher_terms(dx) ]

(Interpolation becomes extrapolation for x outside of the interval 
[x(1),x(n)].)

So it doesn't matter what index is returned by ut_hunt when x falls 
exactly on the lower or upper boundary, since it is clamped anyway.

Klaus



More information about the flash-users mailing list