Skip to content
Snippets Groups Projects
Commit 8a8a7d7e authored by patavirt's avatar patavirt
Browse files

examples: cpr_sns: adjust n-convergence plotting

parent 928fafd1
No related branches found
No related tags found
No related merge requests found
...@@ -112,10 +112,11 @@ def do(W_xi=6, multin=False): ...@@ -112,10 +112,11 @@ def do(W_xi=6, multin=False):
) )
if multin: if multin:
mphi = phi[::2]
res0 = j( res0 = j(
T, T,
h, h,
phi[None, :], mphi[None, :],
alpha_soc=alpha, alpha_soc=alpha,
eta=eta, eta=eta,
L=L[:, None], L=L[:, None],
...@@ -126,7 +127,7 @@ def do(W_xi=6, multin=False): ...@@ -126,7 +127,7 @@ def do(W_xi=6, multin=False):
res1 = j( res1 = j(
T, T,
h, h,
phi[None, :], mphi[None, :],
alpha_soc=alpha, alpha_soc=alpha,
eta=eta, eta=eta,
L=L[:, None], L=L[:, None],
...@@ -142,19 +143,19 @@ def do(W_xi=6, multin=False): ...@@ -142,19 +143,19 @@ def do(W_xi=6, multin=False):
if multin: if multin:
Jx0_mean = np.asarray( Jx0_mean = np.asarray(
[Res(*x).Jx[1:-1].sum(axis=1).mean(axis=0) for x in res0.flat] [Res(*x).Jx[1:-1].sum(axis=1).mean(axis=0) for x in res0.flat]
).reshape(res.shape) ).reshape(res0.shape)
Jx1_mean = np.asarray( Jx1_mean = np.asarray(
[Res(*x).Jx[1:-1].sum(axis=1).mean(axis=0) for x in res1.flat] [Res(*x).Jx[1:-1].sum(axis=1).mean(axis=0) for x in res1.flat]
).reshape(res.shape) ).reshape(res1.shape)
fig, axs = plt.subplots(1, 2, layout="compressed") fig, axs = plt.subplots(1, 2, layout="compressed")
ax = axs[0] ax = axs[0]
ax.plot(phi / pi, Jx_mean.T / abs(Jx_mean).max(axis=1)) ax.plot(phi / pi, Jx_mean.T / abs(Jx_mean).max(axis=1))
if multin: if multin:
ax.plot(phi / pi, Jx0_mean.T / abs(Jx_mean).max(axis=1), "k:") ax.plot(mphi / pi, Jx0_mean.T / abs(Jx0_mean).max(axis=1), "k:", alpha=0.25)
ax.plot(phi / pi, Jx1_mean.T / abs(Jx_mean).max(axis=1), "k:") ax.plot(mphi / pi, Jx1_mean.T / abs(Jx1_mean).max(axis=1), "k:")
ax.set_xlabel(r"$\varphi / \pi$") ax.set_xlabel(r"$\varphi / \pi$")
ax.set_ylabel(r"$I / I_{\mathrm{max}}$") ax.set_ylabel(r"$I / I_{\mathrm{max}}$")
ax.legend(L / xi, title=r"$L/\xi$", loc="lower right") ax.legend(L / xi, title=r"$L/\xi$", loc="lower right")
...@@ -167,7 +168,7 @@ def do(W_xi=6, multin=False): ...@@ -167,7 +168,7 @@ def do(W_xi=6, multin=False):
ax = axs[1] ax = axs[1]
ax.plot(L / xi, 100 * eff(Jx_mean)) ax.plot(L / xi, 100 * eff(Jx_mean))
if multin: if multin:
ax.plot(L / xi, 100 * eff(Jx0_mean), "k:") ax.plot(L / xi, 100 * eff(Jx0_mean), "k:", alpha=0.25)
ax.plot(L / xi, 100 * eff(Jx1_mean), "k:") ax.plot(L / xi, 100 * eff(Jx1_mean), "k:")
ax.set_xlabel(r"$L / \xi$") ax.set_xlabel(r"$L / \xi$")
ax.set_ylabel(r"$\eta$ [%]") ax.set_ylabel(r"$\eta$ [%]")
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment