From 538533c0df24fa8ebff019ab4a389f367a0132a2 Mon Sep 17 00:00:00 2001 From: Klaus <klaus.boehnlein@tu-dresden.de> Date: Wed, 13 Oct 2021 20:42:59 +0200 Subject: [PATCH] recompute right value of a* and adjust dependencies --- Matlab-Programs/symSolveSystem.m | 11 +++++++++-- src/ClassifyMin.py | 8 +++++--- src/Plot_CurvatureLemma1.4.py | 23 ++++++++++++++++------- 3 files changed, 30 insertions(+), 12 deletions(-) diff --git a/Matlab-Programs/symSolveSystem.m b/Matlab-Programs/symSolveSystem.m index 6424b5e7..5dde81e4 100644 --- a/Matlab-Programs/symSolveSystem.m +++ b/Matlab-Programs/symSolveSystem.m @@ -2,9 +2,16 @@ syms q1 q2 q3 q12 b1 b2 a1 a2 -eqn1 = q1 * a1 + (q3 + (q12/2)) *a2 == -2*q1*b1 - q12*b2 +% eqn1 = q1 * a1 + (q3 + (q12/2)) *a2 == -2*q1*b1 - q12*b2 +% +% eqn2 = (q3 + (q12/2)) * a1 + q2 *a2 == -2*q2*b2 - q12*b1 + + +eqn1 = q1 * a1 + (q3 + (q12/2)) *a2 == q1*b1 + (q12*b2/2) + +eqn2 = (q3 + (q12/2)) * a1 + q2 *a2 == q2*b2 + (q12*b1/2) + -eqn2 = (q3 + (q12/2)) * a1 + q2 *a2 == -2*q2*b2 - q12*b1 [A,B] = equationsToMatrix([eqn1, eqn2], [a1,a2]) diff --git a/src/ClassifyMin.py b/src/ClassifyMin.py index 02ce9e28..b4e83546 100644 --- a/src/ClassifyMin.py +++ b/src/ClassifyMin.py @@ -177,7 +177,9 @@ def classifyMin(q1, q2, q3, q12, b1, b2, print_Cases=False, print_Output=False) # is the “exact†solution of the equation. Else, x minimizes the # Euclidean 2-norm || b-Ax ||. If there are multiple minimizing solutions, # the one with the smallest 2-norm is returned. "" - a = np.linalg.lstsq(A, B)[0] # TODO check is this Ok ? + + # a = np.linalg.lstsq(A, B)[0] # TODO check is this Ok ? + a = np.linalg.lstsq(A, -B/2)[0] # TODO check is this Ok ? (UPDATE 13-10-21) print("Solution LGS: a =", a) a1 = a[0] a2 = a[1] @@ -211,9 +213,9 @@ def classifyMin(q1, q2, q3, q12, b1, b2, print_Cases=False, print_Output=False) # ------------------------------------ Elliptic Case ----------------------------------- if determinant >= epsilon: if print_Cases: print('E : elliptic case (determinant greater zero)') - a1_star = -(2*(b1*(q12**2) + 2*b1*q3*q12 - 4*b1*q1*q2 + 4*b2*q2*q3)) / \ + a1_star = (b1*(q12**2) + 2*b1*q3*q12 - 4*b1*q1*q2 + 4*b2*q2*q3) / \ (4*(q3**2) + 4*q3*q12 + (q12**2) - 4*q1*q2) - a2_star = -(2*(b2*(q12**2) + 2*b2*q3*q12 + 4*b1*q1*q3 - 4*b2*q1*q2)) / \ + a2_star = (b2*(q12**2) + 2*b2*q3*q12 + 4*b1*q1*q3 - 4*b2*q1*q2) / \ (4*(q3**2) + 4*q3*q12 + (q12**2) - 4*q1*q2) prod = a1_star*a2_star diff --git a/src/Plot_CurvatureLemma1.4.py b/src/Plot_CurvatureLemma1.4.py index ea4f16b7..f2c9f2de 100644 --- a/src/Plot_CurvatureLemma1.4.py +++ b/src/Plot_CurvatureLemma1.4.py @@ -118,13 +118,14 @@ print('gamma:', gamma) print('----------------------------') -# TODO? : Ask User for Input ... +# Optional - TODO? : +# -Ask User for Input ... # function = input("Enter value you want to plot (y-value):\n") # print(f'You entered {function}') # parameter = input("Enter Parameter this value depends on (x-value) :\n") # print(f'You entered {parameter}') -# Add Option to change NumberOfElements used for computation of Cell-Problem +# -Add Option to change NumberOfElements used for computation of Cell-Problem # --- Define Quantity of interest: @@ -222,7 +223,7 @@ for theta in X_Values: print('angle used') Y_Values.append(type) if yName =='curvature': - print('angle used') + print('curvature used') Y_Values.append(curvature) @@ -230,9 +231,12 @@ print("(Output) Values of " + yName + ": ", Y_Values) idx = find_nearestIdx(Y_Values, 0) -print(' Idx of value closest to 0', idx) +print(' Idx of value closest to 0:', idx) ValueClose = Y_Values[idx] -print('GammaValue(Idx) with mu_gamma closest to q_3^*', ValueClose) +print('GammaValue(Idx) with mu_gamma closest to q_3^*:', ValueClose) +print('Theta(Idx) with curvature closest to 0:', ValueClose) + + @@ -312,14 +316,19 @@ f,ax=plt.subplots(1) # ax.plot(X_Values, Y_Values) -ax.scatter(X_Values, Y_Values) +# ax.scatter(X_Values, Y_Values) # plt.plot(x_plotValues, y_plotValues,'.') # plt.scatter(X_Values, Y_Values, alpha=0.3) # plt.scatter(X_Values, Y_Values) # plt.plot(X_Values, Y_Values,'.') -plt.plot([X_Values[0],X_Values[-1]], [Y_Values[0],Y_Values[-1]]) + + +# plt.plot([X_Values[0],X_Values[-1]], [Y_Values[0],Y_Values[-1]]) # plt.axis([0, 6, 0, 20]) +ax.plot(X_Values[X_Values>jump_xValues[0]], Y_Values[X_Values>jump_xValues[0]] ,'b') +ax.plot(X_Values[X_Values<jump_xValues[0]], Y_Values[X_Values<jump_xValues[0]], 'b') + plt.xlabel(xName) # plt.ylabel(yName) -- GitLab