Potential Vorticity Calculation in MPAS-A

[makinde]

Dear @mgduda,

I hope this message finds you well. I am currently working with MPAS-A simulations to study the dynamics of Cut-off Lows (COLs) and facing challenges in accurately calculating potential vorticity (PV). My goal is to calculate MPAS-A potential vorticity and use it to study the dynamics of COLs in MPAS. First, I embark on testing how to calculate the potential vorticity using the existing reanalysis dataset just to build confidence in the potential vorticity equations. Despite using existing tools like NCL and MetPy, there are discrepancies between my calculated PV and ERA5 reanalysis data.

Methods I have tried:
1. NCL Script (`port_vort_isobaric`):
I utilized the NCL function `port_vort_isobaric` to calculate potential vorticity using global zonal and meridional wind and temperature variables from ERA5 reanalysis. However, the calculated PV differs by approximately > +/ - 3.0 PVU from the ERA5 potential vorticity variable.

2. MetPy Function (`potential_vorticity_baroclinic()`):
I also experimented with MetPy's `potential_vorticity_baroclinic()` function. Unfortunately, the results still deviated from the reanalysis data.


I have prepared figures illustrating the discrepancies between my calculated PV and ERA5 data. These figures highlight the differences at specific altitudes (300 hPa). Figure 1 or the first figure shows the spatial distribution of PV, Static stability, and the potential temperature from the NCL port_vort_isobaric while the second figure compares the calculated PV and the ERA5 dataset, as well as shows the bias in the calculated PV.

Figure 1


View attachment 12624
Figure 2




Specific Questions:
1. Is it possible to obtain the potential vorticity variable directly from MPAS-A without needing to manually calculate it?
2. If not, could you guide me on how to calculate PV using Fortran and incorporate it into MPAS-A?

Incorporating it into MPAS because I am unsure about post-processing MPAS-A results at a full global scale and how `convert_mpas` would handle variable grid results. Most suggestions online state that it is better to calculate potential vorticity at a global scale and I think it may be much easier or better to incorporate potential vorticity calculation directly into MPAS than post-processing variable resolution output from MPAS at a global scale using convert_mpas util. I do not know why I have such a thought.

Your expertise would be immensely valuable in addressing these challenges. Thank you for your time and assistance.

Best regards,
Makinde

 

[AlexLojko]

1. Ertel's potential vorticity should be output by default from MPAS as 'ertel_pv'. The only reason I can think of if it is not being output by default is if you removed the variable from the 'stream_list.atmosphere.output'. If it's not there, add back in ertel_pv.

2. You can check the Fortran code regarding how it's computed. It's best practice to calculate PV on the native model vertical level and then convert to a post-processed level. But this shouldn't be necessary as MPAS already outputs PV. I'm not familiar with computing PV with the NCL scripts, but with Metpy, it doesn't do a good job where there are sharp vertical gradients in PV. So you'll see that where you get closer to the Poles or where there are weather systems present in the mid-latitudes, you will tend to get the largest PV difference errors arising. I'm not sure if this is a particular issue with the finite differencing technique employed in Metpy (or whether it is just an issue of vertical levels being too far apart where PV vertical gradients are large), but this is the likely reason why you get these large PV difference errors around the tropopause.

 

[makinde]

1. Ertel's potential vorticity should be output by default from MPAS as 'ertel_pv'. The only reason I can think of if it is not being output by default is if you removed the variable from the 'stream_list.atmosphere.output'. If it's not there, add back in ertel_pv.

2. You can check the Fortran code regarding how it's computed. It's best practice to calculate PV on the native model vertical level and then convert to a post-processed level. But this shouldn't be necessary as MPAS already outputs PV. I'm not familiar with computing PV with the NCL scripts, but with Metpy, it doesn't do a good job where there are sharp vertical gradients in PV. So you'll see that where you get closer to the Poles or where there are weather systems present in the mid-latitudes, you will tend to get the largest PV difference errors arising. I'm not sure if this is a particular issue with the finite differencing technique employed in Metpy (or whether it is just an issue of vertical levels being too far apart where PV vertical gradients are large), but this is the likely reason why you get these large PV difference errors around the tropopause.

Thank you, AlexLojko.
I sincerely appreciate your valuable advice. Your guidance led me to find Ertel’s potential vorticity (ertel_pv) in the "stream_list.atmosphere.output", effectively resolving my issues.

However, I now have a follow-up query. I aim to compare MPAS’s ertel_pv with observations and/or reanalysis data from sources like CFSR and ERA5 on both pressure and isentropic levels. I need the isentropic levels for doing dynamic tracking and the pressure levels for diagnostics involving other variables. While ERA5 already provides PV on pressure levels, CFSR’s PV remains on isentropic levels. Here are my specific questions:

1. How can I transform PV from pressure levels to isentropic levels and What considerations should I keep in mind during this conversion process?
2. Conversely, how can I convert PV from isentropic levels back to pressure levels?

Given that I need PV data in both isentropic and pressure levels for MPAS, CFSR, and ERA, your expert advice would be immensely helpful.

Thank you once again for your assistance.