Commit 9dac652c authored by ogimgio's avatar ogimgio
Browse files

cant find sol, nelder mead iterating infinitely

parent 92d80037
......@@ -54,12 +54,10 @@ def getCostFunction(controlPoints):
n = len(controlPoints)-1
parametrized_intervalI = []
cost_functional = 0
total_lenght = np.load('total_lenght.npy')
#for every clustered element we need to find its corrispondance interval and do the thing.
for k in range(0,len(clusters)):
print(k)
quad = 0
parametrized_points = []
......@@ -88,6 +86,7 @@ points = np.array([[0.05,-3],[0.1,4],[0.2,10],[0.3,3],[0.4,8],[0.5,0],[0.7,-1],[
#print(getCostFunction(points))
cost_func,curve1 = getCostFunction(points)
curve1 = curve1[curve1[:,2].argsort()]
#points[0] = 3*(points[1]-points[0])
#points[1] = 3*(points[2] -points[1])
......@@ -142,15 +141,13 @@ def auglagrange(fun,h,x0 , lambda0 , tol , kmax):
def fun(x):
print(x)
print(x[0])
return 0.6*x[0]**2+0.5*x[1]*x[0]-x[1]+3*x[0]
def h(x):
return x[0]**2 + x[1]**2 -1
xyposition = np.load('XYpositionElements.npy')
levelSetValuesCentral = np.load('levelSetValues.npy')
levelSetValuesCentral = np.load('levelSetValuesCentral.npy')
def hdistance(x):
slin = np.linspace(0,1,1000)
sxs = []
......@@ -162,7 +159,8 @@ def hdistance(x):
#calcualte d(xi)
d = geo.Point(xyposition[i]).distance(geo.LineString(sxs))
#calculate d(xi) - phi(xi)
err += d - levelSetValuesCentral[i][0]
err += d - abs(levelSetValuesCentral[i][0])
print('error',err)
return err #=> the bigger is the most impact it would have
def secondDerivativeBezier(n,s,x):
......@@ -205,13 +203,14 @@ kmax=500
p=1;#number of equality constraints, so lampbda size depends on h size, is the same
lambda0=np.random.rand(p,1).transpose()
x = auglagrange(fun,h,x0, lambda0, tol, kmax)
print(x)
#print(x)
#test call of final fun
#center of each cluster which i considered -> each cluster has its own velocity!
x0 = np.array(np.load('XYpositionElementsX0.npy'))
#reshape x+y form
x0 = x0.reshape(1,len(x0)*2)
print('xo what is',x0)
print(x0.shape)
res = auglagrange(finalFun,hdistance,x0, lambda0, tol, kmax)
......
......@@ -7,7 +7,7 @@ import random
# Define the solution one wants to load
ProblemName = "Advection_Translation_Bump_2Fluids"
SettingsName = "SecondOrderGalerkin"
Resolution = 80
Resolution = 20
# Load the solution (may take time for coarse meshes)
Mesh, _, Problem, Solution, Settings = DL.Load_Solution(ProblemName, SettingsName, Resolution,0)
......@@ -51,7 +51,7 @@ for el in elementsIndexWithInterface:
i=i+1
#approximate curve lenght
print("Full approximative lenght:", curveLenghtApprox)
#np.save('total_lenght',curveLenghtApprox)
np.save('total_lenght',curveLenghtApprox)
#SET EPSILON
......@@ -110,8 +110,8 @@ def clusterElements(elements,n):
print('Elements considered:', elementConsidered)
clustersOfElements = clusterElements(elementConsidered,1)
print(clustersOfElements)
#np.save('clusters', clustersOfElements)
print('cluster',clustersOfElements)
np.save('clusters', clustersOfElements)
#plot the whole
for el in range(0,len(Mesh.InnerElements)):
......@@ -192,8 +192,8 @@ for i in range(0,len(blueGreenElementsCenters)):
blueGreenElementsNearestRed.append(elementsIndexWithInterface[minI])
print(blueGreenElements[25])
print(blueGreenElementsNearestRed[25])
#np.save('blueGreenElements', blueGreenElements)
#np.save('blueGreenElementsNearestRed', blueGreenElementsNearestRed)
np.save('blueGreenElements', blueGreenElements)
np.save('blueGreenElementsNearestRed', blueGreenElementsNearestRed)
#ordered list of red triangles
orderedElementsIndexWithInterface = []
......@@ -223,7 +223,7 @@ for i in range(0,len(elementsIndexWithInterface) + 1):
break
print(orderedElementsIndexWithInterface)
#np.save("orderedRedElements",orderedElementsIndexWithInterface)
np.save("orderedRedElements",orderedElementsIndexWithInterface)
#setting l1 and l2 for red triangles
intervalList = []
singleInterval = []
......@@ -241,7 +241,7 @@ for i in range(0,len(orderedElementsIndexWithInterface)):
intervalList.append(singleInterval)
print(intervalList)
#np.save('intervalList', intervalList)
np.save('intervalList', intervalList)
#what is l1 and l2 of an element in blue Elements?
print('The l1 and l2 of element %s is the same as the one from the red, and the nearest is %s, so: ' %(blueGreenElements[15],blueGreenElementsNearestRed[15]))
l1 = intervalList[orderedElementsIndexWithInterface.index(blueGreenElementsNearestRed[15])][0]
......@@ -314,13 +314,13 @@ for i in range(0,len(XYpositionElements)):
StateValCenter, LSValCenter = DL.Interpolator(Settings, Mesh, Solution, globalIndexesFullNoDuplicates[i], np.array([XYpositionElements[i]]))
levelSetValuesToExport.append(LSValCenter[0])
np.save('levelSetValues',levelSetValuesToExport)
np.save('levelSetValuesCentral',levelSetValuesToExport)
#'Each global vertices given by globalIndexesFullNoDuplicates
np.save('levelSetValuesALL',levelSetValuesFull)
#i have to delete duplicates because the levelset values are not twice. whats wrong?
"""levelSetValuesFullNoDuplicates = list(dict.fromkeys(np.array(levelSetValuesFull).flatten()))
levelSetValuesFullNoDuplicates = list(dict.fromkeys(np.array(levelSetValuesFull).flatten()))
np.save('levelSetValues',levelSetValuesFullNoDuplicates)
#level set values full: levelSetValuesFullNoDuplicates
......@@ -396,6 +396,4 @@ for i in range(0,len(resultingAlphasLSVPerClusters)):
print('%s Weight(=alpha) is: %s' %(resultingAlphasLSVPerClusters[i][0],alphas_weight_fun(minLSV,maxLSV,abs(resultingAlphasLSVPerClusters[i][0]),0)))
weightedAlphas.append(alphas_weight_fun(minLSV,maxLSV,abs(resultingAlphasLSVPerClusters[i][0]),0))
np.save('weightedAlphas', np.array(weightedAlphas))
"""
\ No newline at end of file
np.save('weightedAlphas', np.array(weightedAlphas))
\ No newline at end of file
No preview for this file type
No preview for this file type
No preview for this file type
bezier_final.png

24.9 KB | W: | H:

bezier_final.png

24.8 KB | W: | H:

bezier_final.png
bezier_final.png
bezier_final.png
bezier_final.png
  • 2-up
  • Swipe
  • Onion skin
No preview for this file type
No preview for this file type
No preview for this file type
No preview for this file type
No preview for this file type
No preview for this file type
No preview for this file type
No preview for this file type
No preview for this file type
No preview for this file type
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment