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.

Sign up for GitHub

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

Jump to bottom

Excessive fall velocity? [BUG/ISSUE] #762

Closed
tscarter opened this issue Jun 25, 2021 · 24 comments
Closed

Excessive fall velocity? [BUG/ISSUE] #762

tscarter opened this issue Jun 25, 2021 · 24 comments
Assignees
@yantosca
Labels
category: Bug Something isn't working

Comments

@tscarter
Copy link

Description of the problem

I have been running v13.0.0 at 2x2.5 fine, but now that I've started to try and run at nested over North America I keep getting an error several hours into my run that states:

GEOS-CHEM ERROR: Excessive fall velocity?
STOP at CALC_FALLVEL, UCX_mod

I haven't messed with any code that I know of that would impact this. Another group member ran into similar issues with her runs at 2x2.5 I believe and was able to resolve that by using an updated restart file.

Description of troubleshooting performed

Because of my labmate's experience, I've tried using different restart files (basic restart file on our server, an updated one from harvard, a regridded 2x2.5 one to nested from the correctly spun up version that I ran for BCs). None of these have helped. I've also had this occur for different days. In this case, I've checked July 1 and 22, 2019 and had the same bug pop up albeit at different time steps. I'm not quite sure what to try next!

GEOS-Chem version

13.0.0

Description of modifications

Describe your modifications here.

Log files

  • Build log (if applicable):

  • Run logs (if applicable):
    log_Andreae_nest.pdf

  • Run script (if applicable):

Software versions

  • CMake version:
  • Compilers (Intel or GNU, and version):
  • NetCDF version:
@tscarter tscarter added the category: Bug Something isn't working label Jun 25, 2021
@yantosca
Copy link
Contributor

Thanks for writing @tscarter. The routine in question is in ucx_mod.F90, here:

geos-chem/GeosCore/ucx_mod.F90

Lines 3092 to 3166 in fdfaef4

!------------------------------------------------------------------------------
! MIT Laboratory for Aviation and the Environment !
!------------------------------------------------------------------------------
!BOP
!
! !IROUTINE: calc_fallvel
!
! !DESCRIPTION: Function CALC\_FALLVEL calculates the terminal velocity of a
! solid particle.
!\\
!\\
! !INTERFACE:
!
FUNCTION CALC_FALLVEL(DENSITY,RADIUS,TCENTER,PCENTER) RESULT(VEL)
!
! !USES:
!
USE ERROR_MOD, ONLY : ERROR_STOP
!
! !INPUT PARAMETERS:
!
REAL(fp),INTENT(IN) :: RADIUS ! Particle radius (cm)
REAL(fp),INTENT(IN) :: DENSITY ! Particle density (kg/m3)
REAL(fp),INTENT(IN) :: TCENTER ! Local temperature (K)
REAL(fp),INTENT(IN) :: PCENTER ! Local pressure (kPa)
!
! !OUTPUT VARIABLES:
!
REAL(fp) :: VEL ! Fall velocity (m/s)
!
! !REMARKS:
! (1) A remark
!
! !REVISION HISTORY:
! 11 Aug 2012 - S. D. Eastham - Initial version
! See https://github.com/geoschem/geos-chem for complete history
!EOP
!------------------------------------------------------------------------------
!BOC
!
! !LOCAL VARIABLES:
!
REAL(fp) :: Vy ! Intermediate velocity (m/s)
REAL(fp),PARAMETER :: eta=6.45e-8_fp ! Constant (kg/(msK))
REAL(fp) :: val_x ! Dimensionless variable
REAL(fp),DIMENSION(3) :: alpha ! Auxiliary variables
REAL(fp) :: PR ! Pressure times radius
!=================================================================
! CALC_FALLVEL begins here!
!=================================================================
DATA ALPHA/1.49e-5_fp,5.02e-6_fp,2.64e-5_fp/
! Sanity check
IF ((RADIUS.le.0.e+0_fp).or.(DENSITY.le.0.e+0_fp)) THEN
VEL=0.e+0_fp
ELSE
! PCENTER (kPa -> Pa) = *1.d3
! RADIUS (cm -> m ) = *1.d-2
! Therefore multiply PR by 10
PR = PCENTER * RADIUS * 10e+0_fp
VAL_X = -1.0e+0_fp*PR/(ALPHA(3)*TCENTER)
VAL_X = ALPHA(2)*TCENTER*EXP(VAL_X)/PR
VAL_X = 1.0e+0_fp + VAL_X + (ALPHA(1)*TCENTER/PR)
Vy = g0*DENSITY*RADIUS*RADIUS*(1.e-4_fp)/(4.5*ETA*TCENTER)
VEL = 0.893e+0_fp * Vy * VAL_X
ENDIF
! Velocities should be of the order of 0.1 m/s
IF (VEL.gt.10.e+0_fp) THEN
CALL ERROR_STOP(' Excessive fall velocity? ', ' CALC_FALLVEL, UCX_mod')
ENDIF
END FUNCTION CALC_FALLVEL

It looks like this the only things that this depends on is the radius & density of the particle, and the temperature and pressure at the grid box. You might want to check to make sure that there isn't any bad input going in (i.e. corrupted met data).

One thing to try is to add some print statements to the IF block right before it stops:

    ! Velocities should be of the order of 0.1 m/s
    IF (VEL.gt.10.e+0_fp) THEN
       print*, '### radius ', radius
       print*, '### density ', density
       print*, '### tcenter ', tcenter
       print*, '### pcenter', pcenter
       print*, '### vel: ', vel
       CALL ERROR_STOP(' Excessive fall velocity? ', ' CALC_FALLVEL, UCX_mod')
    ENDIF

That might help you trace the error. If those all look OK, try addnig more print statements to look at the terms of the equations.

Also tagging @sdeastham here, who is the author of UCX, who may have more ideas.

@sdeastham
Copy link
Contributor

Hi all - I have indeed seen things fall over at this location when there's something more broadly wrong. I've had the fall velocity error check go off when something totally unrelated to the UCX has failed - much like how convection will throw an error (or hundreds of them) when it hits a NaN, even if the NaN didn't originate there. Looking at the restart makes sense, but I'd also suggest looking at the species concentrations immediately prior to the failure in much the same way as Bob is. Given that this is a nested simulation, it's also possible that the transport algorithm is hitting stability issues, although I think that usually happens with the adjoint code.

Overall though I totally agree with Bob's suggestion - this will also highlight if something odd is happening in the meteorology, which is the other way I've seen this get triggered..

@tscarter
Copy link
Author

tscarter commented Jun 25, 2021 via email

@tscarter
Copy link
Author

Following up on this, I believe it was entirely unrelated to UCX as you both guessed. One issue was with the met fields, which this addressed: #667

However, the other issue that I identified was only partially addressed by previously identified issues (b474612). The other issue is that the global run doesn't produce BCs for the full first day of a run and instead creates a file at _0010z.nc4 with the remaining BCs. This only occurs for months at the beginning of a run (i.e., if I run for 3 months, only the first month has this behavior with the BCs).

Screen Shot 2021-06-28 at 2 58 07 PM

Could that be an issue with HISTORY.rc?
HISTORY.txt

Is it related to this? #269

@yantosca
Copy link
Contributor

Thanks for the update @tscarter.

As to the behavior with the boundary conditions, see issue #698. Because the file write time is at the end of a timestep, we had to make a modification in order to keep the file being written out at the end of the last day of a simulation from being overwritten by the first timestep of the next simulation stage.

You should be able to append the 0501 00Z and 0501 0010Z files together with the nco utilities:

ncrcat -hO GEOSChem.BoundaryConditions.20190501_*.nc4 tmp.nc4

Then once you are sure the tmp.nc4 has all the time points, you can do:

mv tmp.nc4 GEOSChem.BoundaryConditions.20190501_0000z.nc4
rm GEOSChem.BoundaryConditions.20190501_0010.nc4

@tscarter
Copy link
Author

tscarter commented Jun 28, 2021 via email •
edited by yantosca

@mattloman
Copy link

Hi, I'm also encountering this issue in an unedited, out of the box 2-year simulation in 13.1.1 (MERRA2 meteorology, 4ºx5º, full chemistry, July 1 2016-2018). I'm using a restart file from http://ftp.as.harvard.edu/gcgrid/geos-chem/10yr_benchmarks/13.0.0/GCClassic/restarts/ but the other restart files I've tried have also given me the same problem.

Here's some of the output from the print statements suggested by @yantosca:

 ### radius   2.638946962287967E-002
 ### density    990.000000017669
 ### tcenter    197.556479136149
 ### pcenter   16.3661502912961
 ### vel:    10.6165913576084

It seems like an issue with the restart file, but it's not clear to me how to fix this based on what's in this thread.

@yantosca
Copy link
Contributor

Thanks for writing @mattloman.

You are just exceeding the threshold (10.61 > 10) by a little bit. I am not sure how the threshold was developed, but I know the UCX development was likely done with a different set of met fields. I think the temperature & pressure of that grid box is OK (it's a stratospheric box).

You could try changing the threshold from 10 to 12 in CALC_FALLVEL and see if that gets you past the error. It might be that the threshold needs refinement.

@sdeastham: any thoughts?

@sdeastham
Copy link
Contributor

@mattloman @yantosca - I'd be worried about that. IIRC the fall velocities are in meters per second, so 10 m/s (i.e. ~800 km/day) was set deliberately to be absurdly high. For reference, a stratospheric particle would be considered to be falling quickly if it dropped at a rate of 2 km/day.

If I'm reading the radius correctly, it's at 0.026 cm, or 260 um. That's absolutely huge - ten times larger than I would expect even under polar conditions. I would take a look to see what the concentrations of H2O and of total NO3 (HNO3 + NIT) are - this suggests to me that one or both of those are blowing up.

@mattloman
Copy link

@sdeastham These values are typical for the locations with high radius:
H2O = 1.838999986312950E-002 total NO3 = 4.000010013560059E-015
Would these be high enough to cause those radius values? I'm not familiar enough with typical stratospheric conditions to be certain.

@sdeastham
Copy link
Contributor

@mattloman - can you confirm the units of those? Is that v/v?

@mattloman
Copy link

@sdeastham Yes, v/v.

@sdeastham
Copy link
Contributor

OK - in that case I think that your water concentration may be much too high. That's 1.8%, or ~18000 ppmv; for the stratosphere I'd expect <10 ppmv pretty much everywhere (I think - it's been a little while since I really thought about this!). I would look into where those very high VMRs are coming from.

@mattloman
Copy link

Thanks! I'll check that out and report back.

@tscarter
Copy link
Author

tscarter commented Jul 7, 2021

Might also be worth trying a different restart file. Our lab group has run into issues using the generic restart files, but using ones from other runs have solved that problem.

@dylanbm
Copy link

dylanbm commented Jul 30, 2021

Hi @mattloman and all - did you end up sorting this out? I'm seeing the same thing in 13.1.1 with GEOS-FP at 4x5 in 08/2016 and wondering if a cause & solution was found. Thanks!

@mattloman
Copy link

@dylanbm I wasn't able to figure out for sure what was causing the high radius - I switched restart files and stopped encountering the issue.

@yantosca Octobox
Copy link
Contributor

yantosca commented Aug 6, 2021

@mattloman @dylanbm I will close out this issue for now. It seems that using a restart file that is well-spun-up seems to avoid this issue. If you are still encountering problems then feel free to open a new issue.

@yantosca yantosca closed this as completed Aug 6, 2021
@dylanbm
Copy link

dylanbm commented Aug 9, 2021 via email

@yantosca yantosca mentioned this issue Oct 12, 2021
Closed
@yangning-code
Copy link

yangning-code commented Nov 1, 2021 •
edited by yantosca

Hi Bob: I'm not sure that's right since I'm using a 2016 restart from the 2010-2019 GEOS-Chem 13.0.0 10-year full-chemistry benchmark.

@tsherwen
Copy link
Contributor

tsherwen commented Jun 17, 2022 •
edited

Hello,

I am getting this error with the 'out of the box' v13.4.1 code for the acid uptake simulation using the shipped Restart file within a couple of timesteps. How did people get around this if using a "well-spun-up" restart file did resolve it?

I have tried using a "well spun up" restart and different dates (inc. 2018/07/01 and 2019/07/01). The runs are using the merra2 meteorology and the reduced atmosphere (47 levels).

@hyungmini
Copy link

Hi, I had the same problem with v13.0.0, full chemistry, APM simulation. In my case, turning off 'Active strat. H2O?' in input.geos made the code going through the error.

@sdeastham
Copy link
Contributor

That is probably an OK workaround if you don't need an active stratosphere. I would keep an eye on your other stratospheric concentrations though - the fact that H2O is getting unreasonable is a bad sign, so simple kludges might just mask rather than fix the issue.

@tessac2
Copy link

tessac2 commented Sep 10, 2022

I have also run into an excessive fall velocity error on my nested grid simulation, after using a spun-up restart file from a run that was at 4x5 resolution (2018-2020), using v13.2.1 full chemistry. I am trying to start my nested grid simulation in July 2019. I was able to successfully run a 1 day nested grid simulation but when I tried to run for 1 month I started to get this issue. I tried to modify the boundary condition files as @yantosca showed in the above example but this did not resolve the issue. I have attached my log file below. I was wondering if anyone else was continuing to have this issue or might have suggestions? Thanks!!

image

[GC.log.201907.txt](https://github.com/geoschem/geos-chem/files/9541360/GC.log.201907.txt)

@CH2-littlefat CH2-littlefat mentioned this issue Feb 13, 2023
Closed
@yantosca yantosca mentioned this issue May 16, 2023
Closed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees

@yantosca yantosca

Labels
category: Bug Something isn't working
Projects
None yet
Milestone
No milestone
Development

No branches or pull requests
9 participants
@hyungmini @tsherwen @yantosca @sdeastham @dylanbm @yangning-code @tscarter @mattloman @tessac2