110 add convection term by HarryMclean · Pull Request #114 · ExeQuantCode/HeatFlow

HarryMclean · 2026-03-26T18:24:44Z

This merge is because while we didn't fully implement the convection term, we did implement, PETsc library's for solving and output file compression.

…lver only works for symmetric Hmatrix

…tFlow into 110-add-convection-term

…tFlow into 110-add-convection-term Fixed merge

…ry output transiant output Added a new input flag _CompressedOutput to param.in to enable binary 'stream' output in mod_output.f90. Also added a USER_MANUAL.md to docs/.

Copilot

Pull request overview

This PR modernizes HeatFlow’s solve and output pipeline in preparation for adding a convection term, by introducing a PETSc-based linear solver, switching the sparse H-matrix build to CSR, and adding optional compressed (binary) output plus per-material velocity input.

Changes:

Add a new PETSc KSP-based solver and wire it into the time evolution.
Build the H matrix directly into CSR arrays (acsr/ia/ja) rather than the previous sparse format path.
Add _CompressedOutput and binary stream output, and add vel to materials/grid cells.

Reviewed changes

Copilot reviewed 16 out of 16 changed files in this pull request and generated 16 comments.

Show a summary per file

File	Description
`src/heatflow/mod_setup.f90`	Builds CSR H-matrix (`acsr/ia/ja`) and adds progress logging.
`src/heatflow/mod_petsc_solver.f90`	New PETSc wrapper to solve CSR systems with configurable preconditioning.
`src/heatflow/mod_evolve.f90`	Switches evolution solve step to PETSc CSR solver; adds diagnostics and CSR integer conversion.
`src/heatflow/mod_output.f90`	Adds optional compressed binary output (`.bin`) and uses global `logname`.
`src/heatflow/mod_inputs.f90`	Adds `_CompressedOutput` flag and parses material `vel` vector.
`src/heatflow/mod_constructions.f90`	Adds `vel` to `heatblock` and `material` types; removes some unused material fields.
`src/heatflow/mod_material.f90`	Extends `material()` to output `vel`; adjusts assignments for input-driven materials.
`src/heatflow/mod_global.f90`	Adds CSR arrays and global `logname`.
`src/heatflow/mod_hmatrix.f90`	Adds (currently commented-out) hooks and helpers for convection conductivity.
`src/heatflow/mod_heating.f90`	Updates heater logic and adds new heater cases; modifies tau usage and comments out radiative term.
`src/heatflow.f90`	Initializes/finalizes PETSc in the main program.
`Makefile`	Reworks build to include PETSc headers/libs and separate compilation/link steps.
`docs/USER_MANUAL.md`	Adds a new user manual covering inputs and run instructions.
`src/.vscode/settings.json`	Adds nested VS Code settings under `src/`.
`README.md`	Fixes install path capitalization (`HeatFlow`).
`src/heatflow/mod_tempdep.f90`	Comments out temperature-dependent property update hook.
`src/heatflow/mod_boundary.f90`	Adds (commented-out) boundary convection term scaffolding.

Comments suppressed due to low confidence (1)

src/heatflow/mod_setup.f90:129

sparse_Hmatrix allocates acsr/ja/ia every time it is called without first deallocating existing allocations. If this routine can be invoked more than once per run (e.g., for temperature-dependent properties or test harnesses), this will leak memory or error on re-allocation. Add if (allocated(...)) deallocate(...) before allocating, or reuse existing arrays when possible.

      if (Periodicx) addit = [addit, (nx-1)]
      if (ny .gt. 1) addit = [addit, nx] ! Add the values to add to the row to get the column
      if ((Periodicy).and.(ny .gt. 1)) addit = [addit, (ny-1)*nx]
      if (nz .gt. 1) addit = [addit, nx*ny]  ! Add the values to add to the row to get the column
      if ((Periodicz).and.(nz .gt. 1)) addit = [addit, (nz-1)*ny*nx]

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.

Copilot · 2026-03-26T18:33:59Z

docs/USER_MANUAL.md

+#### Flags (Logical)
+All flags default to `.False.`. Set to `.True.` (or `T`) to enable.
+- `_Check_Sparse_Full`: Check if simulation is sparse or full.
+- `_Check_Stability`: Perform stability check.
+- `_Check_Steady_State`: Check for steady state convergence.
+- `_WriteToTxt`: Enable writing output to text files.
+- `_Percentage_Completion`: Show progress % in output.
+- `_Test_Run`: Flag for test runs.
+- `_InputTempDis`: Load initial temperature distribution from file.
+- `_FullRestart`: Perform a full restart.
+
+#### Output Control
+Defines the region of the grid to write to output.
+| Keyword | Type | Default | Description |
+| :--- | :--- | :--- | :--- |
+| `write_every` | Integer | `1` | Write output every N steps. |
+| `start_ix`, `end_ix` | Integer | `1`..`Nx` | X-range for output. |
+| `start_iy`, `end_iy` | Integer | `1`..`Ny` | Y-range for output. |
+| `start_iz`, `end_iz` | Integer | `1`..`Nz` | Z-range for output. |
+


The user manual documents many param.in flags, but it doesn't mention the new _CompressedOutput flag or describe the binary stream output format (.bin) produced when it is enabled. Add documentation for _CompressedOutput (default, file naming, and how to read the binary layout: timestamp followed by the subregion temperature array).

Copilot · 2026-03-26T18:33:59Z

src/heatflow/mod_hmatrix.f90

+      ! if not an edge element
+      if ((x_in .ge. 1) .and. (x_in .le. nx) .and. (y_in .ge. 1) .and. &
+          (y_in .le. ny) .and. (z_in .ge. 1) .and. (z_in .le. nz)) then
+
+         rho = grid(x_out,y_out,z_out)%rho
+         CV = grid(x_out,y_out,z_out)%heat_capacity
+         vel_in = grid(x_in,y_in,z_in)%vel
+         vel_out = grid(x_out,y_out,z_out)%vel


calculate_convective_conductivity only bounds-checks the in indices, but then unconditionally reads grid(x_out,y_out,z_out) for rho/CV/vel_out/Length. If this is called for boundary neighbors (e.g., x_out = 0 or nx+1 when convection terms are enabled), it will access out-of-bounds memory. Add bounds checks for the out indices (or route boundary cases through a dedicated boundary-term routine).

Copilot · 2026-03-26T18:33:59Z

src/heatflow/mod_output.f90

+          if (CompressedOutput) then
+             write(*, *) 'Writing Temperature (Compressed) to file'
+             open(logunit,file=logname, status='unknown', access='stream', position='append')
+             write(logunit) real((itime-1)*(time_step), kind=real12)
+             write(logunit) (Temp_cur(start_ix:end_ix, start_iy:end_iy, start_iz:end_iz))
+             close(logunit)


Same issue as the initial write: this OPEN uses access='stream' but doesn't specify form='unformatted' even though the code performs unformatted write(logunit) operations. This will fail at runtime on many compilers unless FORM is set appropriately.

Copilot · 2026-03-26T18:33:59Z

src/heatflow/mod_evolve.f90

+   ! Convert to 32-bit integers for PETSc (only on first call)
+   if (.not. allocated(ia32)) then
+      allocate(ia32(size(ia)), ja32(size(ja)))
+      ia32 = int(ia, kind=kind(ia32))
+      ja32 = int(ja, kind=kind(ja32))
+   end if
+   NA32 = int(NA, kind=kind(NA32))
+
+   call solve_petsc_csr(NA32, ia32, ja32, acsr, S, x, tol, itmax)
+


solve_petsc_csr requires PETSc to be initialized via PetscInitialize, but simulate calls it unconditionally and the module doesn't enforce initialization. This breaks non-main entry points (e.g., unit tests or any caller that invokes simulate without calling petsc_init). Consider initializing PETSc lazily (once) inside the PETSc wrapper, or adding/using an explicit init guard and updating callers/tests accordingly.

Copilot · 2026-03-26T18:34:00Z

src/heatflow/mod_heating.f90

             area = grid(ix,iy,iz)%Length(1)*grid(ix,iy,iz)%Length(2) !???
             !tau divided by time_step squared in setup.f90
-             tau = grid(ix,iy,iz)%tau*(time_step**2.0_real12) 
+             tau = grid(ix,iy,iz)%tau*(time_step) 


grid%tau is stored as tau/time_step^2 in setup (grid%tau = tau*inverse_time*inverse_time), and both hmatrix and the previous heater logic assume that scaling. Multiplying by only time_step here yields tau/time_step rather than the physical tau, which will change the Cattaneo/heater transient terms. This should likely be time_step**2 (or adjust the stored scaling consistently everywhere).

Suggested change

tau = grid(ix,iy,iz)%tau*(time_step)

tau = grid(ix,iy,iz)%tau*(time_step**2)

Copilot · 2026-03-26T18:34:02Z

src/heatflow/mod_inputs.f90

       CALL assignD(buffer,"heat_capacity",dum_mat(i)%heat_capacity,readvarmat(1))! assign heatCapacity
-       CALL assignD(buffer,"h_conv"       ,dum_mat(i)%h_conv       ,readvarmat(2))! assign h_conv
-       CALL assignD(buffer,"kappa"        ,dum_mat(i)%kappa        ,readvarmat(3))! assign kappa
-       CALL assignD(buffer,"kappa3D"      ,dum_mat(i)%kappa3D      ,readvarmat(4))! assign kappa3D
-       CALL assignD(buffer,"rho"          ,dum_mat(i)%rho          ,readvarmat(5))! assign rho
-       CALL assignD(buffer,"sound_speed"  ,dum_mat(i)%sound_speed  ,readvarmat(6))! assign sound_speed
-       CALL assignD(buffer,"tau"          ,dum_mat(i)%tau          ,readvarmat(7))! assign tau
-       CALL assignD(buffer,"em"            ,dum_mat(i)%em            ,readvarmat(8))! assign e
+       CALL assignD(buffer,"kappa"        ,dum_mat(i)%kappa        ,readvarmat(2))! assign kappa
+       CALL assignD(buffer,"rho"          ,dum_mat(i)%rho          ,readvarmat(3))! assign rho
+       CALL assignD(buffer,"tau"          ,dum_mat(i)%tau          ,readvarmat(4))! assign tau
+       CALL assignD(buffer,"em"           ,dum_mat(i)%em           ,readvarmat(5))! assign e
+       CALL assignV(buffer,"vel"          ,dum_mat(i)%vel          ,readvarmat(6)) ! assign velocity
    end do read


With readvarmat now tracking only 6 fields and vel included, every material block in mat.in must specify vel or the parser will stop with an error (any(readvarmat .ne. 1)). This is a breaking input-format change; consider defaulting vel to (0,0,0) and not requiring it in the completeness check, or document that vel is mandatory for all materials.

Copilot · 2026-03-26T18:34:02Z

src/heatflow/mod_setup.f90

+
+      deallocate(row_vals, row_cols)
+      write(*,'(A,I12,A)') " CSR matrix built successfully. Actual nonzeros: ", count, ""
   end subroutine sparse_Hmatrix
 !!!#################################################################################################



The comment says the CSR arrays "will grow if needed", but the code writes into acsr/ja with no bounds check and no growth logic. If nnz_estimate is ever underestimated (e.g., due to stencil/periodic logic changes), this will write past the end of the arrays. Add a count+row_count > size(acsr) check and reallocate/grow safely (or compute nnz exactly up-front).

Copilot · 2026-03-26T18:34:02Z

src/heatflow/mod_output.f90

-          write(logunit,*) real((itime-1)*(time_step)), &
-               (Temp_cur(start_ix:end_ix, start_iy:end_iy, start_iz:end_iz))
+          if (CompressedOutput) then
+             open(logunit,file=logname, status='unknown', access='stream', position='append')


When CompressedOutput is enabled the file is opened as access='stream' but the subsequent write(logunit) statements are unformatted. Without form='unformatted' on the OPEN, this will be a runtime I/O error (default FORM is formatted). Update the OPEN to explicitly use unformatted stream (and keep it consistent in both compressed write sites).

Suggested change

open(logunit,file=logname, status='unknown', access='stream', position='append')

open(logunit,file=logname, status='unknown', access='stream', form='unformatted', position='append')

Copilot · 2026-03-26T18:34:02Z

src/heatflow/mod_evolve.f90

+      write(*,*) "RHS S: min=", minval(S), " max=", maxval(S), " avg=", sum(S)/size(S)
+      write(*,*) "Temp_p: min=", minval(Temp_p), " max=", maxval(Temp_p), " avg=", sum(Temp_p)/size(Temp_p)
+      write(*,*) "Matrix acsr: min=", minval(acsr), " max=", maxval(acsr), " avg=", sum(acsr)/size(acsr)
+      write(*,*) "Matrix size: n=", NA32, " nnz=", size(acsr)


NA32 is printed in the PETSc diagnostics before it is assigned (it is only set a few lines later). This can print a garbage matrix size and makes debugging misleading; assign NA32 before this block or print NA directly.

Suggested change

write(*,*) "Matrix size: n=", NA32, " nnz=", size(acsr)

write(*,*) "Matrix size: n=", NA, " nnz=", size(acsr)

Copilot · 2026-03-26T18:34:03Z

src/.vscode/settings.json

+{
+    "python-envs.defaultEnvManager": "ms-python.python:conda",
+    "python-envs.defaultPackageManager": "ms-python.python:conda",
+    "python-envs.pythonProjects": []
+}


This adds editor-specific VS Code Python environment settings under src/.vscode/. The repo already has a top-level .vscode/settings.json, so this nested settings file is likely accidental and can cause inconsistent workspace behavior. Consider removing src/.vscode/settings.json and keeping editor config only at the repo root (or add it to .gitignore).

Suggested change

{

"python-envs.defaultEnvManager": "ms-python.python:conda",

"python-envs.defaultPackageManager": "ms-python.python:conda",

"python-envs.pythonProjects": []

}

{}

…=1 mpiexec -n 2 ../HeatFlow/bin/ThermalFlow.x"

HarryMclean and others added 27 commits September 2, 2025 13:05

Add ability to read velocity vector from mat.in

51a5af4

Add convection to hmatrix and boundary vector

3bfe37a

Make fix bug, needs knew solver to work, current implementation of so…

09face9

…lver only works for symmetric Hmatrix

Add non symmatrix sparse matrix solver from the intel library mkl

247f862

Fix typo

2fe28bc

Remove not used mat.in variables

0f1ddc0

Merge branch '110-add-convection-term' of github.com:ExeQuantCode/Hea…

7a858d0

…tFlow into 110-add-convection-term

Make with the parallel libraries

31cf45d

Switch to pardiso solver

440e168

Fix bug with output file

96c2ee3

Fix typo

95f9009

Add cattaneo heating

072c957

Adjust heating

5a375f5

Working Bicgstab

3b4da8a

Switch to petsc

2bd3703

Fix petsc finaliser

763ec50

Working petsc, not fully reproducing old yet

ac8835c

Fixed petsc

2e26177

remove debugging prints

61acedf

make more memory efficient

9307f85

Fixed memory issues

b2b58f5

Enable threading

8687aaf

Add ILU, and fix heating case 2

90f4354

Merge branch '110-add-convection-term' of github.com:ExeQuantCode/Hea…

59abd84

…tFlow into 110-add-convection-term Fixed merge

Change blas and lapack to work with threading

83c6bfe

Add ability to change preconditioner and add AMG

2799327

Code keeps filling up my computer! I've implemented a compressed bina…

a1137b5

…ry output transiant output Added a new input flag _CompressedOutput to param.in to enable binary 'stream' output in mod_output.f90. Also added a USER_MANUAL.md to docs/.

HarryMclean requested a review from Copilot March 26, 2026 18:24

HarryMclean self-assigned this Mar 26, 2026

HarryMclean linked an issue Mar 26, 2026 that may be closed by this pull request

Add convection term #110

Open

Copilot started reviewing on behalf of HarryMclean March 26, 2026 18:25 View session

Merge branch 'main' into 110-add-convection-term

b015110

Copilot AI reviewed Mar 26, 2026

View reviewed changes

HarryMclean and others added 3 commits March 26, 2026 19:02

Make PETsc dynamically discovered

c6cf929

Implement MPI solver run with "OMP_NUM_THREADS=1 OPENBLAS_NUM_THREADS…

d922e0e

…=1 mpiexec -n 2 ../HeatFlow/bin/ThermalFlow.x"

Improve fpm compilation

5fe5f03

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

110 add convection term#114

110 add convection term#114
HarryMclean wants to merge 31 commits intomainfrom
110-add-convection-term

HarryMclean commented Mar 26, 2026

Uh oh!

Copilot AI left a comment

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Copilot AI Mar 26, 2026

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

4 participants

	tau = grid(ix,iy,iz)%tau*(time_step)
	tau = grid(ix,iy,iz)%tau(time_step*2)

	open(logunit,file=logname, status='unknown', access='stream', position='append')
	open(logunit,file=logname, status='unknown', access='stream', form='unformatted', position='append')

	write(,) "Matrix size: n=", NA32, " nnz=", size(acsr)
	write(,) "Matrix size: n=", NA, " nnz=", size(acsr)

Conversation

HarryMclean commented Mar 26, 2026

Uh oh!

Copilot AI left a comment

Choose a reason for hiding this comment

Pull request overview

Reviewed changes

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 26, 2026

Choose a reason for hiding this comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

4 participants