Curl Operation

[anthony-meza]

MeshArrays.jl/src/Operations.jl

Lines 84 to 116 in dd3f19a

	"""
	curl(u::MeshArray,v::MeshArray,Γ::NamedTuple)
	Compute curl of a velocity field.
	"""
	function curl(u::MeshArray,v::MeshArray,Γ::NamedTuple)
	uvcurl=similar(Γ.XC)
	fac=exchange(1.0 ./Γ.RAZ,1)
	(U,V)=exchange(u,v,1)
	(DXC,DYC)=exchange(Γ.DXC,Γ.DYC,1)
	[DXC.MA[i].=abs.(DXC.MA[i]) for i in eachindex(U.MA)]
	[DYC.MA[i].=abs.(DYC.MA[i]) for i in eachindex(V.MA)]
	for i in eachindex(U.MA)
	ucur=U.MA[i][2:end,:]
	vcur=V.MA[i][:,2:end]
	tmpcurl=ucur[:,1:end-1]-ucur[:,2:end]
	tmpcurl=tmpcurl-(vcur[1:end-1,:]-vcur[2:end,:])
	tmpcurl=tmpcurl.*fac.MA[i][1:end-1,1:end-1]
	##still needed:
	##- deal with corners
	##- if putCurlOnTpoints
	tmpcurl=1/4*(tmpcurl[1:end-1,2:end]+tmpcurl[1:end-1,1:end-1]+
	tmpcurl[2:end,2:end]+tmpcurl[2:end,1:end-1])
	uvcurl[i]=tmpcurl
	end
	return uvcurl
	end

Should U and V fields be weighted by some distance before computing the curl?

E.g.,

[DXC[i].=abs.(DXC[i]) for i in eachindex(U)]
[DYC[i].=abs.(DYC[i]) for i in eachindex(V)]

#not currently in code
U_DY = U .* DYC #weight by length of opposing side
V_DX = V .* DXC

###compute curl afterward###

[gaelforget]

I understand this question was asked in context of an Ekman pumping calculation.

Can you please try with this?

The underlying function, EkmanTrsp, takes in wind stress and outputs Ekman transport per matter (in m2/s).

To get to Ekman pumping from there you'd calculate the divergence.

Does the results match your expectation for Ekman pumping magnitude?

[anthony-meza]

An Example Climatology of Ekman Pumping Velocity (meters per year) from a lecture I found online

[anthony-meza]

Okay, here are some plots I made of Ekman Transport using different methods.

The "current Implementation" plot shows an Ekman pumping velocity (in meters per year) calculated by taking the curl of the wind stress using the current implementation of curl in MeshArrays.

The "convergence" plot refers to the Ekman pumping velocity calculated using the convergence of the Ekman Transport fields generated by the example you pointed me to. I am using the gradients to calculate this.

The "area-weighted curl" plot refers to the same calculation using the correction to curl that I suggest in this issue.

[gaelforget]

Many thanks @anthony-meza !

I am starting to think that there is a bug in the curl function. Shouldn't be sqrt(RAZ)?

MeshArrays.jl/src/Operations.jl

Line 92 in 45be7ce

fac=exchange(1.0 ./Γ.RAZ,1)

tagging @jm-c to the rescue 👋

[anthony-meza]

I'm not sure actually!

To me, this function is application of Stokes' Theorem. Is this correct? If it is, wouldn't the line integral need to be weighted by a distance? I may be misunderstanding what the code is doing though

Here's an attempt with the sqrt(area) correction. This seems to do pretty well in one tile, but not so much in the others! The smaller scale features are a bit different from what the other methods predict.

!

I think this could just require an adjustment in documentation. We could just say the curl function takes a unit in $m/s$ and the curl operation would return something in $1/(ms)$!