Adjoint models

[RPrudden]

Has anyone given any thought to calculating the adjoints of CliMT models? As far as I can tell there is no easy way to do this, but I'd be delighted to be wrong.

[mcgibbon]

I don't know much about adjoints, but if calculating an adjoint can be done by computing the outputs of sympl components with slight perturbations to the inputs about some reference state, then I would say there should be an easy way to do it, but there isn't an existing way to do it. I don't know that anyone has done this before, maybe @JoyMonteiro has some insight. I think the result would make a neat Python package!

[JoyMonteiro]

@RPrudden you are correct that this is a non-trivial task. A lot of code in climt is in Fortran, and
there aren't really any good free tools (that I know of) that can do automatic differentiation (AD) of Fortran code.

Even if we did get our hands on code like this, ADs have a hard time dealing with model physics because of strong nonlinearities present in the atmosphere in the presence of clouds and water vapour. See this paper for example.

While it would certainly be interesting to see adjoint models for climt components, it is not on my priority list currently.

If you are interested, I think the dynamical core would be a good target to develop an adjoint for.
The core dynamics is contained in a few files that make it quite amenable for such an effort.

[RPrudden]

Many thanks, both! I think computing the adjoint of the dynamical core would be enough to cover the kind of work I have in mind, so I may be able to avoid the difficulties involved in differentiating the physics components.

In terms of the dynamical core, it looks like this directory would be the place to start. I guess you were thinking to use autograd to compute the gradients, is that right?

[JoyMonteiro]

That directory only contains the wrapper. The actual algorithmic part is written in Fortran and resides
here

The easiest way I would think is to translate the relevant parts to python code and then run autograd on it. There are only two core functions:

• advance
  

  climt/climt/_lib/GFS/run.f90

  Line 72 in ddc6530

  subroutine advance(t)

which does the actual time stepping and

• getdyntend
  

  climt/climt/_lib/GFS/dyn_run.f90

  Line 140 in ddc6530

  subroutine getdyntend(dvrtspecdt,ddivspecdt,dvirtempspecdt,&

which calculates the tendencies used by the timestepper. Each will need some helper functions,
but these are your entry points to the code.

There is a subtle issue to keep in mind -- There are parts of the code which do differentiation themselves, calling the shtns library. These will have to be translated by hand to calls to differentiation by shtns (there is a python wrapper to shtns, so not a big issue there).

I think it is doable, but requires some grunt work!

[RPrudden]

Ok, I understand better now. As you say, it looks fairly intensive to translate.. but I am quite interested, so I might try to carve out some time.

[JoyMonteiro]

Great!! We can keep this issue open in case you need any other input. While I can't contribute to the coding effort itself, I can certainly take a look at any code or help you to navigate the fortran code.
It would be really cool if this works out :)

[JoyMonteiro]

A bit of googling around yielded a AD package for Fortran:
http://www-sop.inria.fr/tropics/tapenade.html
which seems to be actively developed.

They seem to have used it for an ocean model (NEMO).