Function for tracers

[ammtD]

Hi,

I would like to force a tracer to a specific value at a specific time and location in my grid. Any creative suggestions on how to do this?
To be more clear, tracer :Rh should be zero everywhere until t = ti, and will be set to 1 at t=ti in part of the grid.

Rh(x, y, z, t) = 0   for t<ti   
and then force  
Rh (x[1:3], y[:], z[1:3] , ti) = 1 

Any help appreciated!

[glwagner]

I think you'll have to chain | logic in a discrete_form forcing function. I tried to make it readable:

@inline t_unforced(t) = t < ti
@inline i_unforced(i) = i < 1 | i > 3
@inline k_unforced(k) = k < 1 | k > 3
@inline unforced(i, k, t) = t_unforced(t) | i_unforced(i) | k_unforced(k)
@inline forcing(i, j, k, grid, clock, fields) = ifelse(unforced(i, k, clock.time), 0.0, 1.0)

c_forcing = Forcing(forcing, discrete_form=true)

I think that will work...

[ammtD]

Finally getting around to follow up on this component of my model construction...

I've tried the above. In principle it seems to do the job, though once the forcing turns on the tracer field does look funny, a sequence of high and low values (see link to animation here)..
The main dispersion of the tracer from the surface does look OK.

Another problem I encounter is when I add a second condition on the first line (aiming to turn of the forcing at a certain time as well):

ti=15minutes;  tf = 20minutes

@inline t_unforced(t) = t < ti | t > tf
@inline i_unforced(i) = i < 1 | i > 3
@inline k_unforced(k) = k < 1 | k > 3
@inline unforced(i, k, t) = t_unforced(t) | i_unforced(i) | k_unforced(k)
@inline forcing(i, j, k, grid, clock, fields) = ifelse(unforced(i, k, clock.time), 0.0, 1.0)
c_forcing = Forcing(forcing, discrete_form=true)

Doesn't work, this model crashes with the error:

TaskFailedException

   nested task error: MethodError: no method matching |(::Float64, ::Float64)
   Closest candidates are:
     |(::Any, ::Any, ::Any, ::Any...) at operators.jl:560
   Stacktrace:
     [1] call
       @ C:\Users\Owner\.julia\packages\Cassette\1lyEM\src\context.jl:456 [inlined]
     [2] fallback
       @ C:\Users\Owner\.julia\packages\Cassette\1lyEM\src\context.jl:454 [inlined]
     [3] _overdub_fallback(::Any, ::Vararg{Any, N} where N)
       @ C:\Users\Owner\.julia\packages\Cassette\1lyEM\src\overdub.jl:586 [inlined]
     [4] overdub(::Cassette.Context{nametype(CPUCtx), Nothing, Nothing, KernelAbstractions.var"##PassType#257", Nothing, Cassette.DisableHooks}, ::typeof(|), ::Float64, ::Float64)
       @ Cassette C:\Users\Owner\.julia\packages\Cassette\1lyEM\src\overdub.jl:586
     [5] t_unforced(::Float64)
       @ .\In[5]:4 [inlined]
     [6] overdub
       @ .\In[5]:4 [inlined]
     [7] unforced(::Int64, ::Int64, ::Float64)
       @ .\In[5]:7 [inlined]
     [8] overdub
       @ .\In[5]:7 [inlined]
     [9] forcing(::Int64, ::Int64, ::Int64, ::RegularRectilinearGrid{Float64, Periodic, 
```  etc 


Any idea why? 

[glwagner]

It's because of order-of-operations rules; | takes precedence over >. So you need to write

@inline t_unforced(t) = (t < ti) | (t > tf)
@inline i_unforced(i) = (i < 1) | (i > 3)
@inline k_unforced(k) = (k < 1) | (k > 3)

It's not too surprising that you see oscillatory errors / noise with a sharp forcing. You could consider smoothing the forcing out with ramps / tanh rather than using sharp cutoffs.

[ammtD]

Oops sorry, OK. fixed that! Works nicely now, thank a lot!
Testing it I was using a 45s timestep; reducing that seems to have solved the noise error.