Diffusion - ROMS fields

[ana-machado]

Hello,
I am trying to add diffusion to my simulations and I am a bit lost on the choices. I am using Parcels with ROMS outputs, and I have a model output that is 3D Diffusion, that would be much better to use as diffusion, instead of constant values. I have already included this variable in my FieldSet (it is a 4D variable). Is there any way I can use it, taking advantage of the diffusion kernels already available?
Thank you very much.

[erikvansebille]

Dear @ana-machado. Thanks for your question. IN my group, we don't use ROMS so we don't have much experience with this. However, perhaps you can find some useful background/ideas in @daanreijnders' preprint on isoneutral subgrid-scale dispersion (implemented in MITgcm)?

[daanreijnders]

Hi Ana, when you say 3D diffusion, does this mean that you have different values of horizontal diffusivity (possibly different in the x- and y- directions) AND vertical diffusivity, or only one of these, but stored in a lon-lat-depth array? Either way, it's definitely possible to write up a kernel so that you could use these in Parcels, and I'll gladly help you on your way.

One caveat is that if the diffusivities are spatially inhomogeneous, you need to know their derivatives in order to implement random walk diffusion correctly (see eq. 28 in van Sebille et al.). Currently, Parcels cannot compute derivatives internally yet (see #970 and #979), so you would either need to pre-process the data by computing spatial derivatives of the diffusivities, or use an in-kernel approximation, as is currently implemented in the advection-diffusion kernels (see also this tutorial).

[ana-machado]

Hi Daan!

Thank you very much for your reply! I don't know why I dind't get the notice of it and just saw it now. Hope you are still available to help me :) What I have are 3D fields, but of horizontal diffusion. Have already added it in my FieldSet, but I am not sure how I should proceed. I was looking at the "acvection-diffusion kernels" (https://github.com/OceanParcels/parcels/blob/202ef087ba8608d4950a096569a6992cf278d022/parcels/application_kernels/advectiondiffusion.py#L12). Which one do you think I should try to implement? In my case I don't have zonal and meridional Kh, but I will probably need to consider the diff3d values for both. Do you think this could me an option? My cells do not have the same length in lat and lon directions, but I can obtain the 2D field of dx and dy values and add them to Fieldset. It this the way?

Some usefull information:
Description of diff3d variable:
float diff3d(time, s_rho, eta_rho, xi_rho) ; long_name = "averaged horizontal diffusivity coefficient" ; standard_name = "ocean_tracer_xy_laplacian_diffusivity" ;

I am using AdvectionRK4, plus some kernels to simulate biological behaviours (my particles are sardine larvae :))
pset.execute(AdvectionRK4 + pset.Kernel(Age) + pset.Kernel(devStage_and_rhoVartime) + pset.Kernel(temperature_mortality) , runtime=delta(days=30), dt=delta(minutes=5), output_file=pset.ParticleFile(name=outFileName, outputdt=delta(minutes=5)),recovery={ErrorCode.ErrorOutOfBounds: DeleteParticle})

Thank you very much, I hope you can help me.

Ana

[daanreijnders]

Hi Ana, sorry for my late reply.

If I understand you correctly, it's a 3D scalar field, right? In that case, I think you may be able to use AdvectionDiffusionM1. You mention you don't have separate meridional and zonal diffusivities, so you can just use the same scalar field for both. The only thing I am not entirely sure of is whether there should be a correction term for the vertical gradient in the horizontal diffusivity fields. By how much does your diffusivity typically change over, say, 1 or 10 meters vertically?

You could also load in the field's horizontal resolutions in the zonal and meridional directions. If you have a look at AdvectionDiffusionM1, you can change dres to load the local field resolution. You'll need a separate resolution in each direction, so I suggest you introduce dres_zonal and dres_meridional.

By the way, if you use the AdvectionDiffusionM1 or AdvectionDiffusionEM kernels, you should not use AdvectionRK4 besides it, since the advection is already taken care of.

Hope this helps. If anything is unclear, just let me know.