[MITgcm-devel] new pkg/matrix checked in

[samar khatiwala]

Hi Jean-Michel

Thanks for taking a look. I can implement most of your suggestions. See my
detailed responses below.

> a) PARAMS.h is not up-to-date (Ed removed seaice_mdsio)

It is hard to keep up :-)
(Its one reason why I am desperate to have this pkg incorporated into the
MITgcm tree.)

>
> c) the other thing is that I don't see clearly why you need to
>  skip the ptracers output (in do_the_model_io) and not
>  just set the proper output frequency to zero.

Because the ptracers output is completely useless when using the
matrix pkg so I just want to avoid the overhead. But now that you mention
it I see there is a variable PTRACERS_dumpFreq that could also be set
thus avoiding changes to the ptracers code. I could simply set this
variable to 0.0 when initializing the matrix package. Would you prefer
I do that?

>
> d) also for ptracers_init.F, what is the reason for setting
>  tIter0=nIter0 and how can we do a correct restart with this ?
>  And what about having tIter0 in data.ptracers ? will this
>  fix your problem ?

The current implementation is such that if the gchem pkg is NOT being
used, tIter0 is set to 0. Thus, if one starts from a pickup file
(nIter0 > 0 in 'data'), ptracers_init attempts to read a pickup_ptracers
file. I have never understood the reason for having this as a default
(and recall raising this issue several months ago). In many (most?)
cases you want to restart the dynamics from a pickup file, but start
the tracer calculation from some prescribed initial condition (if
its an initial value problem), or from zero (if its a forced problem).
The matrix pkg is doing an initial value problem, so I absolutely
don't want to read a ptracers_pickup file. I want to start from a
prescribed initial condition. So I had to change ptracers_init.
Note that the "if (useMATRIX) then" will ensure that ptracers are
correctly restarted when the matrix pkg is NOT being used. That
is, if useMATRIX=.FALSE., the model does exactly what it did before.

Your suggestion of putting tIter0 in data.ptracers would probably work.

>
> e) ptracers_forcing_surf.F : is there a general reason to
>  leave surfaceForcingPtr=0 or is it simply to avoid getting the
>  salinity forcing ? (and in this 2nd case, I would keep this file
>  in matrix_example/code)
>

I want to avoid getting the salinity forcing. In fact, when using the
matrix pkg I don't want any forcing on the ptracers.
I can keep it in my mods directory if you like. Again, IMO its a bit of an
odd default to use the salinity forcing. Why not just set it to zero by
default?

> f) regarding the 3 gad_dst3_adv_?.F , I would leave them in
>  matrix_example/code for now (I would like to have Patrick's
>  feedback on this piece of code)
>

Sure. I also have a question for patrick about this. My understanding is
that dst3 (#30) is a LINEAR advection scheme (we had a long discussion on
the listserv about this last year). That is (unlike dst3fl (#33)) it is
not flux limited. My experience however is that, as currently implemented,
it is not quite linear in the tracer field. (How I arrived at this
conclusion is a long story which I will be happy to relate if anyone is
interested.) The changes I made (and which are only executed if
useMATRIX=.TRUE) ensure linearity. If you look at my bit of the code, you
will notice that it is exactly what appears in the documentation of the
MIT model (without the psi* variables).

Samar