Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Pyrolysis Simulation Time #13187

Open
johodges opened this issue Jul 17, 2024 · 28 comments
Open

Pyrolysis Simulation Time #13187

johodges opened this issue Jul 17, 2024 · 28 comments
Assignees

Comments

@johodges
Copy link
Collaborator

Describe the bug
I was revisiting our internal kinetics optimization routine to integrate optimization of effective battery properties using a similar approach. I noticed the run time of the solid phase only simulations are higher in the latest versions of FDS in the cases used in the optimization. We had a similar discussion on this a few years ago with a difference between FDS 6.2 and 6.7 in #6848. The change made in 6.7.1 in the calculation of temperature dependent material properties brought the runtime down to levels in line with version 6.2.

The attached case runs ~50% slower on the latest FDS than on 6.7.1. On my local machine, the timing for 6.7.1 is ~4 seconds and for the release version of 6.9.1 is ~6 seconds. In my head I was expecting 6.9.1 to run a bit faster since we removed the overhead associated with openmp in the version I was benchmarking. In terms of solid phase only optimization, it may not be too big of a deal, but I am wondering if the difference is also having an impact in HT3D calculations where we are using a more refined grid in the solid phase.

Desktop (please complete the following information):

  • OS: Windows 11
  • Version: FDS 6.9.1 and 6.7.1

Additional context

coneTest2_s6.txt

@johodges
Copy link
Collaborator Author

I am looking into this a bit more. There were some syntax changes to the solid phase pyrolysis routine between versions 6.7.7 and 6.7.8. I am wondering if the changes in syntax are resulting in the difference in runtime since what is being requested by the solid phase solver is not the same. I'll let you know what I find.

@drjfloyd
Copy link
Contributor

We have also made some recent changes in how small the timestep can get as a layer THICKNESS burns away to reduce numerical instablities.

@johodges
Copy link
Collaborator Author

How recent of change are you thinking? Here are the results of a timing study I ran with the input file above:

6.7.4 - 3.5 s
6.7.5 - 3.9 s
6.7.6 - 4.8 s
6.7.7 - 4.8 s
6.7.8 - 5.5 s
6.7.9 - 5.8 s
6.9.1 - 5.8 s

The jump in between 6.7.7 and 6.7.8 could be due to the syntax change, but that is only a part of the difference.

@drjfloyd
Copy link
Contributor

what I was thinking of may have been since 6.9.1.

@rmcdermo
Copy link
Contributor

I recompiled 6.7.1 and compared with latest source on spark. 6.7.1 took 8.236 s and latest took 10.734 s. Nothing obvious to me from the .out files as to why.

@drjfloyd
Copy link
Contributor

For the latest source vs. Jonathan's observation of the release version of 6.9.1 could be changes made since 6.9.1 to timestepping when THICKNESS changes.

@rmcdermo
Copy link
Contributor

I added an output quantity of SUBSTEPS both DEVC and BNDF and only see values of 1 for the latest. So, doesn't that mean we are taking DT_BC=1 with 1 substep for the solid phase?

@johodges
Copy link
Collaborator Author

When I compile the latest master on windows the timing is similar to the pre-compiled 6.9.1 timing for both at ~5.8 s. The 8.2 s on spark verse 5.8 s on my local machine is odd to me. Were you using the build openmp activated? What's the clock speed per core on spark?

@rmcdermo
Copy link
Contributor

Not using openmp.

Intel(R) Xeon(R) Gold 6338 CPU @ 2.00GHz

@gforney
Copy link
Contributor

gforney commented Jul 17, 2024 via email

@johodges
Copy link
Collaborator Author

My laptop is using an Intel Core i7-13800H Processor. The website says max turbo frequency of 5.2 GHz, but task manager shows FDS capping out at ~4.0 GHz.

@marcosvanella
Copy link
Contributor

marcosvanella commented Jul 17, 2024 via email

@mcgratta mcgratta self-assigned this Jul 22, 2024
@mcgratta
Copy link
Contributor

I'll take a look at this today. I want to look at the _cpu.csv file.

@mcgratta
Copy link
Contributor

These are the results I get for the file above on spark:

coneTest2_s6_671.out  Time Stepping Wall Clock Time (s):        7.688
coneTest2_s6_674.out  Time Stepping Wall Clock Time (s):        7.650
coneTest2_s6_675.out  Time Stepping Wall Clock Time (s):        7.774
coneTest2_s6_676.out  Time Stepping Wall Clock Time (s):        9.963
coneTest2_s6_680.out  Time Stepping Wall Clock Time (s):        9.757
coneTest2_s6_new.out  Time Stepping Wall Clock Time (s):        9.823

It appears that there is a significant change between 6.7.5 and 6.7.6.

@mcgratta
Copy link
Contributor

It appears that the increased time in the pyrolysis routine dates back to this

df3b91c975518ceb3e1bf3f38cf97d327d1e7ffd is the first bad commit
commit df3b91c975518ceb3e1bf3f38cf97d327d1e7ffd
Author: drjfloyd <jfloyd@jensenhughes.com>
Date:   Wed Feb 17 13:49:26 2021 -0500

    FDS Source: New wall renoding

 Source/func.f90 |  50 +++++++++++++-----
 Source/read.f90 |  35 +++++++++----
 Source/type.f90 |   6 ++-
 Source/wall.f90 | 155 ++++++++++++++++++++++++++++++++++++++++++++------------
 4 files changed, 186 insertions(+), 60 deletions(-)

@drjfloyd
Copy link
Contributor

This was to reduce instabilities when the wall is renoded. Because we store temperatures at the node center, when we renode, we have to recompute the surface temperature to preserve the correct assumed relationship at boundary. We were getting large swings in the surface temperature which impacted heat transfer, pool evaporation, etc.

@drjfloyd
Copy link
Contributor

You can up the RENODE_DELTA_T on SURF to limit this effect. It defaults to 2 K.

@johodges
Copy link
Collaborator Author

johodges commented Jul 22, 2024

Thanks for looking into it. I am not sure renoding is the issue in this case. The case I uploaded here does not heat very quickly with the surface only heating ~2C over the full simulation of 1800 seconds. I re-ran the case with RENODE_DELTA_T=100 and it did not change the computational time significantly. (ran the comparison with fds 6.9.1 release)

@drjfloyd
Copy link
Contributor

Try adding a PROF for the temperature profile and look to see if there is anything odd with noding.

@mcgratta
Copy link
Contributor

I performed a git bisect using the simple test case at the top. Double check what I did. Run this case with the version that was identified by the git bisect and then run with the version just prior to that . If you see a big change in time, and it is not the renoding, then there might be something else in that commit causing the increased time.

@johodges
Copy link
Collaborator Author

johodges commented Jul 22, 2024

I ran the case with the two bounding commits and did see a difference in computational time.

df3b91c - 4.1 seconds
287e163 - 3.2 seconds

The re-noding commit had the same time with RENODE_DELTA_T=100 and without setting it. I did not see anything odd in the gridding. The attached script will generate the attached video of the time-resolved PROF.

https://github.com/user-attachments/assets/7791f6b9-7824-44b8-a1d5-691ee19a6e92
plot_prof_with_time.txt
coneTest2_s6.txt

Edit: Sorry for the inconsistency in times between this and my earlier timings benchmark post. I had to run this benchmark through my linux boot since my windows machine is busy on another model right now.

@johodges
Copy link
Collaborator Author

johodges commented Jul 23, 2024

Interestingly, most of the computational time is coming from the REMESH_CHECK_IF loop when REMESH_CHECK is True (i.e., the temperature does not exceed RENODE_DELTA_T). If I force the loop to always go through the temperature check (setting RENODE_DELTA_T=1e-10), the computational time in df3b91c drops to 2.2 seconds compared to 4.1 seconds when RENODE_DELTA_T is not set.

The part of the loop that RENODE_DELTA_T bypasses checks if the number of cells in each layer should be decreased if the thickness of the SURF decreased by more than TWO_EPSILON_EB.

Edit: Original comment stated REMESH_CHECK is False, but it should have said True.

@mcgratta
Copy link
Contributor

I do not understand the REMESH_CHECK_IF logic. That is, if the temperature deltas within a layer are less than RENODE_DELTA_T, the number of layer cells is computed, but if at least one delta T is greater, then there is no computation of the layer cells. Also, it seems that the layer is remeshed and temperatures are re-interpolated regardless.

@drjfloyd
Copy link
Contributor

Gave this a little thought. Right now we only store from time step to time step the current wall noding. So if we picked something larger than TWO_EPSILON_EB for the remseh check at line 2291, we would undoubtedly find ourselves in the situation where all the nodes for a wall cell never meet the check for the larger number on any given time step but over many time steps do. If we stored both the current wall noding and the noding as of the last remesh (would only need to allocate this when a cell enters the GET_WALL_NODE_WEIGHTS block at line 1965), then we could set a remesh check based on a wall node having a net decrease over one or more timesteps that is some fraction of the node size since the last renode.

@mcgratta
Copy link
Contributor

If we limit the allocation to PYROLYSIS_MODEL=PYROLYSIS_PREDICTED, it might not use too much memory. However, it will be important to follow through the code how ONE_D%X is allocated and reallocated because that is a big part of the HT3D and variable thickness constructs.

@mcgratta mcgratta assigned drjfloyd and unassigned mcgratta Jul 30, 2024
@firemodels firemodels deleted a comment from drjfloyd Jul 30, 2024
@drjfloyd
Copy link
Contributor

drjfloyd commented Aug 6, 2024

Update. Still working on this. While cleaning up the remeshing logic found some other things that might cause future problems. I think I have a better and somewhat simpler set of logic that should reduce the amount of times we do a full remesh. Doing some testing now. Remains to be seen if after doing this we will still need the RENODE_DELTA_T parameter or not.

@mcgratta
Copy link
Contributor

mcgratta commented Aug 6, 2024

Seems to me that there are always steep gradients in temperature near the surface that will trigger that criterion. But I don't know if these Delta T's are really a problem or not.

@drjfloyd
Copy link
Contributor

drjfloyd commented Aug 9, 2024

I cleaned up the input file - fixing units for SPECIFIC_HEAT and removing TMP_BACK which overrides BACKING='INSULATED'.

With this version much less time is spent in wall; however, you still see the time increase from 6.7.5 to 6.7.6 but it comes down again in 6.9.1. Thew new renoding scheme (plus changes to PYROLYSIS) is about 10 % cheaper than 6.7.5. Timing results below. It is interesting that in 6.9.1 we see time in PART (about 3 %) and this case has no PART inputs or DEVC that invoke particles. Plot shows wall thickness. The new stuff matches 6.9.1 but those differ from 6.7.5 and 6.7.6; however, we have spent a lot of time over the last couple of years focusing on energy and mass conservation for pyrolysis so this isn't surprising.

For renoding walls did a few things:

  • FDS now keeps track of the DX array for walls since the last renoding. There is now a SURF parameter REMESH_RATIO that gives the fractional change in node size needed to trigger a renoding check. Default is 0.05 meaning no check on renoding unless one wall node has changed more then 5 % in size since the last renoding.
  • In stead of RENODE_DELTA_T, there is now REAC_RATE_DELTA on MATL. This defines the fractional step change in reaction rate one is willing to tolerate during a remesh. In READ_MATL, FDS computes for each material with reactions the reaction rate and the reaction rate temperature derivative at the reference temperature. WIth that can estimate how much temperature change can be tolerated in a renode before exceeding REAC_RATE_DELTA. Then in 1DHT, FDS determines the minimum temeperature change based on the materials currently reacting in each layer. It can then use that to compare adjacent nodes in each layer to decide whether or not to allow a full renode. This is to me more physical than the fixed RENODE_DELTA_T.
  • Tried to clean up the logic and simplify it. It is still complex given shrink/swell flags and culling small layers but I think clearer.

It passed firebot on pele, but I will wait for a push until after our dev meeting.

<style> </style>
Version DIVG MASS VELO PRES WALL DUMP PART RADI COMM
675 6.5E-05 0.0E+00 0.0E+00 0.0E+00 4.5E-02 3.6E-02 0.0E+00 6.8E-02 1.6E-03
676 6.5E-05 0.0E+00 0.0E+00 0.0E+00 5.2E-02 3.6E-02 0.0E+00 6.8E-02 1.5E-03
691 8.5E-05 0.0E+00 0.0E+00 0.0E+00 4.5E-02 3.8E-02 3.4E-03 3.7E-02 1.1E-03
New 8.4E-05 0.0E+00 0.0E+00 0.0E+00 4.1E-02 3.6E-02 3.4E-03 3.6E-02 1.1E-03

image

conenew.txt

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

6 participants