OGGM / oggm

Open Global Glacier Model
http://oggm.org
BSD 3-Clause "New" or "Revised" License
223 stars 105 forks source link

Computing Glacier velocity due to ice-deformation #1018

Open ShashankBice opened 4 years ago

ShashankBice commented 4 years ago

I was wondering if there was a direct way one could derive velocity expression via deformation using given estimates of ice-thickness, surface elevation and creep factor A ? And if one could model this over time-steps (how glacier velocity due to deformation evolves) given a mass balance profile ?

fmaussion commented 4 years ago

Thanks for the question!

We output surface velocities of our steady state glacier after inversion if the user wants it, with:

https://github.com/OGGM/oggm/blob/3750ee59fedc1658abf6f22358ab0306b4b520d3/oggm/core/inversion.py#L512

However, this is computed from our estimate of the ice thickness given a mass-balance profile. I don't really know what your question would imply if the ice thickness and the mass-balance profile you want to use are actually from different sources, i.e. non physically consistent.

Regarding the transient runs: we also have a way to compute ice velocities at any time in a model run: https://github.com/OGGM/oggm/blob/74eb713856b64e90e61d65580eb0400190c67776/oggm/core/flowline.py#L1479

But this is not written out per default, but we could add an option for this.

ShashankBice commented 4 years ago

Hi Fabien, Thanks for the reply and thanks for putting up the edu-section. While running the toymodel notebook (oggm-edu) for modelling glacier mass balance, I am able to retrieve the dataframe containing the fluxes and velocities. I was trying to run the compute_velocities function midway in the oggm-edu tutorial for glacier-water-resources. I setup the gdir with preprocessing level of 3, and then ran till this cell


tasks.run_random_climate(gdir, 
                         y0=2000, # climate centered around 2000 (1985-2015)
                         unique_samples=True,  # see documentation (type `tasks.run_random_climate?` in a cell)
                         nyears=300, # simulation time
                         temperature_bias=temp_bias, # add a temperature bias
                         output_filesuffix=file_id, # an identifier for the output file to read it later 
                         store_monthly_step=True, # the default is to store only yearly values
                         mb_elev_feedback='monthly', # the default is to change the MB every year 
                         seed=0);  # this is for predictability of randomness - it's not necessary though

At this point, I ran this cell

tasks.prepare_for_inversion(gdir,add_debug_var=True)

and then tried the compute_velocities function as: compute_velocities(gdir,filesuffix=file_id). I get an keyerror for flux,the exact error being:

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
<ipython-input-23-990de90617d3> in <module>
----> 1 compute_velocities(gdir,filesuffix=file_id)

/srv/conda/envs/notebook/lib/python3.7/site-packages/oggm/utils/_workflow.py in _entity_task(gdir, reset, print_log, **kwargs)
    458                     signal.alarm(cfg.PARAMS['task_timeout'])
    459                 ex_t = time.time()
--> 460                 out = task_func(gdir, **kwargs)
    461                 ex_t = time.time() - ex_t
    462                 if cfg.PARAMS['task_timeout'] > 0:

/srv/conda/envs/notebook/lib/python3.7/site-packages/oggm/core/inversion.py in compute_velocities(gdir, glen_a, fs, filesuffix)
    553 
    554         # this flux is in m3 per second
--> 555         flux = cl['flux']
    556         angle = cl['slope_angle']
    557         thick = cl['thick']

KeyError: 'flux'

Can you please let me know where I am making the obvious mistake ?

I understand what you mean by physically consistent thickness and mass balance profiles. I will use the ice-thickness generated by oggm with whatever mass-balance curve i input, as I play with the functions.

Thanks !

fmaussion commented 4 years ago

Thanks for the question! I think you need to re-compute the ice thickness as well. i.e:

tasks.prepare_for_inversion(gdir, add_debug_var=True)
tasks.mass_conservation_inversion(gdir)
compute_velocities(gdir)

Here are some more information about the inversion in OGGM that you might find useful as well:

Thanks for using OGGM! We hope it will be useful to you

ShashankBice commented 4 years ago

Thanks a lot Fabien ! mass_conservation_inversion did it !

ShashankBice commented 4 years ago

Just 1 more question, I see that there is an option to input user-defined DEMs, can we enter a user-defined thickness grid, if am interested in only playing with surface velocities ? In the code, I see that if the value of fs>0, the velocity components are computed using the Glen's flow law and the flux values are not taken into account. If there is not a built-in method for this, I can try to mimic the inversion_output.pkl with thickness values from other sources and compute the velocities accordingly.. Thanks again and have a great weekend !