AtChem/AtChem2

Constraining Groups of Compounds

AlfredMayhew opened this issue · 5 comments

Hi,
I have a question regarding the possibility of constraining (or holding constant) groups of compounds in a model run. As an example, I'm looking to run a model in which I hold the concentration of NOx constant but would like to allow NO and NO2 to freely vary within that constraint. So for every time step, NO and NO2 can take any value provided the sum of NO and NO2 (the total NOx) is equal to the constrained value. There are some analogues to the RO2 sum here, except I'm looking to constrain the compound group, whether its NOx or RO2.
I can imagine having a new file in the configuration directory where groups of compounds are listed and then these groups could be referenced throughout the model, or just in the constraint files. This could open up a lot of possibilities for what can be done with the models, but I have no idea if this is easily implemented without requiring major changes to the code.

Do you have any thoughts on whether I would be able to adjust the model code to allow for the constraining of groups of compounds in models? And if so, how would I go about this?

I realise this is a very speculative request and so it's absolutely no problem if this would require far too much work, but I figured it was worth asking in case there is some reasonably quick way that this could be implemented.

Thanks,
Alfie

rs028 commented

HI Alfie, I do understand what you mean and I agree it could be a useful option to have.
I am not really sure how to go about it to be honest and I would imagine that it requires some significant coding. Happy to assist if you feel like giving a try! :)

Hi,
Thanks for your reply. No worries, maybe I'll have a crack at it at some point when I'm feeling brave!

Thanks for your help.

As an update to this, I have managed to get a NOx constraint to work by making a few changes. I'll outline what I've done below, as well as the current limitations. I may try making this more generalisable at some point, but this fulfils my needs for now so I'll just put it out there for anyone else who is interested. The changes can be seen on my AtChem2 fork.

  • I created a new function in dataStructures.f90 to return the number corresponding to any species so I can find NO and NO2 easily.
  • The desired NOx concentration is given in the environmentVariables.config file as a new variable named NOx.
  • At every time step when constrained species are updated in the addConstrainedSpeciesToProbSpec subroutine in constraintFunctions.f90, I also check the NOx concentration (sum of NO and NO2). The ratio of NO and NO2 is kept constant while the values are scaled such that the NOx equals the desired value.

As mentioned before, this isn't the best solution and there's a few things that could be changed:

  • This should be made more generalisable, with the ability to specify any number of species groups containing any number of species. It would be best to have a new configuration file in which the species groups are listed.
  • My current implementation now requires NOx to be constrained. Using a NOTUSED value for NOx will set the desired NOx to -1 and is equivalent to removing all nox from the model. Clearly the user should be able to specify whether groups are used or not. (I think that setting NOx to CONSTRAINED does work here though).
  • I have noticed model runs crashing more often than usual, though I suspect this may be because I'm trying to force the chemistry to do unusual things with my new constraints. I'll keep an eye on this and update here if I find a deeper reason for this.

Hopefully that's useful and I'll post again if I do any further work on this.

rs028 commented

Nice work. Thanks for posting this. As you say, it requires a bit more work before we can think of merging it, but it is good to have a workable starting point. Please update if you test it more and/or develop further.

As an update to this, I've developed the branch a bit more to correct the constraining of NOx from a fractional difference (which would lead to run-away partitioning of NOx towards 100% of either NO or NO2) to calculating the absolute difference and apportion it between NO and NO2.
I also found quite a major issue with my modification where using constrained species (since I'm editing the array of unconstrained species using the index of compounds in the array of all species). I'm getting around this by not using constraints, but there still seems to be some odd behaviour which I'm now looking into.
Basically, it looks like this might suit my purposes for now but I wouldn't recommend anyone uses my changes! Just thought I should put this here since I've previously linked my changes!