Feature/implicit forcing #1177

gloopydoop · 2024-03-15T16:54:52Z

Sorry, first pull request, please point out where documentation should go.

For some applications (immersed boundary methods, sponge regions, generation of FST etc) it is common to have a forcing term that is linearly proportional to the velocity, we denote this as f = \chi * u, where \chi is a scalar field.

Imposing this in user_f results in an explicit formulation where f = \chi * u^(n) that can lead to stability issues for large \chi.

Here we impose an implicit formulation f = \chi * u^(n+1)

Predicted problems:

I'm 99% sure if \chi is defined on the boundary, we'll have to adjust the pressure boundary conditions.
Nothing has been tested on GPU's etc yet
not a lot of testing has gone into this, I'm only pull requesting so that it can be discussed

timofeymukha · 2024-03-15T17:25:11Z

src/common/user_intf.f90

+  !> Dummy user (fluid) implicit Brinkman forcing
+  subroutine dummy_implicit_brinkman(chi, t)
+    type(field_t), intent(inout) :: chi
+    real(kind=rp), intent(in) :: t
+    chi%x = 0.0
+  end subroutine dummy_implicit_brinkman
+


This should throw an error, as all the other dummies.

Not all dummies I guess, some are "do nothing" but will not return an error. Like dummy_user_check

MartinKarp

I think this is a good feature. One thing I guess as Tim also mentioned is how we want to treat the brinkman term when there is none. Right now it calls a do-nothing routine, but technically we add a vector of zero each time. If we moved out setting h1,h2 outside of computing the residual for the velocity this might be possible to do more elegantly.

MartinKarp · 2024-03-17T14:49:51Z

YOFAM/README.md

@@ -0,0 +1,3 @@
+#Simple turbulent channel case at Re_tau=180
+The channel starts with a slightly perturbed initial condition which becomes turbulent after around 10 time units.
+If one wants to change the mesh and play around, please do so, by generating a new mesh with genmeshbox. 


Should this be in the main code? If so, please expand on how it differs from usual turb_channel, remove the field files, and rename the folder.

none of this example should be in the main code, it was just the case I was using to test the implicit Brinkman term.
I would envision that the immersed_bunny or something similar would be extended as an example.

MartinKarp · 2024-03-17T14:51:00Z

examples/EXAMPLE_FOR_IMPLICIT_BRINKMAN/README.md

@@ -0,0 +1,3 @@
+#Simple turbulent channel case at Re_tau=180
+The channel starts with a slightly perturbed initial condition which becomes turbulent after around 10 time units.
+If one wants to change the mesh and play around, please do so, by generating a new mesh with genmeshbox. 


Kind of similar comment as for the other example here. Please also avoid ALLCAPS in folder names.

MartinKarp · 2024-03-17T14:51:34Z

examples/EXAMPLE_FOR_IMPLICIT_BRINKMAN/turb_channel.f90

+    type(user_t), intent(inout) :: u
+    u%fluid_user_ic => user_ic
+    u%user_mesh_setup => user_mesh_scale
+	 u%user_implicit_brinkman => bman


Suggested change

u%user_implicit_brinkman => bman

u%user_implicit_brinkman => bman

MartinKarp · 2024-03-17T14:52:52Z

examples/EXAMPLE_FOR_IMPLICIT_BRINKMAN/turb_channel.f90

+	 u%user_implicit_brinkman => bman
+end subroutine user_setup
+
+    ! Harry ---------------------------------------


Suggested change

! Harry ---------------------------------------

!> Insert your comments here

MartinKarp · 2024-03-17T14:57:35Z

src/common/.user_intf.f90.swp

MartinKarp · 2024-03-17T15:00:29Z

src/fluid/bcknd/cpu/pnpn_res_cpu.f90

@@ -147,7 +147,7 @@ subroutine pnpn_vel_res_cpu_compute(Ax, u, v, w, u_res, v_res, w_res, &

    do i = 1, n
       c_Xh%h1(i,1,1,1) = mu
-       c_Xh%h2(i,1,1,1) = rho * (bd / dt)
+       c_Xh%h2(i,1,1,1) = rho * (bd / dt) + chi%x(i,1,1,1)


This is nice and will work I think, but maybe we should also consider moving h1, h2 to the Ax_helm type @njansson @timofeymukha . They do not really belong in coef, but are directly linked to the Helmholtz equation. Maybe something we could issue another PR for?

I think the idea of using a global h1 and h2 that sits somewhere and gets overwritten with right values just prior to a call to an axhelm is unnecessary. I would say we can just pass h1 and h2 to axhelm's compute routine as arguments, or is there some reason not to do that? But I agree that if we should dedicate a components for this in a type, it should be axhelm, and not coef.

If we do that the idea of passing something that computes a linear operator Ax to the ksp solvers will kind of not work as different matrices will have different inputs. Because of that I was leaning a bit towards putting it directly in Ax.

I see, that makes sense.

MartinKarp · 2024-03-17T15:03:54Z

src/fluid/fluid_pnpn.f90

@@ -473,6 +473,8 @@ subroutine fluid_pnpn_step(this, t, tstep, dt, ext_bdf)
    type(ksp_monitor_t) :: ksp_results(4)
    ! Extrapolated velocity for the pressure residual
    type(field_t), pointer :: u_e, v_e, w_e
+    ! Field for implicit Brinkman term
+    ! type(field_t), intent(inout) :: chi_fam


Should this be commented or?

MartinKarp · 2024-03-17T15:05:39Z

src/fluid/fluid_pnpn.f90

@@ -581,10 +583,12 @@ subroutine fluid_pnpn_step(this, t, tstep, dt, ext_bdf)

      ! Compute velocity.
      call profiler_start_region('Velocity residual', 19)
+      ! Potentially compute the Brinkman penalization term
+      call this%userbrinkman%compute(t)


I guess it does not need to be user, I forsee a lot could be captured by typical forcings in the json file. Similar to constant forcing.

This is just a naming suggestion

Agree, can't this be put in the forcing? Same with chi, would be nice to hide these case specific variables under some object/abstraction

…the source terms. NOTE. I have used extra fields in the field list so that I didn't break other source terms, but I suspect it would be more understandable if we passed an additional chi field through all source terms

timfelle

Any chance you could go through and cleanup the places where stuff is no longer changed?
It is kinda difficult to know what is still needed and what is not. Are you still using the changes in user_intf.f90 for instance?

src/source_terms/brinkman_source_term.f90

src/source_terms/source_term.f90

src/source_terms/brinkman_source_term.f90

… sponge. NOTE there is certainly some scaling factor wrong with the current chi implementation!

… passive scalars, test implicit vs explicit, implement pnpn_res on device and sx to include the extrapolated Brinkman terms in the pressure solve

Co-authored-by: Tim Felle Olsen <[email protected]>

…o into feature/implicit_forcing

…orcing

…implicit_forcing

gloopydoop added 3 commits March 15, 2024 11:13

for_martin

f29e4a8

save us Martingit add *

592d19e

first draft of implicit Brinkman penalization

3671bbf

timofeymukha reviewed Mar 15, 2024

View reviewed changes

MartinKarp requested changes Mar 17, 2024

View reviewed changes

MartinKarp added enhancement New feature or request help wanted Extra attention is needed labels Mar 17, 2024

timfelle marked this pull request as draft March 19, 2024 09:07

gloopydoop marked this pull request as ready for review March 19, 2024 10:19

gloopydoop marked this pull request as draft March 19, 2024 10:20

gloopydoop added 3 commits March 21, 2024 09:51

clean up a little

160f50c

clean a little more...

9a3c0bc

timfelle requested changes Mar 21, 2024

View reviewed changes

gloopydoop and others added 16 commits March 21, 2024 16:16

cleaned up user_intf. Still need to implement passive scalars

3d7d948

passive scalars, but no testing

14128f2

added DELETE_fst_sandbox which is a temporary test case for testing a…

99e3107

… sponge. NOTE there is certainly some scaling factor wrong with the current chi implementation!

somewhat happy with the FST example. things to do: test scaling, test…

0bdeaae

… passive scalars, test implicit vs explicit, implement pnpn_res on device and sx to include the extrapolated Brinkman terms in the pressure solve

for Martin to look at

0bd38c2

tempdog

6da63a9

Merge branch 'develop' into feature/implicit_forcing

d11ed60

to test hip

a33353c

Update src/source_terms/brinkman_source_term.f90

63da619

Co-authored-by: Tim Felle Olsen <[email protected]>

Update src/source_terms/brinkman_source_term.f90

371bd97

Co-authored-by: Tim Felle Olsen <[email protected]>

Update src/source_terms/brinkman_source_term.f90

221d578

Co-authored-by: Tim Felle Olsen <[email protected]>

Update src/source_terms/source_term.f90

358b74e

Co-authored-by: Tim Felle Olsen <[email protected]>

Update src/source_terms/source_term.f90

13835fa

Co-authored-by: Tim Felle Olsen <[email protected]>

Update src/source_terms/source_term.f90

b1cad71

Co-authored-by: Tim Felle Olsen <[email protected]>

still hip

fecc8cb

Merge branch 'feature/implicit_forcing' of github.com:ExtremeFLOW/nek…

2bd2164

…o into feature/implicit_forcing

gloopydoop and others added 6 commits April 16, 2024 14:10

if/is implicit

7e546e8

hip

de6eff8

removing devision of rho in vel_res

59a7225

Merge remote-tracking branch 'origin/develop' into feature/implicit_f…

8043c86

…orcing

Make sure the extrapolated velocity is not used.

ec52941

Merge commit '2e61eb39ebfa90d5377223eb869691c23476041b' into feature/…

dea3644

…implicit_forcing

timfelle mentioned this pull request May 24, 2024

Benchmark/cylinder ExtremeFLOW/neko-top#16

Closed

timfelle added 9 commits May 28, 2024 10:45

Remove unused file

d5a7ff1

Add new case instead of modifying the bunny

df934d0

Undo changes that are not needed

2cabd9b

Refactor pressure implicit forcing from this PR

57285e9

Second pass of cleanup

1a95870

Third pass

8cd12c4

Remove unused variable

3bb53bc

Ensure case file follow the 80 character rule.

d1b2af2

Merge commit 'b0232d143f933cc52f28223b23f09c92db856255' into feature/…

0d7194d

…implicit_forcing

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Feature/implicit forcing #1177

Feature/implicit forcing #1177

gloopydoop commented Mar 15, 2024

timofeymukha Mar 15, 2024

MartinKarp Mar 17, 2024

MartinKarp left a comment

MartinKarp Mar 17, 2024

gloopydoop Mar 18, 2024

MartinKarp Mar 17, 2024

MartinKarp Mar 17, 2024

MartinKarp Mar 17, 2024

MartinKarp Mar 17, 2024

MartinKarp Mar 17, 2024

timofeymukha Mar 17, 2024

MartinKarp Mar 17, 2024

timofeymukha Mar 18, 2024

MartinKarp Mar 17, 2024

MartinKarp Mar 17, 2024

MartinKarp Mar 17, 2024

njansson Mar 18, 2024

timfelle left a comment

	u%user_implicit_brinkman => bman
	u%user_implicit_brinkman => bman

	! Harry ---------------------------------------
	!> Insert your comments here

Feature/implicit forcing #1177

Are you sure you want to change the base?

Feature/implicit forcing #1177

Conversation

gloopydoop commented Mar 15, 2024

Choose a reason for hiding this comment

Choose a reason for hiding this comment

MartinKarp left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

timfelle left a comment

Choose a reason for hiding this comment