Page MenuHomec4science

plot_inf_set_3.py
No OneTemporary

File Metadata

Created
Sun, Apr 28, 15:24

plot_inf_set_3.py

import warnings
from copy import deepcopy as copy
import numpy as np
from matplotlib import pyplot as plt
def plotInfSet31FromData(mNone, mus, Z, T, R, E, beta, murange, approx, mu0,
exp = 2., normalizeDen = False):
if hasattr(approx, "mus"):
mu2x = approx.mus(mNone[0]) ** exp
else:
mu2x = mu0[0] ** exp
murangeExp = [[murange[0][mNone[0]] ** exp],
[murange[1][mNone[0]] ** exp]]
mu1 = np.real(np.power([m[mNone[0]] for m in mus], exp))
ZTmin, ZTmax = min(np.min(Z), np.min(T)), max(np.max(Z), np.max(T))
Rmin, Rmax = np.min(R), np.max(R)
Emin, Emax = np.min(E), np.max(E)
if not np.isnan(beta[0]):
eta = R / beta / E
betamin, betamax = np.min(beta), np.max(beta)
etamin, etamax = np.min(eta), np.max(eta)
plt.figure(figsize = (15, 7))
plt.jet()
plt.semilogy(mu1, Z)
plt.semilogy(mu1, T, '--')
mTMP = mus[0]
mTMP[mNone[0]] = None
for l_ in approx.trainedModel.getPoles(mTMP):
plt.plot([np.real(l_ ** exp)] * 2, [ZTmin, ZTmax], 'b:')
plt.plot(mu2x, [ZTmin] * len(mu2x), 'kx')
plt.plot([murangeExp[0][0]] * 2, [ZTmin, ZTmax], 'm:')
plt.plot([murangeExp[1][0]] * 2, [ZTmin, ZTmax], 'm:')
plt.xlim(mu1[0], mu1[-1])
plt.title("|u(mu)|, |u_app(mu)|")
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
plt.semilogy(mu1, R)
for l_ in approx.trainedModel.getPoles(mTMP):
plt.plot([np.real(l_ ** exp)] * 2, [Rmin, Rmax], 'b:')
plt.plot(mu2x, [Rmax] * len(mu2x), 'kx')
plt.plot([murangeExp[0][0]] * 2, [Rmin, Rmax], 'm:')
plt.plot([murangeExp[1][0]] * 2, [Rmin, Rmax], 'm:')
plt.xlim(mu1[0], mu1[-1])
if normalizeDen:
plt.title("|Q(mu)res(mu)|")
else:
plt.title("|res(mu)|")
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
plt.semilogy(mu1, E)
for l_ in approx.trainedModel.getPoles(mTMP):
plt.plot([np.real(l_ ** exp)] * 2, [Emin, Emax], 'b:')
plt.plot(mu2x, [Emax] * len(mu2x), 'kx')
plt.plot([murangeExp[0][0]] * 2, [Emin, Emax], 'm:')
plt.plot([murangeExp[1][0]] * 2, [Emin, Emax], 'm:')
plt.xlim(mu1[0], mu1[-1])
if normalizeDen:
plt.title("|Q(mu)err(mu)|")
else:
plt.title("|err(mu)|")
plt.grid()
plt.show()
if not np.isnan(beta[0]):
plt.figure(figsize = (15, 7))
plt.jet()
plt.plot(mu1, beta)
for l_ in approx.trainedModel.getPoles(mTMP):
plt.plot([np.real(l_ ** exp)] * 2, [betamin, betamax], 'b:')
plt.plot(mu2x, [betamax] * len(mu2x), 'kx')
plt.plot([murangeExp[0][0]] * 2, [betamin, betamax], 'm:')
plt.plot([murangeExp[1][0]] * 2, [betamin, betamax], 'm:')
plt.xlim(mu1[0], mu1[-1])
plt.title("beta(mu)")
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
plt.semilogy(mu1, R / beta)
plt.semilogy(mu1, E, '--')
for l_ in approx.trainedModel.getPoles(mTMP):
plt.plot([np.real(l_ ** exp)] * 2, [Emin, Emax], 'b:')
plt.plot(mu2x, [Emax] * len(mu2x), 'kx')
plt.plot([murangeExp[0][0]] * 2, [Emin, Emax], 'm:')
plt.plot([murangeExp[1][0]] * 2, [Emin, Emax], 'm:')
plt.xlim(mu1[0], mu1[-1])
if normalizeDen:
plt.title("|Q(mu)res(mu)/beta(mu)|")
else:
plt.title("|res(mu)/beta(mu)|")
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
plt.semilogy(mu1, eta)
for l_ in approx.trainedModel.getPoles(mTMP):
plt.plot([np.real(l_ ** exp)] * 2, [etamin, etamax], 'b:')
plt.plot(mu2x, [etamax] * len(mu2x), 'kx')
plt.plot([murangeExp[0][0]] * 2, [etamin, etamax], 'm:')
plt.plot([murangeExp[1][0]] * 2, [etamin, etamax], 'm:')
plt.xlim(mu1[0], mu1[-1])
plt.title("eta(mu)")
plt.grid()
plt.show()
def plotInfSet32FromData(mNone, mus, Ze, Te, Re, Ee, beta, murange, approx,
mu0, exps = [2., 1.], clip = -1,
normalizeDen = False):
if hasattr(approx, "mus"):
mu2x = approx.mus(mNone[0]) ** exps[0]
mu2y = approx.mus(mNone[1]) ** exps[1]
else:
mu2x, mu2y = mu0[mNone[0]] ** exps[0], mu0[mNone[1]] ** exps[1]
murangeExp = [[murange[0][mNone[0]] ** exps[0],
murange[0][mNone[1]] ** exps[1]],
[murange[1][mNone[0]] ** exps[0],
murange[1][mNone[1]] ** exps[1]]]
mu1s = np.unique([m[mNone[0]] for m in mus])
mu2s = np.unique([m[mNone[1]] for m in mus])
mu1 = np.power(mu1s, exps[0])
mu2 = np.power(mu2s, exps[1])
Mu1, Mu2 = np.meshgrid(np.real(mu1), np.real(mu2))
Z = np.log10(Ze)
T = np.log10(Te)
R = np.log10(Re)
E = np.log10(Ee)
ZTmin, ZTmax = min(np.min(Z), np.min(T)), max(np.max(Z), np.max(T))
Rmin, Rmax = np.min(R), np.max(R)
Emin, Emax = np.min(E), np.max(E)
if not np.isnan(beta[0, 0]):
betamin, betamax = np.min(beta), np.max(beta)
if clip > 0:
ZTmax -= clip * (ZTmax - ZTmin)
cmap = plt.cm.bone
else:
cmap = plt.cm.jet
warnings.simplefilter("ignore", category = (UserWarning,
np.ComplexWarning))
plt.figure(figsize = (15, 7))
plt.jet()
p = plt.contourf(Mu1, Mu2, Z, cmap = cmap,
levels = np.linspace(ZTmin, ZTmax, 50))
if clip > 0:
plt.contour(Mu1, Mu2, Z, [ZTmin])
plt.plot(mu2x, mu2y, 'kx')
plt.plot([murangeExp[0][0]] * 2,
[murangeExp[0][1], murangeExp[1][1]], 'm:')
plt.plot([murangeExp[0][0], murangeExp[1][0]],
[murangeExp[1][1]] * 2, 'm:')
plt.plot([murangeExp[1][0]] * 2,
[murangeExp[1][1], murangeExp[0][1]], 'm:')
plt.plot([murangeExp[1][0], murangeExp[0][0]],
[murangeExp[0][1]] * 2, 'm:')
plt.colorbar(p)
plt.title("log10|u(mu)|")
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
p = plt.contourf(Mu1, Mu2, T, cmap = cmap,
levels = np.linspace(ZTmin, ZTmax, 50))
if clip > 0:
plt.contour(Mu1, Mu2, T, [ZTmin])
plt.plot(mu2x, mu2y, 'kx')
plt.plot([murangeExp[0][0]] * 2,
[murangeExp[0][1], murangeExp[1][1]], 'm:')
plt.plot([murangeExp[0][0], murangeExp[1][0]],
[murangeExp[1][1]] * 2, 'm:')
plt.plot([murangeExp[1][0]] * 2,
[murangeExp[1][1], murangeExp[0][1]], 'm:')
plt.plot([murangeExp[1][0], murangeExp[0][0]],
[murangeExp[0][1]] * 2, 'm:')
plt.title("log10|u_app(mu)|")
plt.colorbar(p)
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
p = plt.contourf(Mu1, Mu2, R, levels = np.linspace(Rmin, Rmax, 50))
plt.plot(mu2x, mu2y, 'kx')
plt.plot([murangeExp[0][0]] * 2,
[murangeExp[0][1], murangeExp[1][1]], 'm:')
plt.plot([murangeExp[0][0], murangeExp[1][0]],
[murangeExp[1][1]] * 2, 'm:')
plt.plot([murangeExp[1][0]] * 2,
[murangeExp[1][1], murangeExp[0][1]], 'm:')
plt.plot([murangeExp[1][0], murangeExp[0][0]],
[murangeExp[0][1]] * 2, 'm:')
if normalizeDen:
plt.title("log10|Q(mu)res(mu)|")
else:
plt.title("log10|res(mu)|")
plt.colorbar(p)
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
p = plt.contourf(Mu1, Mu2, E, levels = np.linspace(Emin, Emax, 50))
plt.plot(mu2x, mu2y, 'kx')
plt.plot([murangeExp[0][0]] * 2,
[murangeExp[0][1], murangeExp[1][1]], 'm:')
plt.plot([murangeExp[0][0], murangeExp[1][0]],
[murangeExp[1][1]] * 2, 'm:')
plt.plot([murangeExp[1][0]] * 2,
[murangeExp[1][1], murangeExp[0][1]], 'm:')
plt.plot([murangeExp[1][0], murangeExp[0][0]],
[murangeExp[0][1]] * 2, 'm:')
if normalizeDen:
plt.title("log10|Q(mu)err(mu)|")
else:
plt.title("log10|err(mu)|")
plt.colorbar(p)
plt.grid()
plt.show()
if not np.isnan(beta[0, 0]):
plt.figure(figsize = (15, 7))
plt.jet()
p = plt.contourf(Mu1, Mu2, beta,
levels = np.linspace(betamin, betamax, 50))
plt.plot(mu2x, mu2y, 'kx')
plt.plot([murangeExp[0][0]] * 2,
[murangeExp[0][1], murangeExp[1][1]], 'm:')
plt.plot([murangeExp[0][0], murangeExp[1][0]],
[murangeExp[1][1]] * 2, 'm:')
plt.plot([murangeExp[1][0]] * 2,
[murangeExp[1][1], murangeExp[0][1]], 'm:')
plt.plot([murangeExp[1][0], murangeExp[0][0]],
[murangeExp[0][1]] * 2, 'm:')
plt.title("beta(mu)")
plt.colorbar(p)
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
p = plt.contourf(Mu1, Mu2, R - np.log10(beta),
levels = np.linspace(Emin, Emax, 50))
plt.plot(mu2x, mu2y, 'kx')
plt.plot([murangeExp[0][0]] * 2,
[murangeExp[0][1], murangeExp[1][1]], 'm:')
plt.plot([murangeExp[0][0], murangeExp[1][0]],
[murangeExp[1][1]] * 2, 'm:')
plt.plot([murangeExp[1][0]] * 2,
[murangeExp[1][1], murangeExp[0][1]], 'm:')
plt.plot([murangeExp[1][0], murangeExp[0][0]],
[murangeExp[0][1]] * 2, 'm:')
if normalizeDen:
plt.title("log10|Q(mu)res(mu)/beta(mu)|")
else:
plt.title("log10|res(mu)/beta(mu)|")
plt.colorbar(p)
plt.grid()
plt.show()
plt.figure(figsize = (15, 7))
plt.jet()
p = plt.contourf(Mu1, Mu2, R - np.log10(beta) - E, 50)
plt.plot(mu2x, mu2y, 'kx')
plt.plot([murangeExp[0][0]] * 2,
[murangeExp[0][1], murangeExp[1][1]], 'm:')
plt.plot([murangeExp[0][0], murangeExp[1][0]],
[murangeExp[1][1]] * 2, 'm:')
plt.plot([murangeExp[1][0]] * 2,
[murangeExp[1][1], murangeExp[0][1]], 'm:')
plt.plot([murangeExp[1][0], murangeExp[0][0]],
[murangeExp[0][1]] * 2, 'm:')
plt.title("log10|eta(mu)|")
plt.colorbar(p)
plt.grid()
plt.show()
def plotInfSet3(murange, murangeEff, approx, mu0, nSamples = 200,
marginal = [None, None, .05], exps = [2., 1., 1.], clip = -1,
relative = True, normalizeDen = False, nobeta = False):
mNone = [i for i, m in enumerate(marginal) if m is None]
nNone = len(mNone)
if nNone not in [1, 2]:
raise
if nNone == 1:
exp = exps[mNone[0]]
mu1 = np.linspace(murangeEff[0][mNone[0]] ** exp,
murangeEff[1][mNone[0]] ** exp, nSamples)
mu1s = np.power(mu1, 1. / exp)
mus = []
mBase = copy(marginal)
for m1 in mu1s:
mBase[mNone[0]] = m1
mus += [copy(mBase)]
Z = approx.normHF(mus)
T = approx.normApprox(mus)
R = approx.normRes(mus)
E = approx.normErr(mus)
if relative:
F = approx.normRHS(mus)
R /= F
E /= Z
if normalizeDen:
Qvals = np.abs(approx.trainedModel.getQVal(mus))
R *= Qvals
E *= Qvals
if nobeta:
beta = np.empty(len(mus))
beta[:] = np.nan
else:
beta = approx.HFEngine.stabilityFactor(mus, approx.getHF(mus))
plotInfSet31FromData(mNone, mus, Z, T, R, E, beta, murange, approx,
mu0, exp, normalizeDen)
return mus, Z, T, R, E, beta
if nNone == 2:
exps = [exps[m] for m in mNone]
mu1 = np.linspace(murangeEff[0][mNone[0]] ** exps[0],
murangeEff[1][mNone[0]] ** exps[0], nSamples)
mu2 = np.linspace(murangeEff[0][mNone[1]] ** exps[1],
murangeEff[1][mNone[1]] ** exps[1], nSamples)
mu1s = np.power(mu1, 1. / exps[0])
mu2s = np.power(mu2, 1. / exps[1])
Mu1, Mu2 = np.meshgrid(np.real(mu1), np.real(mu2))
mus = []
mBase = copy(marginal)
for m2 in mu2s:
mBase[mNone[1]] = m2
for m1 in mu1s:
mBase[mNone[0]] = m1
mus += [copy(mBase)]
Ze = approx.normHF(mus).reshape((nSamples, nSamples))
Te = approx.normApprox(mus).reshape((nSamples, nSamples))
Re = approx.normRes(mus).reshape((nSamples, nSamples))
Ee = approx.normErr(mus).reshape((nSamples, nSamples))
if relative:
Fe = approx.normRHS(mus).reshape((nSamples, nSamples))
Re /= Fe
Ee /= Ze
if normalizeDen:
Qvals = np.abs(approx.trainedModel.getQVal(mus).reshape(
(nSamples, nSamples)))
Re *= Qvals
Ee *= Qvals
if nobeta:
betae = np.empty((nSamples, nSamples))
betae[:, :] = np.nan
else:
betae = approx.HFEngine.stabilityFactor(mus, approx.getHF(mus))\
.reshape((nSamples, nSamples))
plotInfSet32FromData(mNone, mus, Ze, Te, Re, Ee, betae, murange,
approx, mu0, exps, clip, normalizeDen)
return mus, Ze, Te, Re, Ee, betae

Event Timeline