lautenberger / elmfire

Eulerian Level set Model of FIRE spread
https://elmfire.io
Eclipse Public License 2.0
23 stars 11 forks source link

Speed up UCB urban model #16

Open dwimjpurnomo opened 1 year ago

dwimjpurnomo commented 1 year ago

I tried with this input and it works, but much slower than using Hamada.

https://drive.google.com/file/d/119VjwXsv0s4xQr8MhSc7dtfCBecbMNPD/view?usp=sharing

lautenberger commented 1 year ago

Thanks Dwi. Are you using this building_fuel_models.csv file?

dwimjpurnomo commented 1 year ago

Yes, I am

lautenberger commented 1 year ago

One thing I noticed right away - in elmfire_spread_rate.f90 line 407 in the subroutine UMD_UCB_BLDG_SPREAD, this line causes division by zero if V_S is zero:

C%LOW = AMIN1((C%VELOCITY_DMS+C%VBACK)/2/V_S,10.0)

I'll change it to something like:

IF (V_S .GT. 1E-4) THEN
   C%LOW = AMIN1((C%VELOCITY_DMS+C%VBACK)/2/V_S,10.0)
ELSE
   C%LOW = 1.0
ENDIF
lautenberger commented 1 year ago

OK, I created a new branch urban-spread-optimization and you can see changes here:

https://github.com/lautenberger/elmfire/compare/main...urban-spread-optimization

These are simple optimizations - A*A instead of A**2, multiplying by reciprocal cellsize instead of dividing by cellsize, etc. I doubt it's going to bring down the run time significantly though.

The computational challenge is that the linked list LB contains thousands of nodes and we're looping over all of those nodes multiple times. Do you think it's possible to limit the number of LB nodes that are looped over using a distance criterion or similar?

dwimjpurnomo commented 1 year ago

The range of LB is only 100m radius, so its 4-ish cells (16 cells area in total). What I know is we can do active cells, I never try it before, but its like only a limited number of cells within the range is going to be calculated.

lautenberger commented 1 year ago

OK, it seems the biggest problem is probably with the LB list. I did some basic profiling yesterday and saw that the LB started off with a few nodes but quickly grew to 10,000+ nodes. So we're looping over 10,000+ LB nodes for every cell with fuel model 91.

Are you running in Visual Studio for debugging? If not, I highly recommend it. If you get stuck I'm happy to help you get up and running. It will save a lot of time in the long run.

dwimjpurnomo commented 1 year ago

Yes, I use Visual Studio now. Maria helped me set it up.

lautenberger commented 1 year ago

Excellent. Are you seeing any improvement in run time with the urban-spread-optimization branch:

https://github.com/lautenberger/elmfire/compare/main...urban-spread-optimization

If the changes look good to you I'd like to merge it with main before going further.

dwimjpurnomo commented 1 year ago

comparisonSpeed

The changes work fine. The time improvement, however, is only around 2%. This is from 10 repetitions of running. I use one of the tutorial (01 - case)

lautenberger commented 1 year ago

That's about what I expected.

To get a significant speed improvement, a change in the algorithm is going to be needed. Rather than a single global building list that may contain tens of thousands of nodes, I think we're going to have to create a local building list for each "tagged" cell that contains at most ~10 ish building nodes. Right now, lots of cpu time is spent analyzing nodes in the LB list where TARGET_R_METERS is > 1 km.

BTW, I'm going to merge the urban-spread-optimization branch into main

dwimjpurnomo commented 1 year ago

Yes, I agree with that.

lautenberger commented 1 year ago

OK - are you familiar enough with the code at this point to create a new branch and make an initial attempt at it?

dwimjpurnomo commented 1 year ago

Yes, I will do it

dwimjpurnomo commented 11 months ago

I have created a pull request for the first attempt of speeding up the simulation. It makes it 5 times faster. For now, it is just using IF command in which if the burning cell is within the predefined neighborhood, then the subsequent commands are executed, otherwise, they are skipped. It means it still goes through all the burning cells, its just that only a limited number of these cells will be thoroughly calculated. I think this can be made faster by limiting the number of elements in the loop, so the time to execute the IF command will be taken back. Still thinking of a way of doing this

dwimjpurnomo commented 11 months ago

This is the pseudocode in MATLAB for the next speed up. What is the best wayto implement it in FORTRAN.

allIX = [LIST_BURNED.IX]; %make an array of all IX coordinate of burning cells and its corresponding index in LIST_BURNED

IX_INCLUDED = (allIX >= IX_LOW_BORDER && allIX <= IX_HIGH_BORDER); % Make all the values outside the limit to be 0 and 1 otherwise

ind = find(IX_INCLUDED); % Find the index number that have 1 value for IX_INCLUDED

for i = 1:numel(ind) % iterate only for the number element in the ind variable DEL_X = C.IX - LIST_BURNED{ind(i)}.IX % only the burning cell that has index stored in ind are considered

and continue the codes......

end

lautenberger commented 11 months ago

I assume you do the same thing for IY as well?

This seems similar (but not identical) to how the narrow band level set method has been implemented in ELMFIRE. Basically LIST_TAGGED tracks all cells that are within BANDTHICKNESS cells of the fire front. That way you calculate gradients and spend cpu time only in a narrow band adjacent to the fire front. The LIST_TAGGED list is generated by calling the subroutine TAG_BAND (in elmfire_level_set.f90) with inputs IXLOC and IYLOC (x and y locations to around which to tag a band), here's the subroutine:

! *****************************************************************************
SUBROUTINE TAG_BAND(NX, NY, IXLOC, IYLOC, T)
! *****************************************************************************

INTEGER, INTENT(IN) :: NX, NY, IXLOC, IYLOC
REAL, INTENT(IN) :: T
INTEGER :: IXTAGSTART, IXTAGSTOP, IYTAGSTART, IYTAGSTOP, IX, IY

IXTAGSTART = MAX(3,    IXLOC - BANDTHICKNESS) 
IXTAGSTOP  = MIN(NX-2, IXLOC + BANDTHICKNESS)
IYTAGSTART = MAX(3,    IYLOC - BANDTHICKNESS) 
IYTAGSTOP  = MIN(NY-2, IYLOC + BANDTHICKNESS)

DO IY = IYTAGSTART, IYTAGSTOP
DO IX = IXTAGSTART, IXTAGSTOP
   IF (ISNONBURNABLE(IX,IY)) CYCLE
   IF (.NOT. TAGGED(IX,IY) .AND. (.NOT. EVERTAGGED(IX,IY)) ) THEN 
      TAGGED    (IX,IY) = .TRUE.
      EVERTAGGED(IX,IY) = .TRUE.
      CALL APPEND(LIST_TAGGED, IX, IY, .FALSE., T)
      NUM_EVERTAGGED = NUM_EVERTAGGED + 1
      EVERTAGGED_IX(NUM_EVERTAGGED) = IX
      EVERTAGGED_IY(NUM_EVERTAGGED) = IY
   ENDIF
ENDDO
ENDDO

! *****************************************************************************
END SUBROUTINE TAG_BAND
! *****************************************************************************

The 2D matrices TAGGED(IX,IY) and EVERTAGGED(IX,IY) are used to avoid adding a tagged cell to the tagged list more than once.

dwimjpurnomo commented 11 months ago

Yes, its for IY as well. So basically I need to declare a new variable that is the same type as LIST_TAGGED?

dwimjpurnomo commented 11 months ago

I use CYCLE. It tremendously speed up the simulation. But if I use it, the building model is not executed. This is my snippet IF (.....) THEN BUILDING_MODEL ELSE

CYCLE

ENDIF

lautenberger commented 8 months ago

This input deck was posted at Issue #35 which has subsequently been closed, re-posting it here:

https://drive.google.com/file/d/1QZNJwqDyKahPp-am9lqUfx1LAvFSs_Qr/view?usp=sharing