[ncl-talk] how to fill an array with an other

[Rashed Mahmood]

Jonathan, I was responding to Laura's question.
I think Dennis has clearly explained the "looping" issue. Thanks Dennis!

Cheers,
Rashed

On Wed, May 6, 2020 at 11:58 AM Buzan, Jonathan <jbuzan at purdue.edu> wrote:

> Hi Rashed,
>
> One way to speed up this process is to use the ncl ndtooned (n-dimension
> to one dimension) function, then loop over the 1d blocks by an amount equal
> to the spread in lats Lons, etc. When your calculations are completed you
> bring your data back to the original n-dimensions. It is complicated to
> code, and takes some effort, but you can substantially speed up your code
> this way.
>
> Cheers,
> -Jonathan
>
>
>
>
> On May 6, 2020, at 8:52 PM, Dennis Shea via ncl-talk <ncl-talk at ucar.edu>
> wrote:
>
> FYI:
>
> In general, looping in an interpreted language [NCL, Matlab, Python, R,
> ...] is 'slow' relative to a compiled language like fortran, etc.
> Sometimes, quite slow. To my knowledge, one 'key' to fortran's speed is
> that fortran does not allow unrestricted pointers.  This allows the
> compiler to do 'data look-ahead' [my words, I am not a computer science
> person]. The compiler can issue 'pre-fetch' instructions tp load (say)
> cache so when the data are needed it is in fast memory. Please correct me
> if I am wrong.
>
> I would say that looping in NCL is reasonably 'quick.'  As example, see
> the code for the Monthly Water Balance Model described at:
> *http://www.ncl.ucar.edu/Applications/crop.shtml*
> <http://www.ncl.ucar.edu/Applications/crop.shtml>
> Code:   *http://www.ncl.ucar.edu/Applications/Scripts/MWBM_DRIVER.ncl*
> <http://www.ncl.ucar.edu/Applications/Scripts/MWBM_DRIVER.ncl>
>
> In the 3-dimension (time,lat,lon) case, the following loops aver all times
> and grid points
>
> The following occurred surprisingly quickly
>
> ;========================================================
> ; Brute Force: Loop over each time step and all 'land' grid points
> ;========================================================
>
>   ntim = dim_data(0)
> [**SNIP**]
>   elseif (rank_prc.eq.3) then
>     nlat = dim_data(1)
>     mlon = dim_data(2)
>
>     remain  := new( (/nlat,mlon/), type_data, type_fill)
>     remain   =  0.0                                   ; 25.4
>     remain at long_name= "Remain Precip Water"
>     remain at units    = "mm"
>
>     prestor := new( (/nlat,mlon/), type_data, type_fill)
>     prestor  =  0.0                                   ; 25.4
>     prestor at long_name= "Remain Precip Water"
>     prestor at units    = "mm"
>
>     if (isscalar(lsmask) .and. lsmask.eq.1) then
>         lsm := new( (/nlat,mlon/), "integer", "No_FillValue")
>         lsm  = 1
>     else
>         lsm  = lsmask
>     end if
>
>     do nl=0,nlat-1
>       do ml=0,mlon-1
>         if (lsm(nl,ml).eq.1) then   ; land only
>           do nt=0,ntim-1
>              prain(nt,nl,ml) = prc(nt,nl,ml)-snow(nt,nl,ml)
>              rodirect(nt,nl,ml)= prain(nt,nl,ml) *directfac           ; Direct Runoff:  Eq (3)
>              prain(nt,nl,ml)   = prain(nt,nl,ml) -rodirect(nt,nl,ml)  ; remaining rain
>              snstor(nt,nl,ml)  = snstor(nt,nl,ml)+snow(nt,nl,ml)      ; psnow accumulates as snstor (page 1)
>              if (snstor(nt,nl,ml).gt.0 .and. tmp(nt,nl,ml).gt.TSNOW) then
>                  if (snstor(nt,nl,ml).le.10) then
>                      smelt(nt,nl,ml)  = snstor(nt,nl,ml)
>                  else
>                      smeltf(nt,nl,ml) = MELTMAX*((tmp(nt,nl,ml)-TSNOW)/((TRAIN-TSNOW)+0.0001))  ; Eq (5)
>                      if (smeltf(nt,nl,ml).gt.MELTMAX) then
>                          smeltf(nt,nl,ml) = MELTMAX
>                      end if
>                      smelt(nt,nl,ml)  = smeltf(nt,nl,ml)*snstor(nt,nl,ml)       ; Eq (6)
>                      snstor(nt,nl,ml) = snstor(nt,nl,ml)-smelt(nt,nl,ml)
>                      if (snstor(nt,nl,ml).lt.0.0) then
>                          snstor(nt,nl,ml) = 0.0
>                      end if
>                 end if
>              end if
>              prain(nt,nl,ml)   = prain(nt,nl,ml) + smelt(nt,nl,ml)
>              pmpe(nt,nl,ml)    = prain(nt,nl,ml) - pet(nt,nl,ml)    ; net water
>
>              if (pmpe(nt,nl,ml).lt.0) then
>                  sms(nt,nl,ml) =  prestor(nl,ml) - abs((pmpe(nt,nl,ml)*prestor(nl,ml))/WHC)
>                  if (sms(nt,nl,ml).lt.0)  then
>                      sms(nt,nl,ml) = 0
>                  end if
>                  delstor = sms(nt,nl,ml)   - prestor(nl,ml)       ; change in soil moisture
>                  aet(nt,nl,ml) = prain(nt,nl,ml) - delstor       ; JAVA: ae = prain + delstor * (-1.0);
>                  prestor(nl,ml) = sms(nt,nl,ml)                   ; save for next iteration
>                  surplus(nt,nl,ml) = 0.0
>             else
>                  aet(nt,nl,ml) = pet(nt,nl,ml)
>                  sms(nt,nl,ml) = prestor(nl,ml) + pmpe(nt,nl,ml)
>                  if (sms(nt,nl,ml).gt.WHC) then
>                      surplus(nt,nl,ml) = sms(nt,nl,ml) - WHC
>                      sms(nt,nl,ml)  = WHC
>                  end if
>                  prestor(nl,ml) = sms(nt,nl,ml)                  ; save for next iteration
>             end if
>
>             deficit(nt,nl,ml) = pet(nt,nl,ml) - aet(nt,nl,ml)
>             runoff(nt,nl,ml) = (surplus(nt,nl,ml)+remain(nl,ml))*runoffFactor
>             remain     = (surplus(nt,nl,ml)+remain(nl,ml))-runoff(nt,nl,ml)
>             if (remain(nl,ml).lt.0) then
>                 remain(nl,ml) = 0
>             end if
>
>             runoff(nt,nl,ml)  = runoff(nt,nl,ml) + rodirect(nt,nl,ml)
>           end do  ; nt
>         end if    ; lsm
>       end do      ; ml
>      end do       ; nl
>
>   end if          ; rank_data
>
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20200507/1659b202/attachment.html>