step3d_t.F bug?

[johnluick]

I find that performance is improved tremendously if the following change is made in step3d_t.F (in the part where mass point sources are applied)...

Code: Select all

FROM:
                DO k=1,N(ng)                
                  FC(i,k)=FC(i,k)+0.5_r8*                              &
     &                    (Qsrc(is,k  )*Tsrc(is,k  ,itrc)+             &
     &                     Qsrc(is,k+1)*Tsrc(is,k+1,itrc))
                END DO
TO:
                DO k=1,N(ng)-1
                  FC(i,k-1)=FC(i,k-1)+0.5_r8*                           &
     &                    (Qsrc(is,k  )*Tsrc(is,k  ,itrc)+              &
     &                     Qsrc(is,k+1)*Tsrc(is,k+1,itrc))
                END DO
                k=n(ng)
                FC(i,k)=FC(i,k)+Qsrc(is,k)*Tsrc(is,k  ,itrc)

Noting that FC is dimensioned as FC(Nsrc,0:N), whereas Qsrc and Tsrc are dimensioned as Qsrc(Nsrc,1:N) and Tsrc(Nsrc,1:N) you see that the limits on the original loop don't look right. Furthermore FC(is,O) is not set but required immediately after this loop.

Unfortunately, this does not completely fix step3d_t.F because if you integrate the concentration (in my case, of a passive tracer dye_01) in ocean_his.nc over the entire volume at each time step it is much less than what is being added in the forcing netcdf file. This note is contributed with the intent of getting someone to fix this properly as I have no further time for it for a while. Thanks! By the way I am testing it with TS_MPDATA and Q_PSOURCE.

[wilkin]

John,

I don't follow what you're trying to do here.

The FC are fluxes are at the "w-points" interfaces between "rho-points" centered divergences. So FC is dimensioned 0:N (bottom to surface) whereas the cell centered divergences are dimensioned 1:N.

ROMS uses a predictor/corrector time step, and the surface/bottom boundary conditions on vertical diffusive fluxes are applied on the predictor step only. So in step3d_t the fluxes FC(0) and FC(N) must be zero.

Your proposed changes would appear to violate FC(0)=FC(N)=0 in step3d_t.

John.

[johnluick]

John,

sorry could not get back to this sooner. Thanks very much for explaining that FC is defined at the w-points. In practical terms, since the fluxes are required at rho points, and this is done by taking a vertical average of the flux above and below each rho-point, the consequence of defining a nonzero FC at the bottom w-point (seabed) is simply to increase the effective horizontal flux of dye at the bottom. Maybe that is why the little hack I did actually made the resulting dye distribution LOOK so much more sensible.

My main concern was that the total amount of passive tracer in the domain (up until the time it reaches the boundary) does not correspond to the amount specified in the input netcdf file. This applied in ROMS as shipped and was not improved by my hack of step3d_t.F. Since TS_MPDATA seems to work very well for points distant from an interior source, I suspected that it had something to do with the way the source is input with Q_PSOURCE. I still think that. But I am hoping someone else does this check on the quantities (or has already done it and can tell me I'm doing something wrong).

Someone in the forum said (or seemed to say) that they abandoned UV_PSOURCE (for river input at the boundary) in favor of a "near-field" model of some sort. I can see how that would be done at the boundary and it makes a lot of sense but I can't see how to implement it for a source in the interior. If someone has any ideas on that I'd like to hear them.

John

[wilkin]

John,

The fluxes are not required at the rho-points. The flux divergence is required at the rho-points, so the FC are not averaged, but differenced. To conserve tracer its essential that FC = 0 at k=0 (z=-h). Otherwise there is an implied flux of tracer in to or out of the seabed.

I have not checked for the exact conservation of the Q SOURCE option, but I do get exact conservation of passive tracers (to 1e-11) using the tracer diagnostics terms and the default surface/bottom passive tracer fluxes of zero.

[johnluick]

John,

thanks for your comments which are always very helpful. I guess I should not have used the term "fluxes" as it was misleading - what I was referring to was the averaging (looks like averaging to me) done in the equation in my first post on this thread (averaging of Qsrc*Tsrc above and below the w point).

John