arbitrary tracer support

[seleneonowe]

Hello! I've just recently started playing around with speedyweather and I love it, definitely so far accomplishing its goals of being easy to use compared to other models I have experience with.

I'm interested in using it to model atmospheres of other solar system bodies, which ties nicely into your to-do list item of wanting an "exoplanet mode" and #462.

For other worlds, different physical processes dominate and so differences in parameterisation schemes are essential. That's fine, as those are generally going to be fine-tuned by and/or written entirely by the user.

Speedy currently has support for a single tracer in the PrimitiveWetModel, water vapour by default - it isn't currently documented, but its properties can be changed by the user. This was discussed in #460 - also mentioned there was that the computational expense of adding new tracers making you hesitant to add support for new ones (in the context of #460, the new ones being water condensates). We can see from the benchmarks the performance hit you get just going from PrimitiveDry to PrimitiveWet.

Tracers in other atmospheres

I'd have to argue from the point of view of generalising to other atmospheres, that support for users adding their own tracers is essential - they should of course be warned about the computational hit of new prognostics, but if this is to be a sandbox for playing with a wide variety of atmospheres rather than just Earth's, then there's really no getting around this. You can just about get away with one tracer on Earth as a first approximation, but on other worlds, that's not necessarily the case.

Mars

At the absolute bare minimum for simulating the Martian atmosphere, I would say you need two tracers: one to keep track of water ice, and the other dust. You'd have to then crudely diagnose the number of particles from only the mass, or vice versa. A slower but more accurate model might have four tracers: mass and number mixing ratios (mmr and nmr) for both, and then parameterisations could assume some properties of the distribution of particle sizes given both of these numbers.

You can't readily neglect either ice or dust, nor can you diagnose one from the other.

NASA's legacy model used two additional tracers (water vapour mmr and dust core mmr) for a total of six.

Summary: definitely can't have less than 2, 4 would be less crude, 6 is great.

Venus

Has complicated chemistry and I'm no Venusian expert. Both water and sulfur based stuff is important and you can't diagnose one from the other, so that's already multiple desirable tracers.

I've seen simplifying assumptions that certain tracers only exist at fixed altitudes and then advecting them in 2d rather than 3d.

I'm not sure what a reasonable number of tracers is or a perfect number, but I suspect one is far too much of an oversimplification.

Titan

Titan has some of the most complex atmospheric chemistry in the solar system. The most important tracer is methane, which can also exist in several phases just like water on earth, so wanting three tracers for it is a very reasonable request, although much like water, you could massively simplify and treat the condensates as diagnostics of the vapour.

You'd also want ethane as a tracer, though, it's also very important.

There is also a haze that is typically treated as a single homogeneous thing rather than a soup of many organic compounds.

You could maybe get away with just two, and diagnose everything from those with some other assumptions about droplet formation etc. I've seen the haze often being treated as a constant or even globally uniform field.

Ideally though the choice of tracer number would be user-specified not baked into the model, since "how few tracers can we get away with while still reasonably approximating the Titanaen atmosphere?" is an interesting research question that's arguably still an open problem; and speedy would be a great tool for exploring that.

Jupiter and Saturn

From what I know about them, personally I'd want to advect two tracers if I had to pick a minimum number: one to represent ammonia and the other methane.

Uranus and Neptune

These you could get away with only one tracer I think? Just methane would probably get you pretty far in simulating these planets.

Triton and Pluto

It's not wildly unreasonable you'd want more than one. There's nitrogen ice but you could diagnose that by altitude maybe. There's methane and carbon monoxide you might want to advect.

Summary

Even just looking at the atmospheres we know a lot about, without going to exoplanet research, support for user-defined tracers is important if we want speedy to be generalised to other atmospheres. I personally feel it's the main thing it's missing - people studying a new atmosphere might well expect they have to write an entirely new paramterisation scheme, but not that they'd have to also know how to modify the dynamical equations to account for new tracers.

Implementation

Not coming from a dynamics background myself I'm not sure exactly how implementation of new tracers should be done. Correct me if this is nonsense but I think "all" that needs doing is:

• add a new tracer transport equation to the primitive equations, with the time derivative of the tracer being the negative wind dot grad tracer (plus some source and sink terms from parameterisations).

• solver would have to be modified to loop over all tracers

• tracers could be passive (in which case that's it) or active (which I'm using to mean the effect of the tracer on density is not negligible) in the latter case the virtual temperature would need to be modified to be the temperature * (1+ sum_over_active_tracers((gas_constant_ratio_for_tracer_i - 1)*qi)).

User interface

Very much up for discussion, but I was envisioning something like having a PrimitiveExoticModel (unsure of naming) in addition to the zero and one tracer Dry and Wet models.

This could be created with a 1d array of AbstractTracer objects.

Each abstract tracer object could be required to have a name (for reference in setting initial conditions or writing output or something), and be either passive or active, and if active then it would also need a gas constant.

The ratio of tracer gas constants to the bulk atmospheric gas constant (dry air on Earth) could then be computed and stored internally since those are the important numbers in the dynamics (I think).

Motivation

This does seem to square with the goals of speedy, since a desire for an exoplanet mode is mentioned on the to-do list and it definitely leans into the "play atmospheric modelling like it's LEGO" desire for flexibility. Unfortunately approximating the atmosphere with a single tracer only sounds reasonable for a subset of known planets, and I think avoiding being too Earth-centric is a noble goal.

This would also assist with the statement about extending speedyweather with new parameterisations: the "we're happy to provide an interface for a complicated cloud scheme you want to attach to speedy, even if it'll live in its own repository rather than here because it's too slow and complex" statement in speedyweather's docs does rather imply that it'd be nice to just say extend AbstractCloudScheme, as opposed to the cloud parameterisation developer having to do some complicated implementation of new prognostics to the dynamical core themselves (a much higher barrier of entry since the people who study one of physics or dynamics tend to be much less familiar with the other).

As an aside, many of these worlds are smaller than Earth, so the same angular grid spacing gives you a higher effective horizontal resolution. You can make up a lot of lost performance from new tracers by the performance gain of simulating a smaller world's atmosphere. Of course the opposite is true for very large worlds unfortunately, but I'm quite interested in the smaller worlds myself!

I am willing to write code to help with this where necessary, though I suspect I'd be more useful helping write some parameterisation implementations than updating the dynamics.

[milankl]

I definitely don't have to be convinced that this is high on the todo list. I was basically just waiting for someone like you laying out the physics that would be needed to be reasonably general. Solar system planets/moons sounds like a great starting point. So many thanks for this overview.

I draft a more specific list that we would need

• SpectralGrid should get a ntracers (default 0) field, as this determines the size to preallocate a bunch of arrays within prognostic/diagnostic variables including the vertical column for parameterizations.
• every model should get a tracer_advection field (default NoTracerAdvection) that contains vectors to denote for every tracer whether it's active or not, name, gas constant and we might even be able to specify 2D advection only (default for our 2D models) which could possibly be with the vertical mean velocity or on a given layer.
• as the last step in the dynamics advect the tracer by looping over 1:ntracers
• every advected quantity requires some diffusion, I hope the default diffusion like we do for humidity is sufficient. It should but we'd need to check that.

I currently have a gut feeling that I would like to keep the tracer advection separate from the dry/wet model distinction because a dry model would also automatically switch off a bunch of parameterizations that are humidity specific and so one can better control "tracers on Earth" (use wet) and "tracers on exoplanets" (use dry). I'm happy to think more whether humidity should always be treated as another tracer or whether it makes sense to let it have its special role. From a physical point of view, yes just a tracer, but this question is more about user interface. E.g. if humidity is the first tracer how does the model know to automatically switch on default parameterizations for the hydrological cycle (wet convection instead of dry, precip, evaporation)? Just by the name of the tracer?
In that sense I currently have the gut feeling that dry/wet should stay and tracer advection is independent of that. But happy to discuss!

I do want to wait for #525 to be merged before actually starting this. But I believe then it's pretty straight forward to get passive tracers without parameterizations in. We would then also need to discuss how to generalise the parameterizations for all those tracers. So question for you, what would be the parameterizations for tracers on all of these planets be? Would a simple tracer_parameterization! function that can then be extended, suffice?

[seleneonowe]

Humid processes

if humidity is the first tracer how does the model know to automatically switch on default parameterizations for the hydrological cycle (wet convection instead of dry, precip, evaporation)? Just by the name of the tracer?
In that sense I currently have the gut feeling that dry/wet should stay and tracer advection is independent of that.

So, regarding default parameterisations, there are actually some "wet" tracers in other planets that mimic the hydrological cycle on Earth that require very similar parameterisation schemes. For example, on Titan you have methane clouds, rain, hail, evaporation, great lakes, surface runoff, etc. People have had success in using schemes we use for water on Earth for methane on Titan with appropriate parameters adjusted (e.g. in the critical relative humidity calculation, maximum droplet sizes, radiation interactions etc) - the reality is quite complicated, the presence of ethane and other factors have big effects, but as a nice simple first approximation you can just say "relative humidity above this threshold produces large scale precipitation" just like you can on Earth.

So it is possible that we might want some appropriate defaults in say the large scale condensation and moist convection schemes for different wet tracers. Regardless of whether we modified the dry/wet distinction in the UI, I think we could just adjust our existing moist schemes to offer a foundation for handling other "wet" tracers.

Most tracers are dry, though.

Dry processes

Would a simple tracer_parameterization! function that can then be extended, suffice?

Yes, I think so. It could always be split up further into new more specific abstract parameterizations where necessary.

Of those dry tracers, things you might want to do with them that I can think of off the top of my head just now:

• deposition removing the tracer from the atmosphere and potentially adding it to the surface if we want to track that (for soluble tracers, this might be split between dry and wet deposition, where e.g. rainfall from other schemes acts as a sink for the tracer). These could extend a tracer_parameterization function just fine I think, although there are bits that you might want to couple to moist schemes and land surface schemes, and I don't know what the ability to pass information between schemes currently is like.

• one of the most likely and most important things our tracers are likely to do is interact with radiation (for example Martian dust absorbs a lot of sunlight; a transparent shortwave scheme on Venus would be a terrible approximation and there's a very strong greenhouse effect for the longwave to deal with; likewise on Titan only about 1% of the top of atmosphere radiation actually reaches the ground). Not sure how to handle that in terms of extending radiation schemes because I'm aware we only have a transparent shortwave scheme right now and a simple longwave one. Since the planned next steps are for a still-simple-but-less-simple radiation scheme that takes gases and albedo into account in some way however idealised, it could well be possible to write such a scheme in a generic fashion, I'll have a think about that. Since more complex radiation schemes are pretty slow there's an argument to be made that keeping complicated schemes planet-specific so you can make as many simplifying assumptions as possible is best; but such things are getting a bit outside the scope of what would be in the main SpeedyWeather repository anyway. Regardless of the nature of such schemes or where they sit, it would be more appropriate to extend the abstract radiation scheme for this than to extend a generic tracer parameterisation. An abstract radiation scheme should thus have access to the tracers somehow.

• interaction with other tracers. This could indeed just extend a tracer_parameterization, and could be a diverse range of processes. Random example, some chemical process where tracer 1 and 2 react together to produce tracer 3 under high enough pressure could be parameterised here.

• sublimation scheme (similar to a surface evaporation scheme) that moves surface ices into atmospheric vapour. Could be generic too.

Overview

We would then also need to discuss how to generalise the parameterizations for all those tracers. So question for you, what would be the parameterizations for tracers on all of these planets be?

There are loads of possibilities. I could go into detail about the best of my knowledge on the schemes needed for specific planets if that would be useful. Being more generic, I feel like there are several general cases that can happen:

• an existing scheme works just as well on the new planet as it does on Earth. No action required.
• an existing scheme works fine on another planet with minor adjustments to account for the different tracers. Straightforward to deal with and could be done in speedyweather.
• existing processes still happen on the new world, but the effect of the new tracers makes them too different to Earth and a new version of an existing type of scheme needs to be made (e.g. a new radiation scheme, a boundary layer scheme that also transports dust from the surface into the atmosphere). Depending on complexity and generalisability, may or may not be appropriate to put in speedyweather.
• an entirely new type of scheme we don't currently have needs written (e.g. dry deposition, submilation of ices; in-atmosphere chemical processes). These could extend a generic tracer parameterisation.
• schemes need to pass information between each other that they didn't previously need to do (e.g. radiation scheme needs to know diagnosed particle size distribution; surface scheme needs to know rainfall or dry deposition rates). Not sure how easily they can currently do that.

As an aside, I know some complicated chemistry and aerosol models that want to keep track of very large numbers of tracers are ran separately at a lower resolution to the main model (because otherwise the number of tracers quickly becomes the limiting factor in runtime, even though it's really not that important to know all of their distributions to such a high resolution) and have a coupling layer to fetch the winds and send information back. Such a model is beyond the scope of what I'm thinking of (I'm considering scenarios where you want two tracers, or eight tracers, not fifty).

Other

every advected quantity requires some diffusion, I hope the default diffusion like we do for humidity is sufficient. It should but we'd need to check that.

Yes, reasonable thing to check. I feel like it's probably okay? Not too sure myself.

Question unrelated to tracers (but still on the topic of simulating other planets)

There is another thing I've thought of which is outside the scope of this issue but I wanted to ask about since my knowledge of how the dynamics are actually solved is quite basic. In speedy we assume an ideal gas law which is actually still fine across a wide range of planetary conditions, even near the Venusian surface I think it's a good enough approximation. It obviously breaks down deeper in the giant planets but at some point you have to give up and couple to a (non-stationary) "surface" anyway. Using the van der Waals equation for your equation of state instead would permit the simulation of very extreme high pressure and low temperature conditions, but it's nonlinear and thus would essentially mean a rewritten (and slower) dynamical core. Since you can simulate such a wide range of environments with just the ideal gas law, and full van der Waals makes things so much more complicated, I wouldn't say there's any urgent need for such a thing.

A question I'm passingly interested in though (just thinking about generalisability to the widest range of planetary environments possible, without modifications so huge they require a totally different model), how possible would it be in the current primitive equations implementation to introduce an altitude-varying gas constant? An approximation between an ideal gas and the van der Waals equation is to introduce a compressibility factor into the ideal gas law, a constant whose deviation from 1 represents the deviation from ideal gas behaviour (higher than 1 where repulsive forces are significant due to a higher collision rate, less than 1 where attractive forces are significant). The value of the constant depends on the conditions, but in GFD you can use height as a proxy for temperature and pressure and set a fixed compressibility factor at each level in an ad-hoc manner as part of your model initialisation. You could then swap your gas constant in the equation of state for a vector of gas constants times these numbers, which I conceptualise as an altitude-varying gas constant. I ask because I'm not actually sure if that's a minor change or a massive change. Depends how much is resting on the assumption that the equation of state is the same everywhere.

[milankl]

@seleneonowe sorry that this took so long but #525 was holding us back. Now this is merged and I've started drafting arbitrary tracer support in #579

[milankl]

We would like to have an interface that's as simple as

add!(simulation, tracer=:co2)
set!(simulation, co2 = (λ,φ) -> φ > 0 ? 1 : 0)

and that allows one to build parameterizations that use or change the tracer concentration, e.g. a tracer source like wildfires, a tracer sink like radioactive decay or an active tracer with a concentration that's used for radiation/precipitation etc.

[milankl]

With the latest commits in #579 the interface would look like

add!(model, Tracer(:co2))

and ntracers = length(model.tracers). The tracers can be added to a model but this has to happen before initialize! to a simulation. That way model.tracer[:co2] is a struct identifying a tracer by its key :co2 and assigns it an id (counts up from 1). simulation.diagnostic_variables.grid.tracers_grid[1]. Alternatively I could see implement tracers as a Dict{Symbol, Spectral/GridVariable2/3D} wherever they are needed, to always identify the through key and not through key that's mapped to an id::Int.

[seleneonowe]

Alternatively I could see implement tracers as a Dict{Symbol, Spectral/GridVariable2/3D} wherever they are needed, to always identify the through key and not through key that's mapped to an id::Int.

I'd personally lean towards having the option to use keys wherever possible (e.g. in the last example something like simulation.diagnostic_variables.grid.tracers_grid[:co2]). In most circumstances the number of tracers being used is likely to be low (if any are in use at all), so it'd typically not be much of a bother to mentally keep track of say which of your two tracers is 1 and which 2, but in general I'd say code would be easier to read and write if keys are being used.

Whatever implementation details you decide to go with, progress so far sounds good!