# Contour plot of streamfunction

I am trying to obtain a graph similar to this one:

D1 = 60;
ρ0 = 1025;
f = 0.545/10^4;
c = ((1 + I)*Pi)/D1;
τx = 0;
τy = 0.07;
τ = τx + I*τy;
ug = 0;
u0 = ((1 - I)*Pi*τ)/(D1*ρ0*f);
Clear[h, ζ, α, vg, ugc, Uek, ψ]
h[x_] := -x/625 + 4
ζ[x_] := (Pi*h[x])/D1
α[x_] := (Cosh[ζ[x]]*Cos[ζ[x]])^2 + (Sinh[ζ[x]]*Sin[ζ[x]])^2
S1[x_] := (Cosh[ζ[x]]*Cos[ζ[x]])/α[x]
S2[x_] := (Sinh[ζ[x]]*Sin[ζ[x]])/α[x]
T1[x_] := (Cosh[ζ[x]]*Sinh[ζ[x]])/α[x]
T2[x_] := (Cos[ζ[x]]*Sin[ζ[x]])/α[x]
vg[x_] := (2*Pi*(((1 - S1[x])*τy)/(T1[x] - T2[x]) + (S2[x]*τx)/(T1[x] - T2[x])))/(ρ0*f*D1) - (ug*(T1[x] + T2[x] - (2*Pi*h[x])/D1))/(T1[x] - T2[x])
ugc[x_] := ug + I*vg[x]
Uek[x_] := τy/(ρ0*f)
ψ[x_, z_] := Re[(1/c)*(u0*(1 - Cosh[c*(z + h[x])]/Cosh[c*h[x]]) + ugc[x]*(Sinh[c*z]/Cosh[c*h[x]]))] - ug[x]*z
With[{x = -40000, z = -300}, {h[x], h[x]/D1, Uek[x], ψ[x, z]}]
ContourPlot[ψ[x, z]/Uek[x]/100, {x, -100000, 0}, {z, -100, 0}]


The output is the following:

The stream function is presented as a series of operations instead of a single value.

What is happening?

Some definitions:

• Is there a minimum working example that works? ... perhaps with the smallest complement of variables? Commented Sep 29, 2019 at 15:18
• I can't get past the error in the φ calculation... Commented Sep 29, 2019 at 15:37

There are a number of typos in your code. As well as correcting them, I recommend writing the code so it doesn't pass functions names as arguments; i.e., writing them so the x are z are only arguments. Also, I would use SetDelayed (:=) in place of Set (=) in your situation. With the typos correcting an the changes made you code looks like this:

Lx = 200;
D1 = 60;
ρ0 = 1025;
f = 0.545/10^4;
hmin = 4;
c = ((1 + I) π)/D1;
τx = 0;
τy = 0.07;
τ = τx + I τy;
ug = 0;
u0 = ((1 - I) π τ)/(D1 ρ0 f);

Clear[h, ζ, α, vg, ugc, Uek, ψ]

h[x_] := (Lx - x)/10^3 + hmin;
ζ[x_] := (π h[x])/D1;
α[x_] := (Cosh[ζ[x]]*Cos[ζ[x]])^2 + (Sinh[ζ[x]]*Sin[ζ[x]])^2

S1[x_] := (Cosh[ζ[x]]*Cos[ζ[x]])/α[x]
S2[x_] := (Sinh[ζ[x]]*Sin[ζ[x]])/α[x]
T1[x_] := (Cosh[ζ[x]]*Sinh[ζ[x]])/α[x]
T2[x_] := (Cos[ζ[x]]*Sin[ζ[x]])/α[x]

vg[x_] :=
(2 π (((1 - S1[x])τy)/(T1[x] - T2[x]) + (S2[x]τx)/(T1[x] - T2[x])))/(ρ0 f D1) -
(ug (T1[x] + T2[x] - (2 π h[x])/D1))/(T1[x] - T2[x])
ugc[x_] := ug + I vg[x]
Uek[x_] :=
(1 - S1[x]) τy/(ρ0 f) + S2[x] τx/(ρ0 f) -
D1/(2 π) ((T1[x] - T2[x]) vg[x] + (T1[x] + T2[x]) ug)


Edited to correct an error pointed out by the OP in a comment below

ψ[x_, z_] :=
Re[1/c (u0 (1 - Cosh[c*(z + h[x])]/Cosh[c h[x]]) + ugc[x] Sinh[c z]/Cosh[c h[x]])] -
ug z


Then you can evaluate for value ranges of x and z like this:

data =
Table[
Quiet @ Check[ψ[x, z]/Uek[x], Indeterminate],
{x, -100, 0, 20}, {z, -100, 0, 20}];


which gives

These values don't look good to me. They certainly won't give a stream plot like the one you show in your question.

Note that

With[{x = -40000, z = -300}, {h[x], h[x]/D1, Uek[x], ψ[x, z]}]


gives me

{221/5, 221/300, 0., 35559.5}

which is different from you show in your edited question. Are you sure your definitions are correct? I am suspicious of the computation of Uek[x].

• Thank you for your detailed reply. I corrected some values and changed some expressions, but also a typo that you made in the stream function expression: the last term of the expression had a ug not ugc multiplying the z variable. Commented Sep 30, 2019 at 17:14
• @ASPVL. Thanks for pointing out my error. I have corrected it. However, there is still a serious problem with code which I believe can only found by someone who understands the science behind the math. I am not such a person, so I don't think I help you further. Commented Sep 30, 2019 at 18:15
• Thank you for your patience, I have added some definitions to the body of the question. I think the formula for Uek (Ekman transport) is written correctly. Though I have to agree with you that the values seem quite wrong... Commented Oct 3, 2019 at 16:18