Note
Go to the end to download the full example code.
DA with Random Initial Conditions on Soil Layers#
This notebook demonstrates how to create an ensemble of models with varying initial conditions for different soil layers. By introducing randomness into the initial conditions, we can better simulate the natural variability in soil properties and improve the robustness of data assimilation (DA) processes.
Estimated time to run the notebook = 2min
nlay = 6
scenario = {
'per_name':['ic'],
# 'per_type': [[None]*nlay],
'per_type': [['multiplicative']*nlay],
'per_nom':[[-15]*nlay],
# 'per_mean':[[-15]*nlay],
'per_mean':[[1]*nlay],
'per_sigma': [[3.75]*nlay],
'per_bounds': [[None]*nlay],
'sampling_type': [['normal']*nlay],
'transf_type':[[None]*nlay],
'listUpdateParm': ['St. var.'],
'pert_control_name': ['layers']*nlay
}
🏁 Initiate CATHY object
🔄 Update hap.in file
🔄 Update dem_parameters file
🔄 Update dem_parameters file
simuWithDA.update_dem_parameters()
simuWithDA.update_veg_map()
simuWithDA.update_soil()
🔄 Update dem_parameters file
───────────────────────────────────────────────────────────────────────────────── ⚠ warning messages above ⚠ ─────────────────────────────────────────────────────────────────────────────────
The parm dictionnary is empty
Falling back to defaults to update CATHYH
This can have consequences !!
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
🔄 Update parm file
🔄 Update soil
homogeneous soil
simuWithDA.NENS = 3
list_pert = perturbate.perturbate(
simuWithDA,
scenario,
simuWithDA.NENS,
)
var_per_dict_stacked = {}
for dp in list_pert:
var_per_dict_stacked = perturbate_parm(
var_per_dict_stacked,
parm=dp,
type_parm = dp['type_parm'], # can also be VAN GENUCHTEN PARAMETERS
mean = dp['mean'],
sd = dp['sd'],
sampling_type = dp['sampling_type'],
ensemble_size = dp['ensemble_size'], # size of the ensemble
per_type= dp['per_type'],
nlayers = nlay,
savefig= os.path.join(simuWithDA.workdir,
simuWithDA.project_name,
simuWithDA.project_name + dp['savefig'])
)
/home/z0272571a@campus.csic.es/.local/share/mamba/envs/pycathy_doc/lib/python3.10/site-packages/scipy/stats/_qmc.py:993: UserWarning: The balance properties of Sobol' points require n to be a power of 2.
sample = self._random(n, workers=workers)
This in normally directly called when using run_DA_sequential()
simuWithDA.ens_valid = list(np.arange(0, simuWithDA.NENS))
simuWithDA._create_subfolders_ensemble()
simuWithDA.apply_ensemble_updates(var_per_dict_stacked,
var_per_dict_stacked.keys(),
cycle_nb=0
)
😟 processor exe not found
🔄 Update ensemble
🍳 gfortran compilation
👟 Run preprocessor
🔄 Update parm file
🛠 Recompile src files [4s]
🍳 gfortran compilation [10s]
✅ Compilation successful!
👟 Run processor
Skip Error Temporary!
0 25.008730776474092
1 -126.16662132728125
2 -31.615134206831023
3 6.790442321853217
4 -1.6144793526819785
5 -34.78058527019098
Skip Error Temporary!
0 2.7741435322026744
1 -62.769218663656986
2 93.22683631327575
3 -39.967779883780054
4 39.40848050456279
5 -19.26234692399839
Skip Error Temporary!
0 -67.15272217842373
1 -8.133053083959016
2 -36.44551685668179
3 92.55597206106512
4 3.1816863506490147
5 -159.11990461772263
Skip Error Temporary!
0 -153.63633060128902
1 94.71334227962858
2 19.475575039373677
3 -51.33852396164443
4 -31.03546447456177
5 0.592970394058403
Skip Error Temporary!
0 -6.939129091446571
1 -21.73957413123815
2 86.01580694567447
3 -109.47691845030182
4 32.573710495866194
5 -64.67435995437594
Skip Error Temporary!
0 -29.862298552196894
1 -5.360808490294345
2 -119.9088469969322
3 78.82459421601052
4 13.20094466946202
5 -49.479903902695604
Skip Error Temporary!
0 81.1925747323183
1 -47.38667439579466
2 -103.90223420331067
3 13.763643346291243
4 -3.1045962702418577
5 -70.06366537389583
Skip Error Temporary!
0 -2.730640261156089
1 -48.9521472220584
2 123.31298248758921
3 -66.73441709067835
4 48.40560255284125
5 -106.59942718547127
Skip Error Temporary!
0 -55.577866925725885
1 20.105357186903753
2 -49.29911029305937
3 24.717586432514977
4 -5.930753062627336
5 -89.72236564099816
Skip Error Temporary!
0 -14.945661275238624
1 -28.290692516976026
2 -52.965561440748935
3 41.785870712505094
4 65.8565234144734
5 -132.42209093157976
Skip Error Temporary!
0 -19.641343129745074
1 38.916659419453055
2 -105.64783969251673
3 11.806187971552863
4 -70.49553258963262
5 -9.405536154002865
Skip Error Temporary!
0 39.09341572091731
1 -75.79983702332217
2 11.827624891999676
3 -36.53894642478697
4 -120.01805669716391
5 72.80497389825122
Skip Error Temporary!
0 -135.90012335057702
1 59.800852601118315
2 22.470817731368044
3 -40.1262699220791
4 -30.27336338582714
5 -10.553245886812388
Skip Error Temporary!
0 25.763786738658712
1 -22.388662964043103
2 -0.1917478944540546
3 -68.22903779964975
4 4.961774994576973
5 -92.00687903270006
Skip Error Temporary!
0 10.348626591430405
1 -41.11157337307887
2 49.1583506213844
3 -77.84969667770252
4 -32.571649980867136
5 1.8533173401687586
Skip Error Temporary!
0 -82.1552707769422
1 67.58420848905988
2 -8.426810108847336
3 -17.21272222963483
4 -46.543870045091154
5 18.569445290352544
Skip Error Temporary!
0 55.38060971719175
1 -71.21394114277172
2 -14.320510660076629
3 -44.61184873871725
4 5.242631356167449
5 -24.76308339750997
Skip Error Temporary!
0 -135.1871874363317
1 30.105282482065846
2 -18.615621638706354
3 15.468736190498465
4 58.1209154313571
5 -69.13512669932494
pl = pv.Plotter(shape=(1,2))
for i, ensi in enumerate([1,3]):
DApath = f'DA_Ensemble/cathy_{ensi}/vtk/'
path = os.path.join(simuWithDA.workdir,
simuWithDA.project_name,
DApath,
simuWithDA.project_name + '.vtk'
)
pl.subplot(0,i)
CTp.show_vtk(path,
'ic_nodes',
ax=pl,
clim = [-25,-5],
#show_scalar_bar=True,
)
_ = pl.add_legend('')
pl.add_title(f'Ensemble nb:{ensi}')
pl.show()
plot ic_nodes
plot ic_nodes
Total running time of the script: (0 minutes 11.812 seconds)