search

lanl / palm_lanl

LANL Contributions to PArallelized Large-eddy simulation Model (PALM)
2 stars 5 forks source link

Issue with 'Revert topography for constant_flux_layer = none' commit 266ccaf93e950b432e41a36188a67bf9de12faaf #67

Closed katsmith133 closed 4 years ago

katsmith133 commented 4 years ago

@cbegeman It seems I may have spoken too soon when saying I was ok with merging #59. After adding commit 266ccaf93e950b432e41a36188a67bf9de12faaf to the PR to fix the non-AMD closure issue, I only re-tested the non-AMD closure bit-for-bit matching and not whether the AMD closure part of the code still worked. Turns out that commit, in particular I've narrowed it down to this part of the commit:
https://github.com/xylar/palm_les_lanl/blob/266ccaf93e950b432e41a36188a67bf9de12faaf/trunk/SOURCE/init_grid.f90#L1895-L1905 broke the AMD closure part of the code. What I mean by that is now when I try to run that same test case that you provided, the one in: /lustre/scratch3/turquoise/cbegeman/palm/jobs/test_ocean_amd_mcphee_onesided, PALM runs with the smallest possible time step and pretty much doesn't get anywhere. Here is beginning of the RUN_CONTROL output:

Screen Shot 2020-07-15 at 3 03 31 PM

opposed to before that commit when AMD was working:

Screen Shot 2020-07-15 at 3 05 15 PM

If I revert back to before that commit, then PALM runs appropriately with the AMD closure, but we get the same bit-for-bit mismatch between what were the master and amd_closure branches before. Hence the whole reason for this commit.

Not sure exactly what the problem is, but I have a question to hopefully help me better understand what is going on in this part of the code: what is the reason for applying the current code implementation? In the original code before any of the AMD commits, there was just this line for flat topography: topo(nzb+1:nzt+1,:,:) = IBSET( topo(nzb+1:nzt+1,:,:), 0 ) Then before commit 266ccaf93e950b432e41a36188a67bf9de12faaf, it was: topo(nzb+1:nzt,:,:) = IBSET( topo(nzb+1:nzt,:,:), 0 ) and then after the commit it was: https://github.com/xylar/palm_les_lanl/blob/266ccaf93e950b432e41a36188a67bf9de12faaf/trunk/SOURCE/init_grid.f90#L1895-L1905 Is there a reason that the original code with just topo(nzb+1:nzt+1,:,:) = IBSET( topo(nzb+1:nzt+1,:,:), 0 ) would not work with AMD or was incorrect?

cbegeman commented 4 years ago

@katsmith133 To respond quickly, on Monday I identified a bugfix that may address the problem that you raised. It's now posted here. Basically, prho doesn't end up getting updated the bottom boundary which then causes the buoyancy-related term in amd. I haven't merged this in because I'd like to identify why bc_h isn't set up correctly so that we don't have to manually enforce neumann conditions at the boundaries, which is what that PR does.

katsmith133 commented 4 years ago

@cbegeman Ok, cool. I will wait for your updates on that PR, but in the mean time I can check to see if it fixes this issue too. Thanks!

cbegeman commented 4 years ago

@katsmith133 Ok, I think that reverting that line in init_grid is fine. It seems that the bottom boundary must be designated as topography for the bc_h surface array to be set up. I updated the new PR.

katsmith133 commented 4 years ago

@cbegeman Ok, I've pulled down the amd_bugfix branch and am testing it now...