From 745772707687d5c548b5788565a27828e315681d Mon Sep 17 00:00:00 2001 From: Klaus <klaus.boehnlein@tu-dresden.de> Date: Fri, 15 Oct 2021 00:30:29 +0200 Subject: [PATCH] Update Plot methods --- src/Plot-Angle-Gamma.py | 220 ++++++++++++--------- src/Plot-Angle-q3.py | 387 +++++++++++++++++++++++++++++++++++++ src/Plot_Angle_Lemma1.4.py | 121 ++++++++---- src/Plotq12.py | 66 +++++-- src/plot-q3-gamma.py | 174 +++++++++++++++++ 5 files changed, 824 insertions(+), 144 deletions(-) create mode 100644 src/Plot-Angle-q3.py create mode 100644 src/plot-q3-gamma.py diff --git a/src/Plot-Angle-Gamma.py b/src/Plot-Angle-Gamma.py index 419f19b1..e4117d35 100644 --- a/src/Plot-Angle-Gamma.py +++ b/src/Plot-Angle-Gamma.py @@ -30,14 +30,26 @@ import pandas as pd 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(2 * value / np.pi)) + N = int(np.round(4 * value / np.pi)) if N == 0: return "0" elif N == 1: - return r"$\pi/2$" + return r"$\pi/4$" elif N == 2: - return r"$\pi$" + return r"$\pi/2$" elif N % 2 > 0: return r"${0}\pi/2$".format(N) else: @@ -45,7 +57,6 @@ def format_func(value, tick_number): - def find_nearest(array, value): array = np.asarray(array) idx = (np.abs(array - value)).argmin() @@ -84,7 +95,6 @@ print('---- Input parameters: -----') alpha = 10.0 # beta = 40.0 # theta = 0.125 - mu1 = 10.0 rho1 = 1.0 # alpha = 2.0 @@ -166,15 +176,13 @@ print('----------------------------') gamma_min = 0.01 gamma_max = 1.0 -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=200) # TODO variable Input Parameters...alpha,beta... print('(Input) Gamma_Values:', Gamma_Values) # mu_gamma = [] # Gamma_Values = '0' -# # --- Options -# # make_Plot = False -make_Plot = True + # Get values for mu_Gamma GetMuGammaVec = np.vectorize(GetMuGamma) @@ -233,111 +241,141 @@ gammaClose = Gamma_Values[detZeroidx] print('gammaClose:', gammaClose) -# Make Plot -if make_Plot: - # plt.figure() - # - # - # - # - # # 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^*$') +# --- Convert to numpy array +Gamma_Values = np.array(Gamma_Values) +angles = np.array(angles) + +# ---------------- Create Plot ------------------- +# plt.figure() + + +mpl.rcParams['text.usetex'] = True +mpl.rcParams["font.family"] = "serif" +mpl.rcParams["font.size"] = "9" +width = 6.28 *0.5 +height = width / 1.618 +fig = plt.figure() +# 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.yaxis.set_major_locator(plt.MultipleLocator(np.pi / 8)) +ax.yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 16)) +ax.yaxis.set_major_formatter(plt.FuncFormatter(format_func)) +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^*$') - f,ax=plt.subplots(1) - # ax.plot(muGammas, angles) - # ax.scatter(muGammas, angles) - ax.plot(Gamma_Values, angles, zorder=3) - ax.scatter(Gamma_Values, angles) +# f,ax=plt.subplots(1) - plt.xlabel("$q_3$") - plt.xlabel("$\gamma$") - plt.ylabel("angle") - ax.grid(True) +# ax.plot(muGammas, angles) +# ax.scatter(muGammas, angles) - # ax.yaxis.set_major_locator(plt.MultipleLocator(np.pi / 2)) - # ax.yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 12)) +ax.plot(Gamma_Values, angles, 'royalblue', zorder=3, ) +# ax.scatter(Gamma_Values, angles) - # 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.set_xlabel(r"$q_3(\gamma)$") +ax.set_xlabel(r"$\gamma$") +ax.set_ylabel(r"angle $\angle$") - # ax.yaxis.set_major_formatter(ticker.FormatStrFormatter('%g $\pi$')) - # ax.yaxis.set_major_locator(ticker.MultipleLocator(base=0.25)) +# plt.xlabel("$q_3$") +# plt.xlabel("$\gamma$") +# plt.ylabel("angle") +# ax.grid(True) - # 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)) +# ax.yaxis.set_major_locator(plt.MultipleLocator(np.pi / 2)) +# ax.yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 12)) - # ---------------------------- 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.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.legend() - # ax.set(xlim=(1.750, 1.880), ylim=(0, math.pi/2.0)) +# ax.yaxis.set_major_formatter(ticker.FormatStrFormatter('%g $\pi$')) +# ax.yaxis.set_major_locator(ticker.MultipleLocator(base=0.25)) - # 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$'] - ax.set_yticks([0, np.pi/8, np.pi/4 ]) - labels = ['$0$',r'$\pi/8$', r'$\pi/4$'] - ax.set_yticklabels(labels) +# 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)) - # Plot Gamma Value that is closest to q3_star - ax.axvline(x = gammaClose, color = 'b', linestyle = 'dashed', label='$\gamma^*$') - ax.axvspan(gamma_min, gammaClose, color='red', alpha=0.5) - ax.axvspan(gammaClose, gamma_max, color='green', alpha=0.5) +# ---------------------------- 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$'] - # plt.xlabel("$q_3(\gamma)$") - plt.xlabel("$\gamma$") - plt.ylabel("angle") - plt.legend(loc='upper center') - plt.show() +# OLD : +ax.set_yticks([0, np.pi/8, np.pi/4 ]) +labels = ['$0$',r'$\pi/8$', r'$\pi/4$'] +ax.set_yticklabels(labels) +# Plot Gamma Value that is closest to q3_star +ax.axvline(x = gammaClose, color = 'midnightblue', linestyle = 'dashed', label='$\gamma^*$') +ax.axvspan(gamma_min, gammaClose, color='red', alpha=0.2) +ax.axvspan(gammaClose, gamma_max, color='green', alpha=0.2) - # 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() - # +ax.legend(loc='upper right') + +# plt.xlabel("$q_3(\gamma)$") +# plt.xlabel("$\gamma$") +# plt.ylabel("angle") +# plt.legend(loc='upper center') + +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() +# diff --git a/src/Plot-Angle-q3.py b/src/Plot-Angle-q3.py new file mode 100644 index 00000000..c44f1f2d --- /dev/null +++ b/src/Plot-Angle-q3.py @@ -0,0 +1,387 @@ +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: -----') +# mu1 = 10.0 +# rho1 = 1.0 +# alpha = 2.56140350877193 #2.56140350877193, 4.0852130325814535 +# beta = 2.0 #5.0 +# theta = 1.0/4.0 +# theta = 1.0/8.0 # 0.5 +# theta = 0.075 # 0.5 +# mu1 = 10.0 +# rho1 = 1.0 +alpha = 10.0 +# beta = 40.0 +# theta = 0.125 +mu1 = 10.0 +rho1 = 1.0 +# alpha = 2.0 +beta = 2.0 #5.0 +# theta = 1.0/4.0 +theta = 1.0/8.0 +# + + + +# mu1 = 10.0 +# rho1 = 1.0 +# alpha = 10.0 +# beta = 2.0 #5.0 +# theta = 1.0/8.0 +# # + + +# mu1 = 10.0 +# rho1 = 1.0 +# # alpha = 10.02021333 +# alpha = 10.0 +# beta = 2.0 +# theta = 0.125 + + + +# mu1 = 10.0 +# rho1 = 1.0 +# # alpha = 10.02021333 +# alpha = 9.0 +# beta = 2.0 +# theta = 0.075 + +# +# mu1 = 10.0 +# rho1 = 1.0 +# alpha = 4.0 +# beta = 10.0 +# theta = 1.0/4.0 + +# alpha = 10 #1.05263158 +# beta = 40.0 #5.0 +# # theta = 1.0/4.0 +# theta = 1.0/8.0 # 0.5 +# +# +# alpha = 2.0 +# beta = 2.0 #5.0 +# theta = 1/4.0 +# rho1 = 1.0 +# mu1 = 10.0 + +# InterestingParameterSet : +# mu1 = 10.0 +# rho1 = 1.0 +# alpha = 10 +# beta = 40.0 + + +# theta = 0.124242 +# gamma = 0.75 +#set gamma either to 1. '0' 2. 'infinity' or 3. a numerical positive value +# gamma = '0' +# gamma = 'infinity' +# gamma = 0.5 +# gamma = 0.25 + +print('mu1: ', mu1) +print('rho1: ', rho1) +print('alpha: ', alpha) +print('beta: ', beta) +print('theta: ', theta) +# print('gamma:', gamma) +print('----------------------------') + +# ---------------------------------------------------------------- + + +gamma_min = 0.01 +gamma_max = 3.0 +Gamma_Values = np.linspace(gamma_min, gamma_max, num=200) # 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) + + +print('ratio(left) =', b2/b1) +print('ratio(right) = ', q1/q3_star) + +# 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, Types, curvature = classifyMin_anaVec(alpha, beta, theta, muGammas, mu1, rho1) + +# _,angles,_,_ = classifyMin_anaVec(alpha, beta, theta, muGammas, mu1, rho1) + +print('angles:', angles) + + + + +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 = np.array(angles) + +# ---------------- Create Plot ------------------- +# plt.figure() + + +mpl.rcParams['text.usetex'] = True +mpl.rcParams["font.family"] = "serif" +mpl.rcParams["font.size"] = "9" +width = 6.28 *0.5 +height = width / 1.618 +fig = plt.figure() +# 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.02)) +ax.xaxis.set_minor_locator(MultipleLocator(0.01)) +ax.yaxis.set_major_locator(plt.MultipleLocator(np.pi / 8)) +ax.yaxis.set_minor_locator(plt.MultipleLocator(np.pi / 16)) +ax.yaxis.set_major_formatter(plt.FuncFormatter(format_func)) +ax.grid(True,which='major',axis='both',alpha=0.3) +ax.set(xlim=(1.760, 1.880), ylim=(-0.1, np.pi/4.0)) + + +# # 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^*$') + + + +# f,ax=plt.subplots(1) + + +# ax.plot(muGammas, angles) +# ax.scatter(muGammas, angles) + +ax.plot(muGammas, angles, 'royalblue', zorder=3, ) +# ax.scatter(Gamma_Values, angles) +ax.set_xlabel(r"$q_3(\gamma)$") +# ax.set_xlabel(r"$\gamma$") +ax.set_ylabel(r"angle $\angle$") + +# 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$'] +ax.set_yticklabels(labels) + + +# Plot Gamma Value that is closest to q3_star +# ax.axvline(x = gammaClose, color = 'midnightblue', linestyle = 'dashed', label='$\gamma^*$') +# ax.axvspan(gamma_min, gammaClose, color='red', alpha=0.2) +# ax.axvspan(gammaClose, gamma_max, color='green', alpha=0.2) + +ax.axvline(x = q1, color = 'forestgreen', linestyle = 'dashed', label='$q_1$') +ax.axvline(x = q2, color = 'darkorange', linestyle = 'dashed', label='$q_2$') + + +ax.legend(loc='upper center') + +# plt.xlabel("$q_3(\gamma)$") +# plt.xlabel("$\gamma$") +# plt.ylabel("angle") +# plt.legend(loc='upper center') + +fig.set_size_inches(width, height) +fig.savefig('Plot-Angle-q3.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() +# diff --git a/src/Plot_Angle_Lemma1.4.py b/src/Plot_Angle_Lemma1.4.py index 18fe424d..de8e97c2 100644 --- a/src/Plot_Angle_Lemma1.4.py +++ b/src/Plot_Angle_Lemma1.4.py @@ -15,9 +15,9 @@ import matplotlib as mpl from matplotlib.ticker import MultipleLocator,FormatStrFormatter,MaxNLocator import pandas as pd -import tikzplotlib -# from pylab import * -from tikzplotlib import save as tikz_save +# import tikzplotlib +# # from pylab import * +# from tikzplotlib import save as tikz_save # Needed ? @@ -45,40 +45,40 @@ mpl.use('pdf') # To set the figure size we construct a function # to convert from pts to inches and to determine # an aesthetic figure height using the golden ratio: -def set_size(width, fraction=1): - """Set figure dimensions to avoid scaling in LaTeX. - - Parameters - ---------- - width: float - Document textwidth or columnwidth in pts - fraction: float, optional - Fraction of the width which you wish the figure to occupy - - Returns - ------- - fig_dim: tuple - Dimensions of figure in inches - """ - # Width of figure (in pts) - fig_width_pt = width * fraction - - # Convert from pt to inches - inches_per_pt = 1 / 72.27 - - # Golden ratio to set aesthetic figure height - # https://disq.us/p/2940ij3 - golden_ratio = (5**.5 - 1) / 2 - - # Figure width in inches - fig_width_in = fig_width_pt * inches_per_pt - # Figure height in inches - fig_height_in = fig_width_in * golden_ratio - - fig_dim = (fig_width_in, fig_height_in) - - return fig_dim - +# def set_size(width, fraction=1): +# """Set figure dimensions to avoid scaling in LaTeX. +# +# Parameters +# ---------- +# width: float +# Document textwidth or columnwidth in pts +# fraction: float, optional +# Fraction of the width which you wish the figure to occupy +# +# Returns +# ------- +# fig_dim: tuple +# Dimensions of figure in inches +# """ +# # Width of figure (in pts) +# fig_width_pt = width * fraction +# +# # Convert from pt to inches +# inches_per_pt = 1 / 72.27 +# +# # Golden ratio to set aesthetic figure height +# # https://disq.us/p/2940ij3 +# golden_ratio = (5**.5 - 1) / 2 +# +# # Figure width in inches +# fig_width_in = fig_width_pt * inches_per_pt +# # Figure height in inches +# fig_height_in = fig_width_in * golden_ratio +# +# fig_dim = (fig_width_in, fig_height_in) +# +# return fig_dim +# @@ -165,7 +165,7 @@ lambda1 = 0.0 gamma = 1.0/4.0 gamma = 'infinity' #Elliptic Setting -# gamma = '0' #Hyperbolic Setting +gamma = '0' #Hyperbolic Setting # gamma = 0.5 @@ -222,7 +222,7 @@ xmax = 0.41 # xmin = 0.01 # xmax = 3.0 numPoints = 200 -# numPoints = 101 +# numPoints = 50 X_Values = np.linspace(xmin, xmax, num=numPoints) print(X_Values) @@ -600,7 +600,7 @@ if gamma == '0': for x in jump_xValues: - plt.axvline(x,ymin=0, ymax= 1, color = 'orange',alpha=0.5, linestyle = 'dashed') + plt.axvline(x,ymin=0, ymax= 1, color = 'orange',alpha=0.5, linestyle = 'dashed', linewidth=1) # plt.axvline(x,ymin=0, ymax= 1, color = 'orange',alpha=0.5, linestyle = 'dashed', label=r'$\theta_*$') # plt.axvline(x_plotValues[1],ymin=0, ymax= 1, color = 'g',alpha=0.5, linestyle = 'dashed') @@ -629,8 +629,47 @@ for x in jump_xValues: +# Find transition point +lastIdx = len(Y_Values)-1 + +for idx, y in enumerate(Y_Values): + if idx != lastIdx: + if abs(y-0) < 0.01 and abs(Y_Values[idx+1] - 0) > 0.05: + transition_point1 = X_Values[idx+1] + print('transition point1:', transition_point1 ) + if abs(y-0.5*np.pi) < 0.01 and abs(Y_Values[idx+1] -0.5*np.pi)>0.01: + transition_point2 = X_Values[idx] + print('transition point2:', transition_point2 ) + if abs(y-0) > 0.01 and abs(Y_Values[idx+1] - 0) < 0.01: + transition_point3 = X_Values[idx+1] + print('transition point3:', transition_point3 ) + +# Add transition Points: +if gamma == '0': + ax.scatter([transition_point1, transition_point2],[np.pi/2,np.pi/2],s=6, marker='o', cmap=None, norm=None, facecolor = 'black', + edgecolor = 'black', vmin=None, vmax=None, alpha=None, linewidths=None, zorder=3) + + ax.text(transition_point1-0.02 , np.pi/2-0.02, r"$1$", size=6, bbox=dict(boxstyle="circle",facecolor='white', alpha=1.0, pad=0.1, linewidth=0.5) + ) + + ax.text(transition_point2+0.012 , np.pi/2-0.02, r"$2$", size=6, bbox=dict(boxstyle="circle",facecolor='white', alpha=1.0, pad=0.1, linewidth=0.5) + ) +else: + ax.scatter([transition_point1, transition_point2, transition_point3 ],[np.pi/2,np.pi/2,0 ],s=6, marker='o', cmap=None, norm=None, facecolor = 'black', + edgecolor = 'black', vmin=None, vmax=None, alpha=None, linewidths=None, zorder=3) + + ax.text(transition_point1-0.02 , np.pi/2-0.02, r"$1$", size=6, bbox=dict(boxstyle="circle",facecolor='white', alpha=1.0, pad=0.1, linewidth=0.5) + ) + + ax.text(transition_point2 +0.015 , np.pi/2-0.02, r"$2$", size=6, bbox=dict(boxstyle="circle",facecolor='white', alpha=1.0, pad=0.1, linewidth=0.5) + ) + + ax.text(transition_point3 +0.005 , 0+0.06, r"$3$", size=6, bbox=dict(boxstyle="circle",facecolor='white', alpha=1.0, pad=0.1, linewidth=0.5) + ) + + fig.set_size_inches(width, height) -fig.savefig('plot.pdf') +fig.savefig('Plot-Angle-Theta.pdf') diff --git a/src/Plotq12.py b/src/Plotq12.py index a2b8fac8..403a7b66 100644 --- a/src/Plotq12.py +++ b/src/Plotq12.py @@ -10,6 +10,11 @@ import matlab.engine # from subprocess import Popen, PIPE #import sys +import matplotlib.ticker as tickers +import matplotlib as mpl +from matplotlib.ticker import MultipleLocator,FormatStrFormatter,MaxNLocator +import pandas as pd + ###################### Plotq12.py ######################### # Input: Lame-Parameter Lambda # compute q12 entry of Q_hom for each Lambda @@ -43,7 +48,7 @@ q12 = [] RUN = True # RUN = False # make_Plot = False -make_Plot = True # vll besser : Plot_muGamma +# make_Plot = True # vll besser : Plot_muGamma if RUN: for lambda1 in Lambda_Values: @@ -56,7 +61,7 @@ if RUN: f.write(filedata) f.close() # Run Cell-Problem - subprocess.run(['./build-cmake/src/dune-microstructure', './inputs/cellsolver.parset'], + subprocess.run(['./build-cmake/src/Cell-Problem', './inputs/cellsolver.parset'], capture_output=True, text=True) #Extract mu_gamma from Output-File with open(OutputFilePath, 'r') as file: @@ -72,15 +77,52 @@ if RUN: # Make Plot -if make_Plot: - plt.figure() - plt.title(r'$q_{12}(\lambda)$-Plot') # USE MATHEMATICAL EXPRESSIONS IN TITLE - plt.plot(Lambda_Values, q12) - plt.scatter(Lambda_Values, q12) - # plt.axis([0, 6, 0, 20]) - plt.xlabel("$\lambda$") - plt.ylabel("$q_{12}$") - plt.legend() - plt.show() +# if make_Plot: + +# --- Convert to numpy array +Lambda_Values = np.array(Lambda_Values) +q12= np.array(q12) + +# ---------------- Create Plot ------------------- +mpl.rcParams['text.usetex'] = True +mpl.rcParams["font.family"] = "serif" +mpl.rcParams["font.size"] = "9" +width = 6.28 *0.5 +height = width / 1.618 +fig = plt.figure() +# 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(2)) +ax.xaxis.set_minor_locator(MultipleLocator(1)) +ax.yaxis.set_major_locator(MultipleLocator(0.1)) +# ax.yaxis.set_minor_locator(MultipleLocator(2)) +ax.grid(True,which='major',axis='both',alpha=0.3) + + +ax.plot(Lambda_Values, q12, 'royalblue', # data +marker='o', # each marker will be rendered as a circle +markersize=4, # marker size +markerfacecolor='orange', # marker facecolor +markeredgecolor='black', # marker edgecolor +markeredgewidth=1, # marker edge width +# linestyle='--', # line style will be dash line +linewidth=1) # line width +ax.set_xlabel(r"$\lambda$") +ax.set_ylabel(r"$q_{12}$") + +fig.set_size_inches(width, height) +fig.savefig('Plot-q12.pdf') + +# plt.figure() +# plt.title(r'$q_{12}(\lambda)$-Plot') # USE MATHEMATICAL EXPRESSIONS IN TITLE +# plt.plot(Lambda_Values, q12) +# plt.scatter(Lambda_Values, q12) +# # plt.axis([0, 6, 0, 20]) +# plt.xlabel("$\lambda$") +# plt.ylabel("$q_{12}$") +# plt.legend() +plt.show() # #--------------------------------------------------------------- diff --git a/src/plot-q3-gamma.py b/src/plot-q3-gamma.py new file mode 100644 index 00000000..9d44bdf7 --- /dev/null +++ b/src/plot-q3-gamma.py @@ -0,0 +1,174 @@ +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 +from HelperFunctions import * +from ClassifyMin import * + +import matplotlib.ticker as tickers +import matplotlib as mpl +from matplotlib.ticker import MultipleLocator,FormatStrFormatter,MaxNLocator +import pandas as pd +# from subprocess import Popen, PIPE +#import sys + +###################### plot-q3-gamma .py ######################### +# ..TODO Documentation ... +# Generalized Plot-Script giving the option to define +# quantity of interest and the parameter it depends on +# to create a plot +# +# Input: Define y & x for "x-y plot" as Strings +# - Run the 'Cell-Problem' for the different Parameter-Points +# (alternatively run 'Compute_MuGamma' if quantity of interest +# is q3=muGamma for a significant Speedup) + +########################################################### + + +# ----- Setup Paths ----- + +InputFile = "/inputs/computeMuGamma.parset" +OutputFile = "/outputs/outputMuGamma.txt" +# path = os.getcwd() +# InputFilePath = os.getcwd()+InputFile +# OutputFilePath = os.getcwd()+OutputFile +# --------- 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: -----') +mu1 = 10.0 +rho1 = 1.0 +# alpha = 2.8 +alpha = 2.0 +beta = 2.0 +theta = 1.0/4.0 +# lambda1 = 0.0 +gamma = 1.0/4.0 +# gamma = 'infinity' + +print('mu1: ', mu1) +print('rho1: ', rho1) +print('alpha: ', alpha) +print('beta: ', beta) +print('theta: ', theta) +print('gamma:', gamma) +print('----------------------------') + +yName = 'q3' +xName = 'gamma' +# --- Define Interval of x-values: +xmin = 0.01 +xmax = 3.0 +numPoints = 30 +X_Values = np.linspace(xmin, xmax, num=numPoints) +print(X_Values) +Y_Values = [] + +q1 = (1.0/6.0)*harmonicMean(mu1, beta, theta) +q2 = (1.0/6.0)*arithmeticMean(mu1, beta, theta) + +for x in X_Values: + # fist replace old parameters in parset + # SetParametersCellProblem(alpha,beta,theta,gamma,mu1,rho1,lambda1,InputFilePath) + SetParametersComputeMuGamma(beta,theta,gamma,mu1,rho1, InputFilePath) + + # replace the sought after x value in the parset + with open(path+"/inputs/computeMuGamma.parset", 'r') as file: + filedata = file.read() + filedata = re.sub('(?m)^'+xName+'=.*',xName+'='+str(x),filedata) + f = open(path+"/inputs/computeMuGamma.parset",'w') + f.write(filedata) + f.close() + #Run Compute_MuGamma since its much faster + print("Run Compute_MuGamma for " + xName + " =" , x) + subprocess.run(['./build-cmake/src/Compute_MuGamma', './inputs/computeMuGamma.parset'], + capture_output=True, text=True) + #Extract mu_gamma from Output-File + with open(path + '/outputs/outputMuGamma.txt', 'r') as file: + output = file.read() + tmp = re.search(r'(?m)^mu_gamma=.*',output).group() # Not necessary for Intention of Program t output Minimizer etc..... + s = re.findall(r"[-+]?\d*\.\d+|\d+", tmp) + Y_Values.append(float(s[0])) + +# ------------end of for-loop ----------------- +print("(Output) Values of " + yName + ": ", Y_Values) +# ----------------end of if-statement ------------- + + + + + + +# --- Convert to numpy array +Y_Values = np.array(Y_Values) +X_Values = np.array(X_Values) + +# ---------------- Create Plot ------------------- +mpl.rcParams['text.usetex'] = True +mpl.rcParams["font.family"] = "serif" +mpl.rcParams["font.size"] = "9" +# width as measured in inkscape +width = 6.28 *0.5 +height = width / 1.618 +fig = plt.figure() +ax = plt.axes((0.15,0.18,0.8,0.8)) +ax.tick_params(axis='x',which='major', direction='out',pad=3) +ax.tick_params(axis='y',which='major', length=3, width=1, direction='out',pad=3) +ax.xaxis.set_major_locator(MultipleLocator(0.5)) +ax.xaxis.set_minor_locator(MultipleLocator(0.25)) +# ax.plot(X_Values, Y_Values) + +ax.plot(X_Values, Y_Values, # data +marker='o', # each marker will be rendered as a circle +markersize=4, # marker size +markerfacecolor='orange', # marker facecolor +markeredgecolor='black', # marker edgecolor +markeredgewidth=1, # marker edge width +# linestyle='--', # line style will be dash line +linewidth=1) # line width + + +ax.set_xlabel(r"$\gamma$") +ax.set_ylabel(r"$q_3(\gamma)$") + + + +# plt.figure() +# # plt.title(r''+ yName + '-Plot') +# plt.plot(X_Values, Y_Values) +# plt.scatter(X_Values, Y_Values) +# # plt.axis([0, 6, 0, 20]) +# plt.xlabel(xName) +# plt.ylabel(yName) + + +# plt.axhline(y = 1.90476, color = 'forestgreen', linestyle = ':', label='$q_1$') +# plt.axhline(y = 2.08333, color = 'darkorange', linestyle = 'dashed', label='$q_2$') +plt.axhline(y = q1, color = 'forestgreen', linestyle = 'dashed', label='$q_1$') +plt.axhline(y = q2, color = 'darkorange', linestyle = 'dashed', label='$q_2$') +ax.legend(loc = 'center right') +# plt.legend() + +fig.set_size_inches(width, height) +fig.savefig('Plot-q3-Gamma.pdf') + + + +plt.show() +# #--------------------------------------------------------------- -- GitLab