How to calculate grid%mu with a corrected surface pressure?

[mirabelle]

I'm working on online WRF assimilation using WRFV3.7.1 , the vertical coordinate is TF.

I'm currently assimilating surface pressure ( grid%psfc) and I had calculated the increment of surface pressure.

Do you know how to use the correction of surface pressure to modify the perturbation dry air mass in column ( grid%mu)?

Dave Gill asked to have this question directly assigned to him

Thanks a lot！

 

[davegill]

Mirabelle,
Thanks for clarifying that this is specifically for the TF vertical coordinate system.

This does not seem like it will be too difficult.

The 3d dry pressure is defined as: p_dry(i,k,j) = eta(k) * mut(i,j) + ptop,
where eta(k) is the WRF vertical coordinate,
mut(i,j) is the dry column pressure,
ptop is the model lid.

If we are interested in the surface pressure, then the full level eta(k=1) is identically defined as 1.0.

So the DRY surface pressure simplifies to p_sfc_dry(i,j) = mut(i,j) + ptop

We know that mut(i,j) = mub(i,j) + mu(i,j),
where mub(i,j) is the reference column pressure (for eaqh (i,j), a time invariant for the simulation),
mu(i,j) is the perturbation column pressure (for each (i,j), mu is time varying)

The total surface pressure p_sfc_tot(i,j) = p_sfc_dry(i,j) + pressure due to the integrated moisture from the model lid all the way to the surface.

Once you compute your total surface pressure increments, then compute the surface pressure increments from the dry pressure by removing the integrated moisture through the column. Then the time difference of the p_sfc_dry(i,j) is defined as the time difference of the mu(i,j).

 

[mirabelle]

Thanks for your quick response!
Cause I don’t know how to calculate the pressure due to the integrated moisture, so
I am wondering how to calculate the p_sfc_dry with p_sfc_tot.

So, I looked up the code for GSI-ENKF, and found some codes in comGSIv3.7_EnKFv1.3/src/enkf/gridio_wrf.f90.

According to comGSIv3.7_EnKFv1.3/src/enkf/gridio_wrf.f90, I write some code below to calculate mu_new:

do j = j_start,j_end
do i = i_start, i_end
qtot=0.
do k=k_start,k_end
do n = PARAM_FIRST_SCALAR,n_moist
qtot = qtot + moist(i,k,j,n)*(grid%php(i,j,k+1)- grid%php(i,j,k))
enddo
enddo
grid%psfc_new(i,j) = grid%psfc(i,j) + increment_psfc(i,j)
p_sfc_dry(i,j) = (grid%psfc_new(i,j) – grid%ptop(i,j))/qtot + grid%ptop(i,j)
grid%mu_new(i,j)= p_sfc_dry(i,j)- (grid%mub(i,j) +grid%ptop(i,j))
enddo


Is this the right way to calculate mu?

 

[davegill]

Mirabelle,
Since you are talking about "increments", are these modifications for use with the WRF DA system?

Once you modify the perturbation column pressure field, note that the other mass fields will likely be inconsistent with the new mu field (geopotential perturbation, pressure perturbation), and the momentum fields will be inconsistent also. The surface pressure increment is a due to the modifications throughout the entire column.

Take a look at dyn_em/module_initialize_real.F. Look for the string "qtot" after the initial declaration. This section of code shows how to compute the integrated moisture, very similar to what you do. This shows some of the initial efforts to balance the mass fields with each other.

As you mention, when you have an increment to the surface pressure, part of that increment is due to mass of the dry air, and part is due to moisture. You want to attribute only the increments of the surface pressure that are from the dry column air mass to the mu field.

psfc_total = psfc_dry + p_at_sfc_due_to_moisture
psfc_total = mu' + mub + ptop + p_at_sfc_due_to_moisture

mu' = psfc_total - p_at_sfc_due_to_moisture - mub - ptop

If the increments are referred to as the old and new values, then we have

mu'_new = psfc_total_new - p_at_sfc_due_to_moisture_new - mub - ptop
mu'_old = psfc_total_old - p_at_sfc_due_to_moisture_old - mub - ptop

We subtract these two.

mu'_new - mu'_old = ( psfc_total_new - psfc_total_old ) - ( p_at_sfc_due_to_moisture_new - p_at_sfc_due_to_moisture_old )

At each (i,j), the increment of psfc_total is a function only of the perturbation fields: mu' and p_at_sfc_due_to_moisture.

mu_new' = mu_old' + psfc_total_increment - ( p_at_sfc_due_to_moisture_new - p_at_sfc_due_to_moisture_old )

An iterative solution is required for mu_new', since the computation of p_at_sfc_due_to_moisture_new is a function of mu_new'.

We know mu_old', psfc_total_increment, and p_at_sfc_due_to_moisture_old.

Assuming an initial estimate of mu_new' (for example, mu_old') we compute p_at_sfc_due_to_moisture_new, which gives us the next estimate of mu_new'. Iterate until some convergence criteria is met or a max number of iterations is exceeded. This will converge quickly.

 

[mirabelle]

Mr. Dave Gill

I would like to add a code to the WRF model to implement online assimilation.
Now I'm trying to add surface pressure assimilation