Introducing Eulerian tracer during a simulation - flushing time

[marcmestres]

Hi all,

I am trying to calculate water e-flushing times in a semienclosed basin according to Jouon et al.'s (2006) definition, which involves seeding the wet nodes of the domain with a Eulerian tracer and seeing how it evolves in time. Everything is fine if I want to introduce the tracer at t = 0, since all I have to do us define it in the initial conditions file. However, in this particular case I want to introduce the tracer at time t = 1 day. Is there a simple way of doing this without using a restart file, or do I need to "customize" the ROMS code? I will be doing this with ROMS version 3.8.

Thanks a lot.

M

[kate]

I'd think it would be easier to use a restart file. Would this be a problem for some reason?

[marcmestres]

Thanks Kate.
Yes, it would. I am adding the flushing time component to an already established operational system that runs four-day long simulations (i.e., one day hindcast+three days forecast) in a large number of domains. The flushing time calculation will only apply to one of these domains. Modifying the existing operational procedure is out of my reach, and I am fairly sure it's not even an option. This is why I am looking for other ways that do not involve the use of the restart file.

M

[wilkin]

I can't think of any way to initialize passive tracers other than via the initial conditions file ... so here a restart.

But you can have floats initialize with a time delay, and you could estimate flushing time from a large float ensemble.

[marcmestres]

John,

Thank you, but if possible I'd rather not use Lagrangian floats. I want to obtain a spatial distribution of the flushing time in the waterbody, and I think a Eulerian approach is better suited for this, Searching through the forum and the cppdefs.h file I see that perhaps I could achieve what I want by

1) defining ANA_PASSIVE to initialize my inert tracer concentration to 0, and

2) defining point sources of tracer using either TRC_PSOURCE or ANA_TRC_PSOURCE. Here I assume that the tracer sources would be specified in a manner equivalent to river sources (e.g., position, Q(t), C(t)), so I could just define a "pulse" of tracer at time t = 1 day.

However, I can't seem to find an ana_trc_psource.h file or any other routine related to trc_psource in my ROMS distribution, although there are plenty of "if define TRC_PSOURCE" and "if define ANA_TRC_PSOURCE" scattered about. Are these options implemented in the code? Should the files be there? Any thought on whether this approach would work?

[kate]

Yes, I think that should work. There are a few ana_trc_psource.h files in https://github.com/kshedstrom/Apps, mostly (all?) for steady tracer sources. They are meant to go with my ROMS fork so I have no idea how incompatible the trunk code is these days.

[wilkin]

No, that won't work because the point sources will also introduce mass (you have to specify discharge in m3/s) and radically change the circulation.

I think your only choice is to stop and restart, or make some radical hack to reset dye_01 concentration at a chosen time (and I'm not sure how to do that).

[rduran]

I know this is not done often, but how about saving the velocity and eddy diffusivity often enough to do the Eulerian computations offline? It will require some coding of course but given the lack of better options ... is there a reason this will not work? (Especially if only interested in the sea-surface 2D problem.)

Another benefit of the offline approach is that Lagrangian Coherent Structures can be likewise computed to visualize the transport barriers that will likely explain the spatial structure of the e-flushing time you'd like to compute. Code for LCS is freely available and computes quickly (probably minutes).

Rodrigo.

[wilkin]

You might be able to do this with a hack to step3d_t.F using the tracer climatology nudging. See the code chuck near line 1559:

Code: Select all

!
!  Nudge towards tracer climatology.
!
        IF (LtracerCLM(itrc,ng).and.LnudgeTCLM(itrc,ng)) THEN
          DO k=1,N(ng)
            DO j=JstrR,JendR
              DO i=IstrR,IendR
                t(i,j,k,nnew,itrc)=t(i,j,k,nnew,itrc)+                  &
     &                             dt(ng)*                              &
     &                             CLIMA(ng)%Tnudgcof(i,j,k,ic)*        &
     &                             (CLIMA(ng)%tclm(i,j,k,ic)-           &
     &                              t(i,j,k,nnew,itrc))
              END DO
            END DO
          END DO
        END IF

You would want to introduce here some factor that switched on the nudging a short interval, say Tf, shortly before the nominal start time of your flushing experiment at, say, time = T0, and then switch off again at exactly time = T0. Then the tracer would quickly build up in place during T0-Tf < t < T0 then evolve only due to circulation. You'll get a little bit of bleeding out of the designated area, but the calculations here are all linear so you can keep the interval very short and the effective nudging rate high.

Tnudgcof(i,j,k,ic) would be set to some finite value (say, 1.0?) in the region you want to tag with passive tracer. This is set either by the file specified by NUDNAME in roms.in (see the template Data/ROMS/CDL/nudgcoef.cdl and variable tracer_NudgeCoef) ... or in Functional/ana_nudgcoef.h where you want to appropriately set CLIMA(ng)%Tnudgcof(i,j,k,itrc) for itrc = 3.

Or, set the spatial region via tclm. It doesn't matter becuase it is the product of Tnudgcoef and tclm that needs to be the region subset.

The switch on and off with time can be accomplished with a clever use of Fortran SIGN function.
tswitch1 = 0.5*(1.0+SIGN(1.0,time-T0-Tf)) ! goes from 0 to 1 at time = T0-Tf
tswitch2 = 0.5*(-1.0+SIGN(1.0,time-T0)) ! goes from 1 to 0 at time = T0

Then tswitch1*tswitch2 should be 1 only in the short interval Tf. (Check my math). You can roll it all into a factor that you add before, say, the dt(ng)* factor in the assignment of t(i,j,k,nnew,itrc).

Remember you'll need to activate ...

Code: Select all

  LnudgeTCLM == F F T                       ! temperature, salinity, inert

in roms.in.

[marcmestres]

I know this is not done often, but how about saving the velocity and eddy diffusivity often enough to do the Eulerian computations offline? It will require some coding of course but given the lack of better options ... is there a reason this will not work? (Especially if only interested in the sea-surface 2D problem.)

Another benefit of the offline approach is that Lagrangian Coherent Structures can be likewise computed to visualize the transport barriers that will likely explain the spatial structure of the e-flushing time you'd like to compute. Code for LCS is freely available and computes quickly (probably minutes).

Rodrigo.

Rodrigo,

that is indeed another option I had not thought about. I will look into it as an alternative to the ROMS "online" approach if I can't get this one to work. Thanks.

Marc

[marcmestres]

Thank you very much for your help, John. I will have a go at modifying step3d_t.F following your indications as soon as I can.

[marcmestres]

Seems I was a bit too optimistic with the "normal" use of the tracer part, i.e., releasing the tracer at t = 0 day. I have been doing (sucessful) preliminary tests with ROMS 3.7, releasing a passive tracer with C = 100 (units) in a specific area of my domain, and C = 0 elsewhere. The tracer then evolves in time in a logical way, and the concentration within my domain ranges between 0 and 100.
Now, when I switch to ROMS 3.8 (as part of the COAWST release, but compiled on its own, without any other model), and I use the same configuration I've been using for the 3.7 tests, the minimum concentration in the domain becomes negative, and the maximum concentration exceeds 100. For these runs I have defined T_PASSIVE, ANA_SPFLUX and ANA_BPFLUX, and provided dye boundary conditions (C=0) via the usual boundary file. At first I thought maybe the problem was related to sharp C gradients, since it occurs close to where the initial C gradients are, but if I replace C=0 with C=90 I get the same behavior. Any idea on why this could be happening?

M

[wilkin]

Look carefully at the logfile output from both runs. ROMS is very complete in documenting what it did, and you may stop some differences in configuration. Consider especially the settings for advection scheme because this can be a source of tracer over/under-shoot.

[jivica]

There are a few options you could explore;
(a) offline BIO ROMS (think Kate managed to build the code a long time ago) with your (pre)saved history file and then play with passive tracer(s), or
(b) use some other advection diffusion reaction code with ROMS currents and diffusion (doesn't matter is it finite difference or unstructured, you can always regird).
As John wrote, the last problem smells as low (odd order) advection scheme, just to be sure not to have over/under shootings you could use some positive definite advection such as MPDATA for example.

Cheers,
Ivica

[kate]

I gave up on offline bio because it fought with wetting and drying and all my domains have wetting and drying. Sorry, I probably left it in some half-baked stage.

[marcmestres]

MPDATA did the trick for the tracer concentration over/undershooting. Thanks, John and Ivica.

[marcmestres]

Hi again,

I finally got round to doing this, following John's instructions in a previous post.

a) I added some code in step3d_t.F that switches on the passive tracer nudging N timesteps before the time T I want it released, trying to keep N as small as possible.

b) I use a nudging.nc file in which I define a Tracer_NudgCoef = 0 in all the cells except the ones I will seed with tracer, where Tracer_NudgCoef = 1.

c) I define a clm.nc file with dye_01 = 0. At the cells I seed with tracer and at the dye_01_time at which I release it, dye_01 = 100

d) I add LnudgeTCLM = F F T in the in file

Doing this, I succeed in releasing the tracer at time T during the simulation, which is what I wanted, and it evolves as expected in time. However, the C values I get are extremely small. For instance, by setting dye_01=100, the initial concentration I get at time T is 0.00077 at all the seeded cells. I have tried changing N (i.e., nudging over a longer/shorter time period) and even different values of the tracer nudging coefficient, and the initial C changes, but I don't seem to get anywhere near the concentration values I set in the clm file. Any idea on why this could be happening?

Thanks

M

[wilkin]

Tracer_NudgCoef is an inverse time scale. So it will interact with DT and your concentration of 100 to approach 100 in discrete steps. The ideal solution will approach 100 as a negative exponential.

You can probably quickly calculate numerically what those discrete steps would be under the assumption of no other physics (diffusion, or advection taking the tracer away). To get a more rapid response, shorten the time scale or increase the concentration.

[marcmestres]

If I set the nudging coefficient to 86400/DT, I get an "initial" concentration at time T that is close enough (approxs. 99.5%) to the value I've specified in the clm file. This is fine for my purposes.
Again, thank you very much for your help!