UEDGE icon indicating copy to clipboard operation
UEDGE copied to clipboard

Possible typo in grid interpolation

Open bendudson opened this issue 2 years ago • 1 comments

There is a suspicious line in this polintp poloidal interpolation function: https://github.com/LLNL/UEDGE/blob/master/bbb/griddubl.m#L1169

        do 30 ixo = ixos, ixof
            if(xo(ixo,iy).gt.xn(ix,iy) .or. ixo.eq.ixof) goto 25
            ixl = ixo
  20    continue
  25    continue

The 30 label is shared with an outer loop, and the 20 label is not used in this function.

The corresponding lines in the radintp radial interpolation function https://github.com/LLNL/UEDGE/blob/master/bbb/griddubl.m#L1108 are:

            do 20 iyo = iyos, iyof
               if(yo(ix,iyo).gt.yn(ix,iy) .or. iyo.eq.iyof) goto 25
               iyl = iyo
  20        continue
  25        continue

I suspect that the radintp version is correct, and the polintp do loop should use 20 rather than 30 or be rewritten as

do ixo = ixos, ixof
   if(xo(ixo,iy).gt.xn(ix,iy) .or. ixo.eq.ixof) exit
   ixl = ixo
end do

Note: I think this bug, if it is a bug, causes extra work to be done but doesn't affect the result: The interpolation is performed for every value of ixo (rather than just one), but the last result computed has the correct value of ixl.

bendudson avatar Jul 19 '23 21:07 bendudson

Thanks Ben. I agree with your assessment that this a bug in polintp. I also like your rewrite of each of these loops (one in polintp and one in radintp), which is more modern syntax that avoids goto statements.

-Tom


Thomas D. Rognlien Email: @.@.> L-440 (B3725) Tel: 925-422-9830 LLNL, 7000 East Ave, P.O. Box 808 Admin support: 925-422-7446 Livermore, CA 94551

From: Ben Dudson @.> Reply-To: LLNL/UEDGE @.> Date: Wednesday, July 19, 2023 at 2:32 PM To: LLNL/UEDGE @.> Cc: Subscribed @.> Subject: [LLNL/UEDGE] Possible typo in grid interpolation (Issue #44)

There is a suspicious line in this polintp poloidal interpolation function: https://github.com/LLNL/UEDGE/blob/master/bbb/griddubl.m#L1169https://urldefense.us/v3/__https:/github.com/LLNL/UEDGE/blob/master/bbb/griddubl.m*L1169__;Iw!!G2kpM7uM-TzIFchu!z0Nxy9U8Q0olqT9HEa2o9oL9WC13FnwPAJuYj-OEXFVCZiJyoFhbC9TkKtzblfAHlhwsch-jDcjPekO2Ptqg_SAVzQ$

    do 30 ixo = ixos, ixof

        if(xo(ixo,iy).gt.xn(ix,iy) .or. ixo.eq.ixof) goto 25

        ixl = ixo

20 continue

25 continue

The 30 label is shared with an outer loop, and the 20 label is not used in this function.

The corresponding lines in the radintp radial interpolation function https://github.com/LLNL/UEDGE/blob/master/bbb/griddubl.m#L1108https://urldefense.us/v3/__https:/github.com/LLNL/UEDGE/blob/master/bbb/griddubl.m*L1108__;Iw!!G2kpM7uM-TzIFchu!z0Nxy9U8Q0olqT9HEa2o9oL9WC13FnwPAJuYj-OEXFVCZiJyoFhbC9TkKtzblfAHlhwsch-jDcjPekO2PtpBH6fPFw$ are:

        do 20 iyo = iyos, iyof

           if(yo(ix,iyo).gt.yn(ix,iy) .or. iyo.eq.iyof) goto 25

           iyl = iyo

20 continue

25 continue

I suspect that the radintp version is correct, and the polintp do loop should use 20 rather than 30 or be rewritten as

do ixo = ixos, ixof

if(xo(ixo,iy).gt.xn(ix,iy) .or. ixo.eq.ixof) exit

ixl = ixo

end do

— Reply to this email directly, view it on GitHubhttps://urldefense.us/v3/__https:/github.com/LLNL/UEDGE/issues/44__;!!G2kpM7uM-TzIFchu!z0Nxy9U8Q0olqT9HEa2o9oL9WC13FnwPAJuYj-OEXFVCZiJyoFhbC9TkKtzblfAHlhwsch-jDcjPekO2PtpDSPSY2g$, or unsubscribehttps://urldefense.us/v3/__https:/github.com/notifications/unsubscribe-auth/AAILAYSDJQ3E4R4FXKX4OT3XRBG7BANCNFSM6AAAAAA2QQHK2M__;!!G2kpM7uM-TzIFchu!z0Nxy9U8Q0olqT9HEa2o9oL9WC13FnwPAJuYj-OEXFVCZiJyoFhbC9TkKtzblfAHlhwsch-jDcjPekO2Pto2bL0aWg$. You are receiving this because you are subscribed to this thread.Message ID: @.***>

trognlien avatar Jul 26 '23 18:07 trognlien