igfuw / UWLCM

University of Warsaw Lagrangian Cloud Model
GNU General Public License v3.0
6 stars 13 forks source link

Add precip_rate diagnostics to blk2m #156

Closed claresinger closed 1 year ago

claresinger commented 1 year ago

This is running, but I don't think it's correct yet. There are negative values in the precip_rate_rr and precip_rate_nr fields. image

@pdziekan Do you see anything obviously wrong with the implementation here? I'm not sure, but I think that what gets returned from rhs_columnwise should be only the terminal velocity advection (not including w or other microphysical sources/sinks of rain), is that correct?

Closes #155

claresinger commented 1 year ago

@pdziekan - My plan was to try and check that the blk1m precipitation flux looked reasonable, but I have encountered a different issue with that route of debugging (#157). If you don't see any obvious problems, what do you suggest I try? Thanks!

pdziekan commented 1 year ago

The erroneous sedimentation fluxes in slvr_blk_2m seem to come from errors in terminal velocity calculations in libcloudphxx/blk_2m.

I've created a libcloudph++ branch in which terminal velocity in blk_2m is calculated using formulas from blk_1m: https://github.com/pdziekan/libcloudphxx/tree/blk_1m_vterm_in_blk_2m

results: Figure_1

I've tried to debug blk_2m velocity formulas, but it came out to be quite complex and I don't have enough time for this right now. I will try to fix them when I'm less busy, but feel free to do it yourself if you feel up to the task :)

The blk_2m velocity code is in: https://github.com/igfuw/libcloudphxx/blob/master/include/libcloudph%2B%2B/blk_2m/terminal_vel_formulae.hpp

Some references:

I think that the blk_2m code was written by Sylwester Arabas (https://github.com/slayoo). He might be able to help.

claresinger commented 1 year ago

Thank you @pdziekan! The figure above is very promising. I had also been looking at terminal velocity and found a few bugs so far in the blk2m formula (https://github.com/igfuw/libcloudphxx/pull/443). But I haven't found the most important bug yet it seems because the precipitation flux field still is very discontinuous.

Screenshot 2023-05-25 at 09 34 04
claresinger commented 1 year ago

@pdziekan can we merge this too?

pdziekan commented 1 year ago

nancheck(ijk) needs to be done by each thread, so it should not be inside an if (rank==0) clause. I will fix that after merging.