Plot-Angle-GammaV2_SubPlots.py

import numpy as np
import matplotlib.pyplot as plt
import sympy as sym
import math
import os
import subprocess
import fileinput
import re
import matlab.engine
import sys
from ClassifyMin import *
from HelperFunctions import *
# from CellScript import *
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.cm as cm
from vtk.util import numpy_support
from pyevtk.hl import gridToVTK
import time
import matplotlib.ticker as ticker

import matplotlib as mpl
from matplotlib.ticker import MultipleLocator,FormatStrFormatter,MaxNLocator
import pandas as pd

# from matplotlib import rc
# rc('text', usetex=True) # Use LaTeX font
#
# import seaborn as sns
# sns.set(color_codes=True)


def format_func(value, tick_number):
    # # find number of multiples of pi/2
    # N = int(np.round(2 * value / np.pi))
    # if N == 0:
    #     return "0"
    # elif N == 1:
    #     return r"$\pi/2$"
    # elif N == 2:
    #     return r"$\pi$"
    # elif N % 2 > 0:
    #     return r"${0}\pi/2$".format(N)
    # else:
    #     return r"${0}\pi$".format(N // 2)
    # find number of multiples of pi/2
    N = int(np.round(4 * value / np.pi))
    if N == 0:
        return "0"
    elif N == 1:
        return r"$\pi/4$"
    elif N == 2:
        return r"$\pi/2$"
    elif N % 2 > 0:
        return r"${0}\pi/2$".format(N)
    else:
        return r"${0}\pi$".format(N // 2)



def find_nearest(array, value):
    array = np.asarray(array)
    idx = (np.abs(array - value)).argmin()
    return array[idx]


def find_nearestIdx(array, value):
    array = np.asarray(array)
    idx = (np.abs(array - value)).argmin()
    return idx


InputFile  = "/inputs/computeMuGamma.parset"
OutputFile = "/outputs/outputMuGamma.txt"
# --------- Run  from src folder:
path_parent = os.path.dirname(os.getcwd())
os.chdir(path_parent)
path = os.getcwd()
print(path)
InputFilePath = os.getcwd()+InputFile
OutputFilePath = os.getcwd()+OutputFile
print("InputFilepath: ", InputFilePath)
print("OutputFilepath: ", OutputFilePath)
print("Path: ", path)

print('---- Input parameters: -----')
alpha = 10
mu1 = 1.0
rho1 = 1.0
beta = 2.0  #5.0
theta = 1.0/8.0
#

alpha = -0.5
beta = 40.0
theta= 1/8.0



# # INTERESTING! from pi/2:
alpha = -0.5
beta = 40.0
theta= 1/8.0
#
# # # INTERESTING! from pi/2:
# alpha = -0.2
# beta = 25.0
# theta= 1/2

# INTERESTING!:
# alpha = -0.5
# beta = 5.0
# theta= 1/30



# INTERESTING!:
# alpha = -0.25
# beta = 10.0
# theta= 3/4


# # INTERESTING!:
alpha = -0.25
beta = 10.0
theta= 1/8

#
# INTERESTING!:
# alpha = -0.25
# beta = 5.0
# theta= 1/8
#


# # INTERESTING!:
alpha = -0.5
beta = 10.0
theta= 1/8



alpha_1 = -1.0
alpha_2 = -0.75
alpha_3 = -0.70

# alpha_1 = -1
# alpha_2 = -0.5
# alpha_3 =  -0.25

angles_1 = []
angles_2 = []
angles_3 = []

beta = 2.0
theta= 0.25
# beta = 10.0
# theta= 0.5




print('mu1: ', mu1)
print('rho1: ', rho1)
print('alpha_1: ', alpha_1)
print('alpha_2: ', alpha_2)
print('alpha_3: ', alpha_3)
print('beta: ', beta)
print('theta: ', theta)
# print('gamma:', gamma)
print('----------------------------')

# ----------------------------------------------------------------


gamma_min = 0.01
gamma_max = 1.5
# Gamma_Values = np.linspace(gamma_min, gamma_max, num=100)    # TODO variable Input Parameters...alpha,beta...
Gamma_Values = np.linspace(gamma_min, gamma_max, num=40)    # TODO variable Input Parameters...alpha,beta...
print('(Input) Gamma_Values:', Gamma_Values)
# mu_gamma = []

# Gamma_Values = '0'



# Get values for mu_Gamma
GetMuGammaVec = np.vectorize(GetMuGamma)
muGammas = GetMuGammaVec(beta,theta,Gamma_Values,mu1,rho1, InputFilePath ,OutputFilePath )
print('muGammas:', muGammas)

q12 = 0.0
q1 = (1.0/6.0)*harmonicMean(mu1, beta, theta)
q2 = (1.0/6.0)*arithmeticMean(mu1, beta, theta)
print('q1: ', q1)
print('q2: ', q2)
b1 = prestrain_b1(rho1, beta, alpha,theta)
b2 = prestrain_b2(rho1, beta, alpha,theta)
q3_star = math.sqrt(q1*q2)
print('q3_star:', q3_star)

# TODO these have to be compatible with input parameters!!!
# compute certain ParameterValues that this makes sense
# b1 = q3_star
# b2 = q1
print('b1: ', b1)
print('b2: ', b2)

# return classifyMin(q1, q2, q3, q12,  b1, b2,  print_Cases, print_Output)



# classifyMin_anaVec = np.vectorize(classifyMin_ana)
# G, angles, Types, curvature = classifyMin_anaVec(alpha, beta, theta, muGammas,  mu1, rho1)
classifyMin_anaVec = np.vectorize(classifyMin_ana)
G, angles_1, Types, curvature = classifyMin_anaVec(alpha_1, beta, theta, muGammas,  mu1, rho1)
G, angles_2, Types, curvature = classifyMin_anaVec(alpha_2, beta, theta, muGammas,  mu1, rho1)
G, angles_3, Types, curvature = classifyMin_anaVec(alpha_3, beta, theta, muGammas,  mu1, rho1)

# _,angles,_,_ = classifyMin_anaVec(alpha, beta, theta, muGammas,  mu1, rho1)

print('angles_1:', angles_1)
print('angles_2:', angles_2)
print('angles_3:', angles_3)


idx = find_nearestIdx(muGammas, q3_star)
print('GammaValue Idx closest to q_3^*', idx)
gammaClose = Gamma_Values[idx]
print('GammaValue(Idx) with mu_gamma closest to q_3^*', gammaClose)



determinantVec = np.vectorize(determinant)

detValues = determinantVec(q1,q2,muGammas,q12)
print('detValues:', detValues)


detZeroidx = find_nearestIdx(detValues, 0)
print('idx where det nearest to zero', idx)
gammaClose = Gamma_Values[detZeroidx]
print('gammaClose:', gammaClose)


# --- Convert to numpy array
Gamma_Values = np.array(Gamma_Values)


angles_1 = np.array(angles_1)
angles_2 = np.array(angles_2)
angles_3 = np.array(angles_3)

# ---------------- Create Plot -------------------
# plt.figure()


mpl.rcParams['text.usetex'] = True
mpl.rcParams["font.family"] = "serif"
mpl.rcParams["font.size"] = "9"
width = 6.28
height = width / 1.618
height = width / 2.5
fig = plt.figure()

fig,ax = plt.subplots(nrows=1,ncols=3,figsize=(width,height)) # more than one plot
# fig,ax = plt.subplots(nrows=1,ncols=3,figsize=(width,height),sharey=True) # Share Y-axis


fig = plt.figure()
gs = fig.add_gridspec(1,3, hspace=0.2, wspace=0.1)
ax = gs.subplots(sharey=True)

# ax = plt.axes((0.15,0.21 ,0.75,0.75))
# ax = plt.axes((0.15,0.21 ,0.8,0.75))
# ax.tick_params(axis='x',which='major', direction='out',pad=5)
# ax.tick_params(axis='y',which='major', length=3, width=1, direction='out',pad=3)
# ax.xaxis.set_major_locator(MultipleLocator(0.1))
# ax.xaxis.set_minor_locator(MultipleLocator(0.05))
ax[0].yaxis.set_major_locator(plt.MultipleLocator(np.pi / 8))
ax[0].yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 16))
ax[0].yaxis.set_major_formatter(plt.FuncFormatter(format_func))
ax[1].yaxis.set_major_locator(plt.MultipleLocator(np.pi / 8))
ax[1].yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 16))
ax[1].yaxis.set_major_formatter(plt.FuncFormatter(format_func))
ax[2].yaxis.set_major_locator(plt.MultipleLocator(np.pi / 8))
ax[2].yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 16))
ax[2].yaxis.set_major_formatter(plt.FuncFormatter(format_func))

ax[0].grid(True,which='major',axis='both',alpha=0.3)
ax[1].grid(True,which='major',axis='both',alpha=0.3)
ax[2].grid(True,which='major',axis='both',alpha=0.3)
# ax.grid(True,which='major',axis='both',alpha=0.3)


# # plt.title(r'angle$-\mu_\gamma(\gamma)$-Plot')
# # plt.title(r'angle$-\gamma$-Plot')
# plt.plot(Gamma_Values, angles)
# plt.scatter(Gamma_Values, angles)
# plt.plot(muGammas, angles)
# plt.scatter(muGammas, angles)
# # plt.axis([0, 6, 0, 20])
# # plt.axhline(y = 1.90476, color = 'b', linestyle = ':', label='$q_1$')
# # plt.axhline(y = 2.08333, color = 'r', linestyle = 'dashed', label='$q_2$')
# plt.axvline(x = q1, color = 'b', linestyle = ':', label='$q_1$')
# plt.axvline(x = q2, color = 'r', linestyle = 'dashed', label='$q_2$')
# # plt.axvline(x = q3_star, color = 'r', linestyle = 'dashed', label='$\gamma^*$')

ax[0].plot(Gamma_Values, angles_1, 'royalblue', zorder=3, )
ax[1].plot(Gamma_Values, angles_2, 'royalblue', zorder=3, )
ax[2].plot(Gamma_Values, angles_3, 'royalblue', zorder=3, )

ax[0].set_xlabel(r"$\gamma$")
ax[0].set_ylabel(r"angle  $\alpha$")
ax[0].xaxis.set_minor_locator(MultipleLocator(0.25))
ax[0].xaxis.set_major_locator(MultipleLocator(0.5))

ax[1].set_xlabel(r"$\gamma$")
# ax[1].set_ylabel(r"angle  $\alpha$")
ax[1].xaxis.set_minor_locator(MultipleLocator(0.25))
ax[1].xaxis.set_major_locator(MultipleLocator(0.5))

ax[2].set_xlabel(r"$\gamma$")
# ax[2].set_ylabel(r"angle  $\alpha$")
ax[2].xaxis.set_minor_locator(MultipleLocator(0.25))
ax[2].xaxis.set_major_locator(MultipleLocator(0.5))
# Labels to use in the legend for each line
line_labels = [r"$\theta_\mu  = 1.0$", r"$\theta_\mu  = 2.0$",  r"$\theta_\mu  = 5.0$", r"$\theta_\mu  = 10.0$"]


# ax.plot(Gamma_Values, angles, 'royalblue', zorder=3, )
# ax.set_xlabel(r"$\gamma$")
# ax.set_ylabel(r"angle $\alpha$")

# plt.xlabel("$q_3$")
# plt.xlabel("$\gamma$")
# plt.ylabel("angle")
# ax.grid(True)


# ax.yaxis.set_major_locator(plt.MultipleLocator(np.pi / 2))
# ax.yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 12))

# ax.yaxis.set_major_locator(plt.MultipleLocator(np.pi / 4))
# ax.yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 12))
# ax.yaxis.set_major_formatter(plt.FuncFormatter(format_func))

# ax.yaxis.set_major_formatter(ticker.FormatStrFormatter('%g $\pi$'))
# ax.yaxis.set_major_locator(ticker.MultipleLocator(base=0.25))


# ax.yaxis.set_major_formatter(ticker.FuncFormatter(
# lambda val,pos: '{:.0g}$\pi$'.format(2*val/np.pi) if val !=0 else '0'))
# ax.yaxis.set_major_locator(ticker.MultipleLocator(base=0.5*np.pi))

# ---------------------------- show pi values ------------------------------------
# ax.axvline(x = q1, color = 'b', linestyle = ':', label='$q_1$')
# ax.axvline(x = q2, color = 'r', linestyle = 'dashed', label='$q_2$')
# ax.legend()
# # ax.set(xlim=(1.750, 1.880), ylim=(0, math.pi/2.0))
# ax.set(xlim=(1.760, 1.880), ylim=(-0.1, np.pi/4.0))
# # ax.set_yticks([0,  np.pi/4 ,np.pi/2])
# # labels = ['$0$', r'$\pi/4$', r'$\pi/2$']
# ax.set_yticks([0, np.pi/8, np.pi/4 ])
# labels = ['$0$',r'$\pi/8$', r'$\pi/4$']
# ax.set_yticklabels(labels)
# ---------------------------------------------------------------


# ax.set(xlim=(1.750, 1.880), ylim=(0, math.pi/2.0))

# ax.set(xlim=(1.760, 1.880), ylim=(-0.1, np.pi/4.0))
# ax.set(xlim=(1.760, 1.880), ylim=(-0.1, np.pi/4.0))
# ax.set_yticks([0,  np.pi/4 ,np.pi/2])
# labels = ['$0$', r'$\pi/4$', r'$\pi/2$']


# OLD :
# ax.set_yticks([0, np.pi/8, np.pi/4 ])
# labels = ['$0$',r'$\pi/8$', r'$\pi/4$']

labels = ['$0$',r'$\pi/8$', r'$\pi/4$' ,r'$3\pi/8$' , r'$\pi/2$']
ax[0].set_yticks([0, np.pi/8, np.pi/4, 3*np.pi/8 , np.pi/2, ])
ax[1].set_yticks([0, np.pi/8, np.pi/4, 3*np.pi/8 , np.pi/2 ])
ax[2].set_yticks([0, np.pi/8, np.pi/4, 3*np.pi/8 , np.pi/2 ])

ax[0].set_yticklabels(labels)
ax[1].set_yticklabels(labels)
ax[2].set_yticklabels(labels)

for i in range(3):
    ax[i].set_ylim([0-0.1, np.pi/2+0.1])

# Plot Gamma Value that is closest to q3_star
l1 = ax[0].axvline(x = gammaClose, color = 'midnightblue', linestyle = 'dashed', linewidth=1, label='$\gamma^*$')
l2 = ax[1].axvline(x = gammaClose, color = 'midnightblue', linestyle = 'dashed', linewidth=1, label='$\gamma^*$')
l3 = ax[2].axvline(x = gammaClose, color = 'midnightblue', linestyle = 'dashed', linewidth=1, label='$\gamma^*$')
# color elliptic/hyperbolic region
# ax.axvspan(gamma_min, gammaClose, color='red', alpha=0.2)
# ax.axvspan(gammaClose, gamma_max, color='green', alpha=0.2)


# ax[0].legend(loc='upper right')
# ax[2].legend(l3,
#             labels=['$\gamma^*$'],
#             loc='upper left',
#             bbox_to_anchor=(1,1))

# ax[2].legend(
#             loc='upper left',
#             bbox_to_anchor=(1,1))

# fig.legend([l1],     # The line objects
#            # label='$\gamma^*$',   # The labels for each line
#            loc="center right",   # Position of legend
#            borderaxespad=0.15    # Small spacing around legend box
#            # title="Legend Title"  # Title for the legend
#            )


line_labels = [r"$\gamma^*$"]

fig.legend([l1], [r"$\gamma^*$"],
            # bbox_to_anchor=[0.5, 0.92],
            bbox_to_anchor=[0.5, 0.94],
            loc='center', ncol=3)



# plt.subplots_adjust(wspace=0.4, hspace=0.0)
# plt.tight_layout()

# Adjust the scaling factor to fit your legend text completely outside the plot
# (smaller value results in more space being made for the legend)
# plt.subplots_adjust(right=0.9)
plt.subplots_adjust(bottom=0.2)

fig.set_size_inches(width, height)
fig.savefig('Plot-Angle-Gamma.pdf')

plt.show()




# plt.figure()
# plt.title(r'angle$-\mu_\gamma(\gamma)$-Plot')
# plt.plot(muGammas, angles)
# plt.scatter(muGammas, angles)
# # plt.axis([0, 6, 0, 20])
# # plt.axhline(y = 1.90476, color = 'b', linestyle = ':', label='$q_1$')
# # plt.axhline(y = 2.08333, color = 'r', linestyle = 'dashed', label='$q_2$')
# plt.axvline(x = 1.90476, color = 'b', linestyle = ':', label='$q_1$')
# plt.axvline(x = 2.08333, color = 'r', linestyle = 'dashed', label='$q_2$')
# plt.xlabel("$\mu_\gamma$")
# plt.ylabel("angle")
# plt.legend()
# plt.show()
#