Model top instability after DART

[Giovanni Conti]

The simulation runs fine, but when we activate DART data assimilation some disturbance is excited and transported on the top of the model and the patterns of temperature, winds are degraded.

This happens also if we try to assimilate observations only in the lower troposphere and limit the impact of the observations to remain inside the troposphere (differences between forecast and analysis is zero above the level we set). Seems we produce some gravity waves, maybe.

We can improve the situation only reducing a lot the physical and dynamical time step ATM_NCPL=288 NSPLIT=256

If we turn off Data Assimilation we can recover a smoother field after 20 days.

I've used all the damping flag, second order divergence, fourth order divergence, fourth order divergence and diffusion. Nothing solved the issues. The second order divergence (default ) gave better results.

Does anyone have any ideas of what might be behind this behavior and any ways to control it?
-------------------------------------
The model is CMCC CM3 based on cesm2.1.1-exp17 (69af836c8a857ccac1b36efc04b0008770e5970d) with a different ocean component (NEMO 3.6).
(Unable to provide manage_externals output right now due to inability to access NCAR svn servers)
Code modification is substitution of NEMO for POP which is not relevant for this simulation which is AMIP.

From README.case
2022-04-21 15:18:35: /users_home/csp/gc02720/cesm2_1_intel20_1/cime/scripts/create_clone --case /users_home/csp/gc02720/cesm-exp/conv/conv_0001 --clone /users_home/csp/gc02720/cesm-exp/conv/case_template
---------------------------------------------------
2022-04-21 15:18:35: Compset longname is HIST_CAM60%WCSC_CLM50%SP_CICE%PRES_DOCN%DOM_MOSART_SGLC_SWAV
---------------------------------------------------
2022-04-21 15:18:35: Compset specification file is /users_home/csp/gc02720/cesm2_1_intel20_1/cime/../components/cam//cime_config/config_compsets.xml
---------------------------------------------------
2022-04-21 15:18:35: Pes specification file is /users_home/csp/gc02720/cesm2_1_intel20_1/cime/../components/cam//cime_config/config_pes.xml
---------------------------------------------------
*** original clone README follows ****
2022-04-21 14:58:42: /users_home/csp/gc02720/cesm2_1_intel20_1/cime/scripts/create_newcase --case /users_home/csp/gc02720/cesm-exp/conv/case_template --machine zeus --res f09_f09_mg17 --project R000 --queue p_short --walltime 2:00 --pecount 36x1 --ninst 1 --compset HIST_CAM60%WCSC_CLM50%SP_CICE%PRES_DOCN%DOM_MOSART_SGLC_SWAV --run-unsupported
---------------------------------------------------
2022-04-21 14:58:42: Compset longname is HIST_CAM60%WCSC_CLM50%SP_CICE%PRES_DOCN%DOM_MOSART_SGLC_SWAV
---------------------------------------------------
2022-04-21 14:58:42: Compset specification file is /users_home/csp/gc02720/cesm2_1_intel20_1/components/cam//cime_config/config_compsets.xml
---------------------------------------------------
2022-04-21 14:58:42: Pes specification file is /users_home/csp/gc02720/cesm2_1_intel20_1/components/cam//cime_config/config_pes.xml
---------------------------------------------------
2022-04-21 14:58:42: Forcing is Historic transient
---------------------------------------------------
2022-04-21 14:58:42: Using None coupler instances
---------------------------------------------------
2022-04-21 14:58:42: Component ATM is CAM cam6 physics:WACCM specified chemistry:
---------------------------------------------------
2022-04-21 14:58:42: ATM_GRID is 0.9x1.25
---------------------------------------------------
2022-04-21 14:58:42: Component LND is clm5.0:Satellite phenology:

This information is from one of the clones used in the data assimilation experiment. Only initial conditions are changed. I do not think the compilers or machine are relevant here because the simulation works when not using data assimilation.

 

[Giovanni Conti]

This is how the Temperature field appears at the WACCM top if we activate the assimilation. We can see the disturbance already after 6H of forecast (so pure evolution after assimilation)

This is the best we obtain decreasing both the physical and dynamical time step, it is better but still with some disturbance (but it take really too long the run with this configuration to be useful)

 

[raeder]

Hi Gio,
how did you create your initial ensemble?

If it's not close to the observations, then there will be large adjustments at the start,
which may propagate upward and grow, as you suggested.
It might help to increase the obs error variance in the early days,
so that they have a smaller impact, but still pull the ensemble toward the obs.
But that involves making modified copies of the observation files.

You might get a similar effect by assimilating a subset of observations,
either by removing some types from the assimilate_these_obs list,
or, again, by using obs_seq files with fewer observations of all types, or selected types.

A third option is to create your initial ensemble by giving DART (filter) a single model state
and allowing it (input.nml:model_nml) to add small perturbations to it for each member.
Then many of the observations will be rejected at the start (due to small ensemble spread)
but the ensemble will gradually be pulled to the observations, a higher fraction of obs
will be assimilated, and eventually you'll have an ensemble that can accept all of
the observations without sending big waves to the model top. That's the theory and hope.
WACCM is sensitive.

 

[Giovanni Conti]

Hi Kevin,
Sorry for the late reply but we had some unscheduled electrical black-out and we did not finished our experiments for the understanding.

Back to your question (how did you create your initial ensemble?):

- At the beginning we created the initial ensemble starting from a previous analysis, that already showed some bubbles at the top, by means
Of perturb_single_inst program, perturbing T and US with perturbation_amplitude=0.2, 0.5.
Witout DA the model is able to remove the bubbles for all the members (after around 20 days), but if we activate the DA then these bad features come back, or are maintained, also if we set
model_damping_ends_at_level=38, so no obs influence above this level. Is the perturbation too big? Perturbing the whole column do we break the hydrostatic balance?
What we do not understand is why removing DA the model is able to adjust the state, considering that also with active DA we do not introduce disturbances above a certain level.


-Then we tried to understand if some set of observations could create some problem. We run experiment
1 conventional (ACARS, AIRCRAFT, RADIOSONDES) +AMSUA
2 conventional+gpsro
3 conventional+satwind
4 conventional
But in all the experiments we maintain the bubbles feature at the top. NCPL=288, NSPLIT=16

-We tried to assimilate only conventional on the first few levels near the surface, in this way the model adjust the top. If we restore the assimilation of conventional until the mid troposphere then the bubbles come back again.


-Then we tried to create other sets of I.C. for new experiments:
1 the first 28 levels, up to 35 hPa, that are the same for CAM6 and WACCM are taken from NCAR reanalysis (for 30 members), above this level all the members are equal to a WACCM climatology. We assimilate everything until 35 hPa and no bubbles are formed at the top of the model. We used NCPL=288, NSPLIT=60…so really slow evolution…
2 the same as above but starting from the same WACCM climatology (land, cice are different for every members as above)
3 the same as above but with 32 level (so basically the IC for the atmosphere is the NCAR reanalysis, the boundary for the 32 levels also required less chemical species)





Top left:
IC differences between the 32 level case and the 70 level case. Clearly the only differences are in a thin band between 35 hPa and 3.6hPa where the WACCM pressure levels and reanalysis are not the same


Top right (1):
6hr zonal mean temperature changes starting from the full column WACCM IC

Bottom left(2):
6hr zonal mean temperature changes starting from the 32 level IC (NCAR reanalysis)



Bottom right(3):
6hr zonal mean temperature changes starting from the hybrid IC

It is interesting that 1 and 3 in just 6 hr looks complementary at the top, and further there is this huge intrusion in 3 that start above the SH.

We started to do some test also trying to turn off the gw, just to see if they have some effect. However setting to false
use_gw_convect_dp, and use_gw_front it gives some error, because the model requires some variable that we did not provided. Then we tried setting to zero the efficiency effgw_cm and effgw_beres_dp. We are analyzing the results now.

 

[raeder]

The largest T perturbation I've used is 0.1 K, and I haven't perturbed US.
My feeling is that your perturbations are larger than is ideal,
but I can't say that they are definitely a problem.
The goal is to use perturbations that are small enough to prevent model instabilities
but large enough that it doesn't take weeks to spin up. Since WACCM is sensitive,
it could be worth trying smaller, if you haven't tried that yet.

"At the beginning"; am I right that you started from a single analysis model state
(e.g. an ensemble mean) and not from an ensemble?
It's highly likely that WACCM was developed to do 2 things; not create
unrealistic fluctuations and to damp out any that occur,
using the minimum damping that will do the job.
The DA increments seem to be frequently pushing the model beyond its ability to
damp the non-physical aspects of the increments.
Even when the increments are only in the troposphere, the small gravity waves
that they create propagate upward and grow in amplitude as the density drops.

That's an interesting result, that assimilating obs only near the surface
doesn't lead to problems in the upper atmosphere. That might be because
there are small numbers of observations there, so the model state is not changed much.
You can check that by comparing the increments (analysis - hindcast) there
with the increments in the mid and upper troposphere in the case where you
assimilate all tropospheric obs.
That would be consistent with my idea (previous post) that assimilating more gently
could reduce the high fluctuations.

To be sure that we're looking at the same atmosphere, are the levels in the pictures
labeled with 1 near the surface? That's the opposite of how WACCM labels them.
And does "6hr change" include assimilation at hour 6?
And it's a single, 6 hour, assimilation cycle?

In the upper left figure, those are big, spatially "sudden" differences in ICs.
Does the color scale even show how big they are? They seem to be saturated at
the maximum values on the scale.
But the sudden transition from CAM analysis to WACCM climatology around level 28
doesn't seem to cause a problem, when assimilation is done below level 28,
so that's interesting (and hopeful?).

I'm suspicious about the lower left figure. My interpretation is that
above level 32 the T changes during 6 hours are all < 0.1 K, which seems unlikely.
It's especially unlikely because the WACCM climatology is not a natural WACCM
model state, so it should be somewhat unbalanced, which would lead to rapid changes.
But this is the zonal mean, which I don't look at much, so maybe it's possible.

That's an excellent observation, that the 2 pictures on the right are so complementary.
I assume that you've checked that the differences were calculated in the same order
(hour 6 - hour 0). The upper atmosphere is not affected by the assimilation at hour 6,
so those complementary patterns must be the result of the different ICs below level 32;
WACCM climatology versus CAM analysis. In a nonlinear system I don't know how we could
even construct ICs that evolve to opposite patterns like that, much less have that happen
by accident.

Am I right in thinking that the lower right figure corresponds to
"other sets of I.C." item 1? I can't see from the zonal mean whether there are
or aren't "bubbles" at the top. If there aren't in the lower right figure,
then there aren't in the upper right figure either, since they have essentially
the same pattern.

 

[Giovanni Conti]

You right Kevin the levels are labeled in the reverse order in the picture above.
Yes that was a single assimilation cycle.

We have some new update about our problems with WACCM+DART.
We created an ensemble extracting 30 member of NCAR reanalysis and then above 35 hPa, where the CAM6 and WACCM levels start to be different, we just plug the WACCM climatology.
After a few weeks of free run by using fv_div24del2flag=2 (default) NCPL=288 NSPLIT=64, NSPLTRAC=16,NSPLTVRM=16 we started to assimilate observations up to ~2 hPa and damping the DA correction above 0.7 hPa. With this setup, although really slow, the assimilation of all the NCAR reanalysis observation +AMSUA (aqua satellite ch 8,9,10,11,12,13,14) does not disturb as before the top that remain smooth (not 100% perfect, there is still some little disturbance as in the second plot of the second post, but way better than before).

The problem is that some ripples at the poles, in particular in the South Pole and along the Antartic coast, started to raise in the VS field and after some days also temperature, zonal wind and surface pressure where affected by that. After that the meridional velocity at the south pole suddenly diverged (400 m/s) and the simulation crashed. In order to get rid of the ripple we used fv_div24del2flag=42 and the laplacian diffusion removed completely the problem at the surface . . . but the top of the model started to behave bad again (if the assimilation is not turned off).
Furthermore we started to have some warning about nucleation and optical depth unreasonably too high..and model crashed again.


So, we were wondering if in order to maintain some internal equilibrium in WACCM is not possible to have an hybrid damping e.g. 0.5*div4del2 + 0.5*div2 in a way that we can cure at the same time the top and the bottom. In the code I see a comment that says del2 and div2 cannot coexist.

Do you have any suggestions about this?


Thanks a lot

Cheers
/gio