Skip to content

P
Python

m0_59476899 / Python
与 Fork 源项目一致

Fork自 inscode / Python

通知 1
Star 0
Fork 0
P
Python
提交 20a22b00 编写于 作者: 6 62ac1e676314f46ab57773b8
浏览文件
操作

Auto commit

上级 df352508
展开全部 隐藏空白更改
内联 并排
Showing with 1847 addition and 1 deletion
.inscode
浏览文件 @ 20a22b00
run = "pip install -r requirements.txt;python3 main.py"
run = "pip install -r requirements.txt;cd optimization;python3 manage.py runserver 0.0.0.0:8080 --noreload"
[env]
VIRTUAL_ENV = "/root/${PROJECT_DIR}/venv"
......
optimization/.idea/encodings.xml 0 → 100644
浏览文件 @ 20a22b00
<?xml version="1.0" encoding="UTF-8"?>
<project version="4">
<component name="Encoding" addBOMForNewFiles="with NO BOM" />
</project>
\ No newline at end of file
optimization/.idea/misc.xml 0 → 100644
浏览文件 @ 20a22b00
<?xml version="1.0" encoding="UTF-8"?>
<project version="4">
<component name="JavaScriptSettings">
<option name="languageLevel" value="ES6" />
</component>
<component name="ProjectRootManager" version="2" project-jdk-name="Python 3.7" project-jdk-type="Python SDK" />
</project>
\ No newline at end of file
optimization/.idea/modules.xml 0 → 100644
浏览文件 @ 20a22b00
<?xml version="1.0" encoding="UTF-8"?>
<project version="4">
<component name="ProjectModuleManager">
<modules>
<module fileurl="file://$PROJECT_DIR$/.idea/optimization.iml" filepath="$PROJECT_DIR$/.idea/optimization.iml" />
</modules>
</component>
</project>
\ No newline at end of file
optimization/.idea/optimization.iml 0 → 100644
浏览文件 @ 20a22b00
<?xml version="1.0" encoding="UTF-8"?>
<module type="PYTHON_MODULE" version="4">
<component name="FacetManager">
<facet type="django" name="Django">
<configuration>
<option name="rootFolder" value="$MODULE_DIR$" />
<option name="settingsModule" value="optimization/settings.py" />
<option name="manageScript" value="$MODULE_DIR$/manage.py" />
<option name="environment" value="&lt;map/&gt;" />
<option name="doNotUseTestRunner" value="false" />
<option name="trackFilePattern" value="migrations" />
</configuration>
</facet>
</component>
<component name="NewModuleRootManager">
<content url="file://$MODULE_DIR$" />
<orderEntry type="inheritedJdk" />
<orderEntry type="sourceFolder" forTests="false" />
</component>
<component name="TemplatesService">
<option name="TEMPLATE_CONFIGURATION" value="Django" />
<option name="TEMPLATE_FOLDERS">
<list>
<option value="$MODULE_DIR$/../optimization\templates" />
</list>
</option>
</component>
<component name="TestRunnerService">
<option name="PROJECT_TEST_RUNNER" value="Unittests" />
</component>
</module>
\ No newline at end of file
optimization/.idea/workspace.xml 0 → 100644
浏览文件 @ 20a22b00
此差异已折叠。
optimization/caculation/admin.py 0 → 100644
浏览文件 @ 20a22b00
from django.contrib import admin
# Register your models here.
optimization/caculation/algorithm/KKT.py 0 → 100644
浏览文件 @ 20a22b00
from sympy import *
def main():
#构造拉格朗日等式
if a == "min":
L = -f + l1 * f1 + l2 * f2
g1 = diff(L, x1)
g2 = diff(L, x2)
print(
f"KKT条件：\n{g1}=0 ； {g2}=0 ; \n{l1 * f1}=0 ; {l2 * f2}=0 ; {f1}>=0 ; {f2}>=0 ; \n x1>=0 ; x2>=0 ; l1>=0 ; l2>=0 ; l3>=0")
print("__________开始计算_______________")
y = 10 ** 50
#λ1≠0；λ2≠0；λ3＝0
x110 = solve([f1, f2], [x1, x2])
y110 = f.subs(x110)
if (f1.subs(x110) >= 0) & (f2.subs(x110) >= 0) & (f3.subs(x110) >= 0):
print(f"λ1≠0；λ2≠0；λ3＝0时：解得{x110}，此时y={y110},此解符合约束条件")
if y > y110:
x = x110
y = y110
else:
print(f"λ1≠0；λ2≠0；λ3＝0时：解得{x110}，此时y={y110},此解不符合约束条件")
#λ1≠0；λ2＝0；λ3≠0
x101 = solve([f1, f3], [x1, x2])
y101 = f.subs(x101)
if (f1.subs(x101) >= 0) & (f2.subs(x101) >= 0) & (f3.subs(x101) >= 0):
print(f"λ1≠0；λ2＝0；λ3≠0时：解得{x101}，此时y={y101},此解符合约束条件")
if y > y101:
x = x101
y = y101
else:
print(f"λ1≠0；λ2＝0；λ3≠0时：解得{x110}，此时y={y101},此解不符合约束条件")
#λ1≠0；λ2＝0；λ3＝0
x100 = solve([f1, g1.subs([(l2, 0), (l3, 0)]), g2.subs([(l2, 0), (l3, 0)])], [l1, x1, x2])
y100 = f.subs(x100)
if (f1.subs(x100) >= 0) & (f2.subs(x100) >= 0) & (f3.subs(x100) >= 0) :
print(f"λ1≠0；λ2＝0；λ3＝0时：解得{x100}，此时y={y100},此解符合约束条件")
if y > y100:
x = x100
y = y100
else:
print(f"λ1≠0；λ2＝0；λ3＝0时：解得{x100}，此时y={y100},此解不符合约束条件")
#λ1＝0；λ2≠0；λ3≠0
x011 = solve([f2, f3], [x1, x2])
y011 = f.subs(x011)
if (f1.subs(x011) >= 0) & (f2.subs(x011) >= 0) & (f3.subs(x011) >= 0):
print(f"λ1＝0；λ2≠0；λ3≠0时：解得{x011}，此时y={y011},此解符合约束条件")
if y > y011:
x = x011
y = y011
else:
print(f"λ1＝0；λ2≠0；λ3≠0时：解得{x011}，此时y={y011},此解不符合约束条件")
#λ1＝0；λ2＝0；λ3≠0
x001 = solve([f3, g1.subs([(l1, 0), (l2, 0)]), g2.subs([(l1, 0), (l2, 0)])], [l3, x1, x2])
y001 = f.subs(x001)
if (f1.subs(x001) >= 0) & (f2.subs(x001) >= 0) & (f3.subs(x001) >= 0):
print(f"λ1＝0；λ2＝0；λ3≠0时：解得{x001}，此时y={y001},此解符合约束条件")
if y > y001:
x = x001
y = y001
else:
print(f"λ1＝0；λ2＝0；λ3≠0时：解得{x001}，此时y={y001},此解不符合约束条件")
#λ1＝0；λ2≠0；λ3＝0
x010 = solve([f2,g1.subs([(l1, 0), (l3, 0)]), g2.subs([(l1, 0), (l3, 0)])], [l2, x1, x2])
y010 = f.subs(x010)
if (f1.subs(x010)>=0) & (f2.subs(x010)>=0) & (f3.subs(x010)>=0):
print(f"λ1＝0；λ2≠0；λ3＝0时：解得{x010}，此时y={y010},此解符合约束条件")
if y > y010:
x = x010
y = y010
else:
print(f"λ1＝0；λ2≠0；λ3＝0时：解得{x010}，此时y={y010},此解不符合约束条件")
elif a == "max":
L = f + l1 * f1 + l2 * f2 + l3 * f3
g1 = diff(L, x1)
g2 = diff(L, x2)
print(
f"KKT条件：\n{g1}=0 ； {g2}=0 ; \n{l1 * f1}=0 ; {l2 * f2}=0 ; {l2 * f3}=0 ; \n{f1}>=0 ; {f2}>=0 ; {f3}>=0 ;\n x1>=0 ; x2>=0 ; l1>=0 ; l2>=0 ; l3>=0")
print("__________开始计算_______________")
y = -10 ** 50
#λ1≠0；λ2≠0；λ3＝0
x110 = solve([f1, f2], [x1, x2])
y110 = f.subs(x110)
if (f1.subs(x110) >= 0) & (f2.subs(x110) >= 0) & (f3.subs(x110) >= 0):
print(f"λ1≠0；λ2≠0；λ3＝0时：解得{x110}，此时y={y110},此解符合约束条件")
if y < y110:
x = x110
y = y110
else:
print(f"λ1≠0；λ2≠0；λ3＝0时：解得{x110}，此时y={y110},此解不符合约束条件")
#λ1≠0；λ2＝0；λ3≠0
x101 = solve([f1, f3], [x1, x2])
y101 = f.subs(x101)
if (f1.subs(x101) >= 0) & (f2.subs(x101) >= 0) & (f3.subs(x101) >= 0):
print(f"λ1≠0；λ2＝0；λ3≠0时：解得{x101}，此时y={y101},此解符合约束条件")
if y < y101:
x = x101
y = y101
else:
print(f"λ1≠0；λ2＝0；λ3≠0时：解得{x110}，此时y={y101},此解不符合约束条件")
#λ1≠0；λ2＝0；λ3＝0
x100 = solve([f1, g1.subs([(l2, 0), (l3, 0)]), g2.subs([(l2, 0), (l3, 0)])], [l1, x1, x2])
y100 = f.subs(x100)
if (f1.subs(x100) >= 0) & (f2.subs(x100) >= 0) & (f3.subs(x100) >= 0) :
print(f"λ1≠0；λ2＝0；λ3＝0时：解得{x100}，此时y={y100},此解符合约束条件")
if y < y100:
x = x100
y = y100
else:
print(f"λ1≠0；λ2＝0；λ3＝0时：解得{x100}，此时y={y100},此解不符合约束条件")
#λ1＝0；λ2≠0；λ3≠0
x011 = solve([f2, f3], [x1, x2])
y011 = f.subs(x011)
if (f1.subs(x011) >= 0) & (f2.subs(x011) >= 0) & (f3.subs(x011) >= 0):
print(f"λ1＝0；λ2≠0；λ3≠0时：解得{x011}，此时y={y011},此解符合约束条件")
if y < y011:
x = x011
y = y011
else:
print(f"λ1＝0；λ2≠0；λ3≠0时：解得{x011}，此时y={y011},此解不符合约束条件")
#λ1＝0；λ2＝0；λ3≠0
x001 = solve([f3, g1.subs([(l1, 0), (l2, 0)]), g2.subs([(l1, 0), (l2, 0)])], [l3, x1, x2])
y001 = f.subs(x001)
if (f1.subs(x001) >= 0) & (f2.subs(x001) >= 0) & (f3.subs(x001) >= 0):
print(f"λ1＝0；λ2＝0；λ3≠0时：解得{x001}，此时y={y001},此解符合约束条件")
if y < y001:
x = x001
y = y001
else:
print(f"λ1＝0；λ2＝0；λ3≠0时：解得{x001}，此时y={y001},此解不符合约束条件")
#λ1＝0；λ2≠0；λ3＝0
x010 = solve([f2,g1.subs([(l1, 0), (l3, 0)]), g2.subs([(l1, 0), (l3, 0)])], [l2, x1, x2])
y010 = f.subs(x010)
if (f1.subs(x010)>=0) & (f2.subs(x010)>=0) & (f3.subs(x010)>=0):
print(f"λ1＝0；λ2≠0；λ3＝0时：解得{x010}，此时y={y010},此解符合约束条件")
if y < y010:
x = x010
y = y010
else:
print(f"λ1＝0；λ2≠0；λ3＝0时：解得{x010}，此时y={y010},此解不符合约束条件")
print(f"___________计算结果________________")
print(f"x点位：{x}；此时y:{y}")
if __name__ == '__main__':
x1, x2, l1, l2, l3 = symbols('x1,x2,l1,l2,l3')
#f直接写进来，不用转化
f = x1**2-2*x1+2*x2**2
a="min"
#约束条件应该是f*>=0
f1 = x1-x2**2
f2 = x2-x1
main()
optimization/caculation/algorithm/OLS.py 0 → 100644
浏览文件 @ 20a22b00
from sympy import *
import numpy as np
def OLS(A,b):
output=''
x1, x2= symbols('x1, x2')
A=np.array(A)
b=np.array(b)
F=0*x1
for i in range(np.shape(A)[0]):
F=F+(A[i][0]*x1+A[i][1]*x2-b[i])**2
output+=f'令F(x)={F},\n则原问题化为minF(x)(x∈R)\n记A={A}，x=[x1,x2]^T,b={b},\n则F(x)=(Ax-b)^T(Ax-b)。\n'
output+=f'{A.T}'
ATA=np.dot(A.T,A)
ATb=np.dot(A.T,b)
x=np.dot(np.linalg.inv(ATA),ATb)
output+=f'∴A^TA:\n{ATA}\n A^Tb:\n{ATb}\n (A^TA)^-1:\n{np.linalg.inv(ATA)}\n∴近似解x=(A^TA)^-1 * A^Tb={x}^T'
return output
if __name__ == '__main__':
#等式左边的系数
A=[[1,2],
[3,-1],
[2,-3]]
#等式右边的常数
b=[5,
3,
0]
output=OLS(A,b)
print(output)
optimization/caculation/algorithm/requirements.txt 0 → 100644
浏览文件 @ 20a22b00
Flask==2.2.3
numpy==1.24.3
sympy==1.11.1
optimization/caculation/algorithm/web.py 0 → 100644
浏览文件 @ 20a22b00
from flask import Flask
app = Flask(__name__)
@app.route("/show/info")
def index():
return "中"
if __name__=='__main__':
app.run()
\ No newline at end of file
optimization/caculation/algorithm/一维度/二次逼近法.py 0 → 100644
浏览文件 @ 20a22b00
from sympy import *
import numpy as np
def main():
global init,step
n=0
print(f'1.给定初始区间{init}，给定精度{dx}，设f(x1)>f(x2),f(x2)<f(x3)，令k=1，x_{n}=x2')
x1=init[0]
x3=init[1]
x2 = (x1+x3)/2
x_k_1 = x2
k=1
a, b, c = symbols('a,b,c')
F = a + b * x + c * x ** 2
while 1:
n=n+1
print(f'\n_______________第{n}次迭代___________________________________________')
abc = solve([F.subs(x,x1) - f.subs(x,x1), F.subs(x,x2) - f.subs(x,x2), F.subs(x,x3) - f.subs(x,x3)], [a, b, c])
x_k = solve(diff(F, x).subs(abc), x)[0]
print(f'2.设φ(x)=a+bx+cx^2,令φ(xi)=f(xi) (i=1,2,3)\n'
f' 解得{abc}，则φ(x)={abc[a]}+{abc[b]}x+{abc[c]}x^2\n'
f' 解得φ(x)极小点x_{n}=-b/2c={x_k}')
x_all = np.array([x1, x2, x3, x_k])
y_all = np.array([f.subs(x, x1), f.subs(x, x2), f.subs(x, x3), f.subs(x, x_k)])
if abs(f.subs(x,x_k) - f.subs(x,x_k_1)) >= dx:
print(f' {abs(f.subs(x,x_k) - f.subs(x,x_k_1))}=|f(x_k)-f(x_k-1)| >= dx,转3')
x_k_1=x_k
x2=x_all[np.argsort(y_all)[0]]
x1=x_all[np.argsort(y_all)[1]]
x3=x_all[np.argsort(y_all)[2]]
print(f'3.计算x1,x2,x3,x_{n}函数值:\n'
f' f(x1)={f.subs(x,x1)}\n'
f' f(x2)={f.subs(x,x2)}\n'
f' f(x3)={f.subs(x,x3)}\n'
f' f(x1)={f.subs(x,x_k)}\n'
f' 其中x1,x2,x3,x_{n}中函数值最小的点设为x2={x_all[np.argsort(y_all)[0]]},左右两点设为x1,x3={x_all[np.argsort(y_all)[1]],x_all[np.argsort(y_all)[2]]}')
else:
print(f' {abs(f.subs(x,x_k) - f.subs(x,x_k_1))}=|f(x_k)-f(x_k-1)| < dx,结束算法')
print(f'\n++++++++++++++++++++++++++++++结果+++++++++++++++++++++++++++++++++++++++\n'
f'极值点在{x2}附近')
break
if __name__ == '__main__':
x ,x1,x2,x3= symbols('x,x1,x2,x3')
f = x**2+x # 设置初始函数
init = [0,4] # 初始区间
dx = 0.01 # 迭代精度
m = "min"
main()
optimization/caculation/algorithm/一维度/进退法.py 0 → 100644
浏览文件 @ 20a22b00
from sympy import *
def main():
global init,step
x1=init
k=0
print(f'1.给定初始点{init},初始步长step={step}>0,计算f(x1)={f.subs(x,x1)},并令k=0')
n=0
while 1:
n=n+1
print(f'_______________第{n}次迭代___________________________________')
x3 = x1 + step
k = k + 1
print(f'2.令x3=x1+step={x3}，计算f(x3)，令k=k+1={k}')
if f.subs(x,x1) > f.subs(x,x3):
print(f'3.此时x1,x2,x3={x1, x2, x3},f(x3)={f.subs(x, x3)} < {f.subs(x, x1)}=f(x1)，转4')
x2 = x1
x1 = x3
step = 2 * step
print(f'4.令x2=x1,x1=x3,step=2*step={step},转2')
else:
if k == 1:
step = -step
x2 = x3
print(f'3.此时x1,x2,x3={x1, x2, x3},f(x3)={f.subs(x, x3)} >= {f.subs(x, x1)}=f(x1)，转5')
print(f'5.若k=1,令x2=x3,step=-step={step},转2')
else:
print(f'3.此时x1,x2,x3={x1, x2, x3},f(x3)={f.subs(x,x3)} > {f.subs(x,x1)}=f(x1)，转5\n'
f'4.若k≠1，停止计算')
break
print(f"\n+++++++++++++++++++++结果++++++++++++++++++++++++++++++")
if x2<x3:
print(f"区间[{x2},{x3}]即为包含极小点的区间")
else:
print(f"区间[{x3},{x2}]即为包含极小点的区间")
if __name__ == '__main__':
x ,x1,x2,x3= symbols('x,x1,x2,x3')
f = (x - 1) ** 2 # 设置初始函数
init = 30 # 初始值
step = 1 # 初始步长
m = "min"
main()
optimization/caculation/algorithm/一维度/黄金分割法.py 0 → 100644
浏览文件 @ 20a22b00
from sympy import *
def main(a,b):
if m=="max":
i = 0
while True:
i += 1
x1 = b - 0.618 * (b - a)
x2 = a + 0.618 * (b - a)
if f.subs(x,x1) > f.subs(x,x2):
b = x2
print(f'____________________第{i}次迭代_______________________________________\n'
f'λ{i} = {a}+0.382*({b}-{a})={x1}\n'
f'μ{i} = {a}+0.618*({b}-{a})={x2}\n'
f'f(λ1)={x1}^2={f.subs(x,x1)} f(μ1)={x2}^2={f.subs(x,x2)}\n'
f'f(λ1)>f(μ1) 新的搜索区间[{a},{b}]')
elif f.subs(x,x1) <= f.subs(x,x2):
a = x1
print(f'____________________第{i}次迭代_______________________________________\n'
f'λ{i} = {a}+0.382*({b}-{a})={x1}\n'
f'μ{i} = {a}+0.618*({b}-{a})={x2}\n'
f'f(λ1)={x1}^2={f.subs(x,x1)} f(μ1)={x2}^2={f.subs(x,x2)}\n'
f'f(λ1)<f(μ1) 新的搜索区间[{a},{b}]')
DX = abs(b - a)
if DX <= dx:
print("+++++++++++++++++++++结果+++++++++++++++++++++++++++++++++++++++++++++++++++")
print(f'达到精度需迭代{i}次')
print(f'最终的搜索区间为：[{min(a, b)},{max(a, b)}]')
print(f'确定最大值的两端函数值为：{f.subs(x,a)}, {f.subs(x,b)}，函数最大值：{f.subs(x,(a + b) / 2)}')
break
else:
pass
if m=="min":
i = 0
while True:
i += 1
x1 = a + 0.382 * (b - a)
x2 = a + 0.618 * (b - a)
if f.subs(x,x1) > f.subs(x,x2):
a = x1
print(f'____________________第{i}次迭代_______________________________________\n'
f'λ{i} = {a}+0.382*({b}-{a})={x1}\n'
f'μ{i} = {a}+0.618*({b}-{a})={x2}\n'
f'f(λ1)={x1}^2={f.subs(x,x1)} f(μ1)={x2}^2={f.subs(x,x2)}\n'
f'f(λ1)>f(μ1) 新的搜索区间[{a},{b}]')
elif f.subs(x,x1) <= f.subs(x,x2):
b = x2
print(f'____________________第{i}次迭代_______________________________________\n'
f'λ{i} = {a}+0.382*({b}-{a})={x1}\n'
f'μ{i} = {a}+0.618*({b}-{a})={x2}\n'
f'f(λ1)={x1}^2={f.subs(x,x1)} f(μ1)={x2}^2={f.subs(x,x2)}\n'
f'f(λ1)<f(μ1) 新的搜索区间[{a},{b}]')
DX = abs(b - a)
if DX <= dx:
print("+++++++++++++++++++++结果+++++++++++++++++++++++++++++++++++++++++++++++++++")
print(f'达到精度需迭代{i}次')
print(f'最终的搜索区间为：[{min(a, b)},{max(a, b)}]')
print(f'确定最大值的两端函数值为：{f.subs(x,a)}, {f.subs(x,b)}，函数最大值：{f.subs(x,(a + b) / 2)}')
print(f'第k次迭代后区间长度为{b-a}(初始区间长度的((sqrt(5)-1)/2)^k倍)')
break
else:
pass
if __name__ == '__main__':
x = symbols('x')
f = x**2-4*x+2# 设置初始函数
a = -1 # 区间下限
b = 6 # 区间上限
dx = 0.01 # 迭代精度
m = "min"
main(a,b)
optimization/caculation/algorithm/可行方向法/Zoutendijk.py 0 → 100644
浏览文件 @ 20a22b00
import numpy as np
from sympy import *
def graphic(fd,Fd):
a=0
bujiru=0
for i in range(len(Fd)):
for j in range(i+1,len(Fd),1):
if len(solve([Fd[i],Fd[j]],[d1,d2]))==0:
continue
else:
d_new = np.array([[float(solve([Fd[i], Fd[j]], [d1, d2])[d1]), float(solve([Fd[i], Fd[j]], [d1, d2])[d2])]])
for Fdi in Fd:
if Fdi.subs([(d1,d_new[0][0]),(d2,d_new[0][1])])<0:
bujiru=1
if bujiru==0:
if a == 0:
d = d_new
else:
d = np.concatenate((d, d_new), axis=0)
a = a + 1
bujiru=0
y=[[]for i in range(len(d))]
for m in range(len(d)):
y[m]=fd.subs([(d1,d[m][0]),(d2,d[m][1])])
d_min=d[np.argmin(np.array(y))]
return d_min
def search(di,step_max):
new_x1 = init[0] + t * di[0]
new_x2 = init[1] + t * di[1]
f_new = f.subs([(x1, new_x1), (x2, new_x2)]) # 新函数
grad_new = diff(f_new, t) # 对t偏导
t_val = nsolve(grad_new, t, 0)
if t_val>step_max:
print(f'使用一维搜索计算λ={t_val}>λmax\n'
f'最终取λ={step_max}')
return step_max
else:
print(f'使用一维搜索计算λ={t_val}<λmax\n'
f'最终取λ={t_val}')
return t_val
def main():
global a,init,d
g1=diff(f,x1)
g2=diff(f,x2)
a=np.array(a)
x_num=np.shape(a)[1]-1
init=np.array(init)
print(f'目标函数的梯度∇f(x)=[{g1}\n'
f' {g2}]')
num=0
while 1:
num=num+1
print(f"___________________第{num}次迭代__________________")
g1_new=g1.subs(x1,init[0])
g2_new=g2.subs(x2,init[1])
print(f'1.在点x^({num})={init},∇f(x)=[{g1_new},{g2_new}]')
A1_num=0
A2_num=0
for i in range(len(a)):
if sum(a[i][:-1]*init)+a[i][-1]==0:
A1_num=A1_num+1
else:
A2_num=A2_num+1
A1=np.ones([A1_num,x_num])
A2=np.ones([A2_num,x_num])
b1=np.ones(A1_num)
b2=np.ones(A2_num)
A1_num=0
A2_num=0
for i in range(len(a)):
if sum(a[i][:-1]*init)+a[i][-1]==0:
A1[A1_num]=a[i][:-1]
b1[A1_num]=a[i][-1]
A1_num=A1_num+1
else:
A2[A2_num]=a[i][:-1]
b2[A2_num]=a[i][-1]
A2_num=A2_num+1
b1=-b1
b2=-b2
print(f'起作用约束(紧约束)系数矩阵,约束右端\nb1={b1}\nA1={A1}\n'
f'不起作用约束(紧约束)系数矩阵，约束右端\nb2={b2}\nA2={A2}')
Fd=[[]for i in range(A1_num+4)]
fdi=0*d1
for n in range(A1_num):
for m in range(x_num):
fdi=fdi+A1[n][m]*d[m]
Fd[n]=fdi
fdi=0*d1
Fd[n+1]=d1+1
Fd[n+2]=-d1+1
Fd[n+3]=d2+1
Fd[n+4]=-d2+1
fd=g1_new*d[0]+g2_new*d[1]
di=graphic(fd,Fd)
print(f'2.求在点x^({num})=init处下降可行方向d\n'
f'即解min∇f(x^{num})^Td={fd}\n'
f's.t. A1*d>=0\n'
f' |dj|<=1(j=1,2))\n'
f'通过图解法求得d^{num}={di}')
dd=np.dot(np.array(A2),np.array(di).T)
bb=np.array(b2)-np.dot(np.array(A2),np.array(init).T)
for i in range(len(dd)):
if (dd[i]<0) & ((dd[i]/bb[i])<=(dd[0]/bb[0])):
step_max = dd[i]/bb[i]
print(f'3.求步长：\n'
f'∵d=A2*d{num}={dd},b=b2-A2x{num}={bb}\n'
f'∴λmax=min(d/b|d<0)={step_max}')
t = search(di,step_max)
init_new = init+t*di
print(f'目标点位{init_new}')
if (init_new==init).all():
print(f'\n_________________________结果______________________\n'
f'最终点位{init_new},最终函数值{f.subs([(x1,init_new[0]),(x2,init_new[1])])}')
break
else:
init=init_new
if __name__ == '__main__':
x1, x2 ,d1,d2,t= symbols('x1,x2,d1,d2,t')
f=x1**2+x2**2-2*x1-4*x2+6
d=[d1,d2]
a=[[-2,1,1],[-1,-1,2],[1,0,0],[0,1,0]]
init=[0,0]
main()
optimization/caculation/algorithm/可行方向法/图解法.py 0 → 100644
浏览文件 @ 20a22b00
import numpy as np
from sympy import *
def graphic(F):
for i in range(len(F)):
for j in range(i+1,len(F),1):
x_new=np.array([[float(solve([F[i],F[j]],[x1,x2])[x1]),float(solve([F[i],F[j]],[x1,x2])[x2])]])
if (i==0)&(j==1):
x=x_new
else:
x=np.concatenate((x,x_new),axis=0)
y=[[]for i in range(len(x))]
for m in range(len(x)):
y[m]=f.subs([(x1,x[m][0]),(x2,x[m][1])])
x_min=x[np.argmin(np.array(y))]
return x_min
if __name__=="__main__":
x1, x2 = symbols('x1,x2')
f=5*x1-6*x2
f1=10-x1-2*x2
f2=5-2*x1+x2
f3=4-x1+4*x2
f4=x1
f5=x2
F=[f1,f2,f3,f4,f5]
x_min=graphic(F)
\ No newline at end of file
optimization/caculation/algorithm/各种解.py 0 → 100644
浏览文件 @ 20a22b00
# -*- coding: utf-8 -*-
# @File : 外点法.py
# @Time : 2021/4/14
# @Author : jiayang xie
from sympy import *
def get_dic(x_2):
if isinstance(x_2,list):
a={}
a[x1]=x_2[0][0]
a[x2]=x_2[0][1]
return a
else:
return x_2
def main():
for xi in x:
F_i=F+[xi]
x_val=solve(F_i)
print(x_val)
print(f"基本解有{len(x)}个，分别是{x}")
if __name__ == '__main__':
x1, x2 ,x3,x4 = symbols('x1,x2,x3,x4')
f1 = -3*x1+x2+x3-4
f2 = x1-x2-2
F = [f1,f2]
x=[x1,x2,x3]
main()
\ No newline at end of file
optimization/caculation/algorithm/惩罚函数/内点法.py 0 → 100644
浏览文件 @ 20a22b00
# -*- coding: utf-8 -*-
# @File : 外点法.py
# @Time : 2021/4/14
# @Author : jiayang xie
from sympy import *
def get_dic(x_2):
if isinstance(x_2,list):
a={}
a[x1]=x_2[0][0]
a[x2]=x_2[0][1]
return a
else:
return x_2
def main():
#用-ln(f1)求
g = f -r*ln(f1)
x = solve([diff(g,x1),diff(g,x2)],[x1,x2])
print(f"F(x,r)=f(x)+r*P(x)={f}+r*P(x)\n"
f"P(x)=-r*ln({f1})\n"
f"F(x,r)={f}-r*ln({f1})\n"
f"________解析法求驻点_________________________________________")
print(f"0=dg(x,r)/dx1 = {diff(g,x1)}")
print(f"0=dF(x,r)/dx2 = {diff(g,x2)}")
print(f'______________结果______________________________')
print(x)
if len(x)==0:
print(f"当{f1}<0时，无解")
else:
x=get_dic(x)
if f1.subs(x).subs(r,10**-50) > 0:
print(f"r->0+，此时解得{x},符合条件\n")
print('+++++++++++++++++++++++++++最好别用第二种±sqrt(r)代码写不出来++++++++++++++++++++++++++++++++++++++++++++++++\n')
#用r/f1求
g = f + r/f1
x = solve([diff(g,x1),diff(g,x2)],[x1,x2])
print(f"F(x,r)=f(x)+r*P(x)={f}+r*P(x)\n"
f"P(x)=-r*ln({f1})\n"
f"F(x,r)={f}+r/({f1}) \n"
f"________解析法求驻点_________________________________________")
print(f"0=dg(x,r)/dx1 = {diff(g,x1)}")
print(f"0=dF(x,r)/dx2 = {diff(g,x2)}")
print(f'______________结果______________________________')
if len(x)==0:
print(f"当{f1}<0时，无解")
else:
x=get_dic(x)
print(f"r->0+，此时解得{x},符合条件\n")
if __name__ == '__main__':
x1, x2 ,r = symbols('x1'), symbols('x2'), symbols('r')
f = x1**2+x2**2-4*x2 # 设置初始函数
##约束条件要都>=0
f1 = -x1**2-x2+1
main()
\ No newline at end of file
optimization/caculation/algorithm/惩罚函数/外点法.py 0 → 100644
浏览文件 @ 20a22b00
# -*- coding: utf-8 -*-
# @File : 外点法.py
# @Time : 2021/4/14
# @Author : jiayang xie
from sympy import *
import numpy as np
def get_dic(x_2):
if isinstance(x_2,list):
a={}
a[x1]=x_2[0][0]
a[x2]=x_2[0][1]
return a
else:
return x_2
def main():
#f1>=0时
F1 = f
x_1 = solve([diff(F1,x1),diff(F1,x2)],[x1,x2])
#f1<0时
F2 = f + Mk * f1 * f1
x_2 = solve([diff(F2,x1),diff(F2,x2)],[x1,x2])
print(f"F(x,Mk)=f(x)+Mk*P(x)={f}+Mk*P(x)\n"
f"P(x)=(min({f1},0))^2)\n"
f"F(x,Mk)={f}+Mk*(min({f1},0))^2)\n"
f"F(x,Mk)={f}+Mk*({f1})^2) ({f1}>=0)\n"
f" {f} ({f1}<0)\n"
f"________解析法求驻点_________________________________________")
print(f"0=dF(x,Mk)/dx1={diff(F1,x1)} ({f1}>=0)\n"
f" ={diff(F2,x1)} ({f1}<0)\n")
print(f"0=dF(x,Mk)/dx2={diff(F1,x2)} ({f1}>=0)\n"
f" ={diff(F2,x2)} ({f1}<0)\n")
print(f'______________求解______________________________')
if len(x_1)==0:
print(f"当{f1}>=0时，无解")
else:
if f1.subs(x_1) >= 0:
print(f"当{f1}>=0时，此时解得{x_1},符合条件")
else:
print(f"当{f1}>=0时，此时解得{x_1},矛盾")
if len(x_2)==0:
print(f"当{f1}<0时，无解")
else:
x_2=get_dic(x_2)
if f1.subs(x_2).subs(Mk,10**50) < 0:
print(f"当{f1}<0时，Mk->正无穷，此时解得{x_2},符合条件")
else:
print(f"当{f1}<0时，Mk->正无穷，此时解得{x_2},矛盾")
if __name__ == '__main__':
x1, x2 ,Mk = symbols('x1'), symbols('x2'), symbols('Mk')
f = x1**4-2*x2+2*x2**2 # 设置初始函数
##约束条件要都>=0
f1 = 2*x1**2+5-x1-x2**2
main()
\ No newline at end of file
optimization/caculation/algorithm/梯度/共轭梯度法(FR).py 0 → 100644
浏览文件 @ 20a22b00
# -*- coding: utf-8 -*-
# @File : 最速下降法.py
# @Time : 2021/4/14
# @Author : laipinyan
from sympy import *
import numpy as np
def main():
global init
for i in range(nub_iter):
print(f'__________________第{i+1}次迭代_________________________________')
#算共轭梯度
g1_new = diff(f, x1).subs([(x1, init[0]),(x2,init[1])]) # 对x偏导
g2_new = diff(f, x2).subs([(x1, init[0]),(x2,init[1])]) # 对y偏导
g_new = np.array([float(g1_new),float(g2_new)]).T
if i==0:
d_new=-g_new
print(f'∵∇f(x)=({diff(f, x1)},{diff(f, x2)})^T\n'
f'∴d^{i + 1}=-∇f(x)=({d_new[0]},{d_new[1]})^T\n'
f'∴x^{i + 1}+λd^{i + 1}=({x1}-{d_new[0]}λ,{x2}-{d_new[1]}λ)^T')
d=d_new
g=g_new
else:
d_new=-g_new+(np.linalg.norm(g_new))**2/(np.linalg.norm(g))**2*d
print(g_new)
print((np.linalg.norm(g_new))**2)
print((np.linalg.norm(g)) ** 2)
print(f'∵∇f(x)=({diff(f, x1)},{diff(f, x2)})^T\n'
f'∴β{i}=||g{i+1}||^2/||g{i}||^2={(np.linalg.norm(g_new))**2/(np.linalg.norm(g))**2}\n'
f'∴d^{i + 1}=-g{i+1}+β{i}d^{i}=({d_new})^T\n'
f'∴x^{i + 1}+λd^{i + 1}=({x1}-{float(d_new[0])}λ,{x2}-{float(d_new[1])}λ)^T')
d=d_new
g=g_new
new_x1 = init[0] + t*(d_new[0])
new_x2 = init[1] + t*(d_new[1])
f_new = f.subs([(x1, new_x1), (x2, new_x2)]) # 新函数
grad_new = diff(f_new, t) # 对t偏导
t_val = nsolve(grad_new, t,0) # 求t
init_new = (new_x1.subs(t, t_val), new_x2.subs(t, t_val)) # 新点
f_value = f.subs([(x1, init_new[0]), (x2, init_new[1])]) # 目标函数值
print(f'令φ(λ)=f(x^{i+1}+λd^{i+1})={f_new}\n'
f'求解minφ(λ)\n'
f'令φ´(λ)={grad_new}=0\n'
f'解得λ{i+1}={t_val}\n'
f'∴x^{i+2}=x^{i+1}+λd^{i+1}={init_new}\n'
f' f(x^{i+2})={f_value}')
print(f"{i+1}次迭代总结：原始点位{init}，共轭方向{d},最优步长{t_val},目标点位{init_new},目标函数值{f_value}")
init=init_new
if __name__ == '__main__':
x1, x2, t = symbols('x1'), symbols('x2'), symbols('t')
f = x1**2-4*x1+2*x2**2 # 设置初始函数
init = (1, -1) # 设置初始点
nub_iter = 2 # 设置迭代次数
main()
optimization/caculation/algorithm/梯度/最速下降法(梯度法().py 0 → 100644
浏览文件 @ 20a22b00
# -*- coding: utf-8 -*-
# @File : 最速下降法.py
# @Time : 2021/4/14
# @Author : laipinyan
from sympy import *
import numpy as np
def main():
global init
for i in range(nub_iter):
print(f'__________________第{i+1}次迭代_________________________________')
#算梯度
g1 = diff(f, x1).subs([(x1, init[0]),(x2,init[1])]) # 对x偏导
g2 = diff(f, x2).subs([(x1, init[0]),(x2,init[1])]) # 对y偏导
print(f'∵∇f(x)=({diff(f, x1)},{diff(f, x2)})^T\n'
f'∴d^{i+1}=-∇f(x)=(-{g1},-{g2})^T\n'
f'∴x^{i+1}+λd^{i+1}=({x1}-{g1}λ,{x2}-{g2}λ)^T')
new_x1 = init[0] + t*(-g1)
new_x2 = init[1] + t*(-g2)
f_new = f.subs([(x1, new_x1), (x2, new_x2)]) # 新函数
grad_new = diff(f_new, t) # 对t偏导
t_val = nsolve(grad_new, t,0) # 求t
init_new = (new_x1.subs(t, t_val), new_x2.subs(t, t_val)) # 新点
f_value = f.subs([(x1, init_new[0]), (x2, init_new[1])]) # 目标函数值
print(f'令φ(λ)=f(x^{i+1}+λd^{i+1})={f_new}\n'
f'求解minφ(λ)\n'
f'令φ´(λ)={grad_new}=0\n'
f'解得λ{i+1}={t_val}\n'
f'∴x^{i+2}=x^{i+1}++λ{i+1}d^{i+1}={init_new}\n'
f' f(x^{i+2})={f_value}')
print(f"{i+1}次迭代总结：原始点位{init}，最速下降方向{np.array([-g1,-g2])},最优步长{t_val},目标点位{init_new},目标函数值{f_value}")
init=init_new
if __name__ == '__main__':
x1, x2, t = symbols('x1'), symbols('x2'), symbols('t')
f = 2*x1**2+x2**2-4*x1 # 设置初始函数
init = (0, 1) # 设置初始点
nub_iter = 2 # 设置迭代次数
main()
optimization/caculation/algorithm/梯度/牛顿法.py 0 → 100644
浏览文件 @ 20a22b00
#----牛顿法求根-----#
import numpy as np
from sympy import *
def main():
global init
L=0
while 1: #采用残差来判断
g1 = diff(f, x1).subs([(x1, init[0]), (x2, init[1])]) # 对x偏导
g2 = diff(f, x2).subs([(x1, init[0]), (x2, init[1])]) # 对y偏导
g = np.array([float(g1), float(g2)]).T
g11 = diff(diff(f, x1), x1).subs([(x1, init[0]), (x2, init[1])])
g12 = diff(diff(f, x1), x2).subs([(x1, init[0]), (x2, init[1])])
g21 = diff(diff(f, x2), x1).subs([(x1, init[0]), (x2, init[1])])
g22 = diff(diff(f, x2), x2).subs([(x1, init[0]), (x2, init[1])])
G = np.linalg.inv(np.array([[float(g11), float(g12)], [float(g21), float(g22)]]))
print(g)
if (np.linalg.norm(g)) > e:
init_new = init - np.dot(G, g)
L = L + 1 # 统计迭代次数
f_value = f.subs([(x1, init_new[0]), (x2, init_new[1])]) # 目标函数值
print(f'第{L}次迭代，原始点位为{init}，牛顿方向为{-np.dot(G, g)},目标点位：{init_new}，目标函数值：{f_value}')
init=init_new
else:
print("_____________over__________________")
print(f"最终点位x：{init}") # 输出数值解
print(f"最终f(x):{f_value}") # 验证解的正确性
print(f"最终迭代次数{L}") # 输出迭代次数
break
if __name__ == '__main__':
e = 10 ** (-9) # 误差要求
L = 0 # 初始化迭代次数
x1, x2 = symbols('x1'), symbols('x2')
f = 2*x1**2-2*x1*x2+x2**4+x2 # 设置初始函数
# f = x1**2+4*x2**2+x1*x2
init = (-2,1) # 设置初始点
main()
optimization/caculation/algorithm/梯度/高斯牛顿法.py 0 → 100644
浏览文件 @ 20a22b00
#----牛顿法求根-----#
import numpy as np
from sympy import *
def gradiant(FF,init,d,num):
print(f"原始点x{num}={init}")
new_x1 = init[0] + t * d[0]
new_x2 = init[1] + t * d[0]
f_new = FF.subs([(x1, new_x1), (x2, new_x2)]) # 新函数
grad_new = diff(f_new, t) # 对t偏导
t_val = nsolve(grad_new, t, 0) # 求t
init_new = [new_x1.subs(t, t_val), new_x2.subs(t, t_val)] # 新点
f_value = FF.subs([(x1, init_new[0]), (x2, init_new[1])]) # 目标函数值
print(f'令x^{num + 1}=x^{num}+λd^{num}，求解minφ(λ)\n'
f'解得λ={t_val}\n'
f'新点：{init_new}')
return t_val,init_new,f_value
def main():
global init
A=np.empty([len(F),len(x)])
f_val=[[]for i in range(len(F))]
FF=0*x1
num=0
while 1: #采用残差来判断
num=num+1
print(f"___________________第{num}次迭代______________________________________")
for i in range(len(F)):
for j in range(len(x)):
A[i][j]=diff(F[i],x[j]).subs([(x1, init[0]), (x2, init[1])])
f_val[i]=F[i].subs([(x1, init[0]), (x2, init[1])])
FF=FF+(F[i])**2
A=np.array(A)
f_val=np.array(f_val)
r=np.linalg.matrix_rank(A)
if r==len(x):
d=-np.dot(np.dot(np.linalg.inv(np.dot(A.T,A)),A.T),f_val)
elif r<len(x):
d=-np.dot(A.T,f_val)
t_val,init_new,f_value=gradiant(FF,init,d,num)
if np.linalg.norm(np.dot(A.T,f_val).astype('float')) >= e:
print(f'第{num}次迭代，原始点位为{init}，高斯牛顿方向为{d},目标点位：{init_new}，目标函数值：{f_value}')
init=init_new
else:
print("_____________over__________________")
print(f"最终点位x：{init_new}") # 输出数值解
print(f"最终f(x):{f_value}") # 验证解的正确性
print(f"最终迭代次数{num}") # 输出迭代次数
break
if __name__ == '__main__':
# e = 10 ** (-9) # 误差要求
x1, x2= symbols('x1,x2')
# x=[x1,x2]
# f1=2*x1+2*x2-3
# f2=x1-2*x2-1
# f3=x1+4*x2-3
# F=[f1,f2,f3]
# init=[0,1] ##定义初始点位
t=symbols("t") ##定义λ变量
# main()
FF=x1**2-2*x1+3*x2**2
init=[2,-3]
d=[4,1]
num=0
gradiant(FF, init, d, num)
optimization/caculation/algorithm/直接法/原始Powell方法.py 0 → 100644
浏览文件 @ 20a22b00
# -*- coding: utf-8 -*-
# @File : 最速下降法.py
# @Time : 2021/4/14
# @Author : laipinyan
from sympy import *
import numpy as np
def main():
global init,D
num=0
init_new=init
D = np.array(D)
while 1:
num=num+1
print(f'__________________第{num}次迭代_____________________________________________________________________________')
#算梯度
for i in range(len(D)):
new_x1 = init_new[0] + t * D[i][0]
new_x2 = init_new[1] + t * D[i][1]
f_new = f.subs([(x1, new_x1), (x2, new_x2)]) # 新函数
grad_new = diff(f_new, t) # 对t偏导
t_val = nsolve(grad_new, t, 0) # 求t
print(f'令x^({num},{i})=x^({num-1})={init_new}^T，从x^({num},{i})出发沿着方向d^({num},{i+1})={D[i]}做一维搜索，即求解minf(x^({num},{i})+td^({num},{i+1}))\n'
f'∵x^({num},{i})+td^({num},{i+1})=({new_x1},{new_x2})\n'
f'∴令φ(λ)=f(x^({num},{i})+td^({num},{i+1}))={f_new}\n'
f'令φ´(λ)={grad_new}=0\n'
f'解得λ{i + 1}={t_val}')
new_x1 = new_x1.subs(t, t_val)
new_x2 = new_x2.subs(t, t_val)# 新点
init_new = [new_x1, new_x2]
print(f'∴x^({num},{i+1})=x^({num},{i})+td^({num},{i+1})={new_x1,new_x2}^T')
print('______')
n,m=D.shape
r=np.linalg.matrix_rank(D.astype(np.float))
if r<n:
print(f'\n________________________________结果______________________________________________________\n'
f'两个搜索方向线性相关，迭代后x^({num},{i+1})=x^({num},{i})+td^({num},{i+1})={new_x1,new_x2}^T')
break
else:
d = np.array([[new_x1 - init[0], new_x2 - init[1]]])
D = np.concatenate((D, d), axis=0)
new_x1 = init_new[0] + t * d[0][0]
new_x2 = init_new[1] + t * d[0][1]
f_new = f.subs([(x1, new_x1), (x2, new_x2)]) # 新函数
grad_new = diff(f_new, t) # 对t偏导
t_val = nsolve(grad_new, t, 0) # 求t
print(
f'令x^({num},{i + 1})=x^({num - 1})={init_new}^T，从x^({num},{i})出发沿着方向d^({num},{i + 2})=x({num},{i + 1})-({num},0)={d}做一维搜索，即求解minf(x^({num},{i})+td^({num},{i + 2}))\n'
f'∵x^({num},{i + 1})+td^({num},{i + 1})=({new_x1},{new_x2})\n'
f'∴令φ(λ)=f(x^({num},{i + 1})+td^({num},{i + 2}))={f_new}\n'
f'令φ´(λ)={grad_new}=0\n'
f'解得λ{i + 2}={t_val}')
new_x1 = new_x1.subs(t, t_val)
new_x2 = new_x2.subs(t, t_val) # 新点
init_new = [new_x1, new_x2]
print(f'∴x^({num},{i + 2})=x^({num},{i})+td^({num},{i + 2})={new_x1, new_x2}^T')
if np.linalg.norm(
np.array([float(init_new[0]) - float(init[0]), float(init_new[1]) - float(init[1])])) >= dx:
D = D[1:]
init = init_new
else:
print(f'\n________________________________结果______________________________________________________\n'
f'最终点位：x^({num},{i + 1})=x^({num},{i})+td^({num},{i + 1})={new_x1, new_x2}^T')
break
if __name__ == '__main__':
x1, x2, t = symbols('x1'), symbols('x2'), symbols('t')
f = (x1+x2)**2+(x1-1)**2 # 设置初始函数
init = [2, 1] # 设置初始点
D=[[1,0],[0,1]] # 设置初始搜索方向
dx = 0.01 # 设置迭代精度
main()
optimization/caculation/algorithm/直接法/模式搜索法.py 0 → 100644
浏览文件 @ 20a22b00
# -*- coding: utf-8 -*-
# @File : 最速下降法.py
# @Time : 2021/4/14
# @Author : laipinyan
from sympy import *
import numpy as np
def main():
global init,step
num=1
print(f'给定初始点x{num},初始步长{step},加速因子{a}>=1,缩减率{b}∈(0,1),精度{dx}>0\n')
y = init
while 1:
print(f'__________________第{num}次迭代_________________________________')
print(f'令y{1}=x{num},k=1,j=1')
#轴向搜索
e=step*np.eye(n)
for i in range(n):
if f.subs([(x1,(y+e[i])[0]),(x2,(y+e[i])[1])])<f.subs([(x1,y[0]),(x2,y[1])]):
y_new = y + e[i]
print(f'∵f(y^{i+1}+ξe{i+1})={f.subs([(x1,(y+e[i])[0]),(x2,(y+e[i])[1])])} < {f.subs([(x1,y[0]),(x2,y[1])])}=f(y^{i+1})')
else:
print(f'∵f(y^{i+1}+ξe{i+1})={f.subs([(x1,(y+e[i])[0]),(x2,(y+e[i])[1])])} > {f.subs([(x1,y[0]),(x2,y[1])])}=f(y^{i+1})')
if f.subs([(x1,(y-e[i])[0]),(x2,(y-e[i])[1])])<f.subs([(x1,y[0]),(x2,y[1])]):
y_new = y - e[i]
print(f'∵f(y^{i+1}-ξe{i + 1})={f.subs([(x1,(y-e[i])[0]),(x2,(y-e[i])[1])])} < {f.subs([(x1,y[0]),(x2,y[1])])}=f(y^{i+1})')
else:
print(f'∵f(y^{i+1}-ξe{i + 1})={f.subs([(x1,(y-e[i])[0]),(x2,(y-e[i])[1])])} < {f.subs([(x1,y[0]),(x2,y[1])])}=f(y^{i+1})')
y_new = y
if f.subs([(x1,(y+e[i])[0]),(x2,(y+e[i])[1])])<f.subs([(x1,y[0]),(x2,y[1])]):
print(f'∴y^{i+2}=y^{i+1}+ξe{i+1}={y+e[i]}')
elif f.subs([(x1,(y-e[i])[0]),(x2,(y-e[i])[1])]) < f.subs([(x1,y[0]),(x2,y[1])]):
print(f'∴y^{i+2}=y^{i+1}-ξe{i+1}={y-e[i]}')
else:
print(f'∴y^{i + 2}=y^{i + 1}={y}')
y=y_new
if f.subs([(x1,y[0]),(x2,y[1])])<f.subs([(x1,init[0]),(x2,init[1])]):
print(f'∵f(y^{n+1})<f(x^{num})\n'
f'∴令x^{num+1}=y^{n+1}\n'
f'取加速方向d^{num}=x^{num+1}-x^{num}={y-init}\n'
f'模式搜索得到目标点位{y}')
init_new = y
y=init_new+a*(init_new-init)
else:
if step>dx:
step=b*step
y=init
init_new=init
print(f'∵f(y^{n + 1})>=f(x^{num})且未达到规定精度\n'
f'∴令step=β*step\n'
f'取加速方向d^{num}=x^{num + 1}-x^{num}={y - init}\n'
f'模式搜索得到目标点位{init_new}')
else:
break
init = init_new
num=num+1
if __name__ == '__main__':
x1, x2 = symbols('x1'), symbols('x2')
f = x1**2+x2**2 # 设置初始函数
n = 2 #维度数，即有多少个Xi
init = [1, 1] # 设置初始点
step = 0.25 #初始步长
a = 1 #加速因子
b = 0.2 #缩减率
dx = 0.01 #精度
main()
optimization/caculation/algorithm/线性/线性规划.py 0 → 100644
浏览文件 @ 20a22b00
from sympy import *
import numpy as np
def main(i,j):
B=np.array([A[:,i],A[:,j]]).T
cb=np.array([c[i],c[j]])
B_1=np.linalg.inv(B)
xb=np.dot(B_1,b)
f_val=np.sum(cb*xb)
print(f'______________令B=[p{i+1}p{j+1}]__________________')
print(f'B={B}\n'
f',cB={cb}\n'
f'B^-1={B_1}\n'
f'xB=B^-1 b={xb}\n'
f'得到相应解f={f_val}')
return f_val,xb
if __name__ == '__main__':
x1, x2, x3, x4, x5,x6,x7 = symbols('x1, x2, x3, x4, x5,x6,x7')
f=x1-x2
c=[1,-1,0,0,0]
A=[[1,1,1,1,0],
[-1,1,2,0,1]]
b=[5,
6]
A=np.array(A)
b=np.array(b)
aa=0
for i in range(np.shape(A)[1]):
for j in range(i+1,np.shape(A)[1],1):
f_,x_=main(i,j)
if (x_>0).all():
if aa == 0:
x_fin = x_
f_fin = f_
i_fin = i
j_fin = j
aa=1
if f_fin >= f_:
x_fin = x_
f_fin = f_
i_fin = i
j_fin = j
x=np.zeros(len(c))
x[i_fin]=x_fin[0]
x[j_fin]=x_fin[1]
print(f"\n________________结论________________________\n"
f"最优解：{x},最优值fmin={f_fin}")
optimization/caculation/apps.py 0 → 100644
浏览文件 @ 20a22b00
from django.apps import AppConfig
class CaculationConfig(AppConfig):
default_auto_field = 'django.db.models.BigAutoField'
name = 'caculation'
optimization/caculation/models.py 0 → 100644
浏览文件 @ 20a22b00
from django.db import models
# Create your models here.
optimization/caculation/tests.py 0 → 100644
浏览文件 @ 20a22b00
from django.test import TestCase
# Create your tests here.
optimization/caculation/views.py 0 → 100644
浏览文件 @ 20a22b00
from django.shortcuts import render
from django.shortcuts import HttpResponse
from .algorithm.OLS import OLS
# Create your views here.
def home(request):
return render(request,'home.html')
def page1(request):
if request.method == 'POST':
form = request.POST
A = [[float(request.POST.get('x1_1')), float(request.POST.get('x2_1'))],
[float(request.POST.get('x1_2')), float(request.POST.get('x2_2'))],
[float(request.POST.get('x1_3')), float(request.POST.get('x2_3'))]]
b = [float(request.POST.get('c_1')),
float(request.POST.get('c_2')),
float(request.POST.get('c_3'))]
output = OLS(A, b)
return render(request,'page1.html',{'output':output,'form':form})
else:
return render(request,'page1.html')
\ No newline at end of file
optimization/manage.py 0 → 100644
浏览文件 @ 20a22b00
#!/usr/bin/env python
"""Django's command-line utility for administrative tasks."""
import os
import sys
def main():
"""Run administrative tasks."""
os.environ.setdefault('DJANGO_SETTINGS_MODULE', 'optimization.settings')
try:
from django.core.management import execute_from_command_line
except ImportError as exc:
raise ImportError(
"Couldn't import Django. Are you sure it's installed and "
"available on your PYTHONPATH environment variable? Did you "
"forget to activate a virtual environment?"
) from exc
execute_from_command_line(sys.argv)
if __name__ == '__main__':
main()
optimization/optimization/asgi.py 0 → 100644
浏览文件 @ 20a22b00
"""
ASGI config for optimization project.
It exposes the ASGI callable as a module-level variable named ``application``.
For more information on this file, see
https://docs.djangoproject.com/en/4.2/howto/deployment/asgi/
"""
import os
from django.core.asgi import get_asgi_application
os.environ.setdefault('DJANGO_SETTINGS_MODULE', 'optimization.settings')
application = get_asgi_application()
optimization/optimization/settings.py 0 → 100644
浏览文件 @ 20a22b00
import os
"""
Django settings for optimization project.
Generated by 'django-admin startproject' using Django 4.2.
For more information on this file, see
https://docs.djangoproject.com/en/4.2/topics/settings/
For the full list of settings and their values, see
https://docs.djangoproject.com/en/4.2/ref/settings/
"""
from pathlib import Path
# Build paths inside the project like this: BASE_DIR / 'subdir'.
BASE_DIR = Path(__file__).resolve().parent.parent
# Quick-start development settings - unsuitable for production
# See https://docs.djangoproject.com/en/4.2/howto/deployment/checklist/
# SECURITY WARNING: keep the secret key used in production secret!
SECRET_KEY = 'django-insecure-e_jsrlgm4cle*t1zjqmrqnj$#rh84n0zvz-%*$av*j^_t(&yz8'
# SECURITY WARNING: don't run with debug turned on in production!
DEBUG = True
ALLOWED_HOSTS = []
# Application definition
INSTALLED_APPS = [
'django.contrib.admin',
'django.contrib.auth',
'django.contrib.contenttypes',
'django.contrib.sessions',
'django.contrib.messages',
'django.contrib.staticfiles',
]
MIDDLEWARE = [
'django.middleware.security.SecurityMiddleware',
'django.contrib.sessions.middleware.SessionMiddleware',
'django.middleware.common.CommonMiddleware',
'django.middleware.csrf.CsrfViewMiddleware',
'django.contrib.auth.middleware.AuthenticationMiddleware',
'django.contrib.messages.middleware.MessageMiddleware',
'django.middleware.clickjacking.XFrameOptionsMiddleware',
]
ROOT_URLCONF = 'optimization.urls'
TEMPLATES = [
{
'BACKEND': 'django.template.backends.django.DjangoTemplates',
'DIRS': [os.path.join(BASE_DIR, 'templates')]
,
'APP_DIRS': True,
'OPTIONS': {
'context_processors': [
'django.template.context_processors.debug',
'django.template.context_processors.request',
'django.contrib.auth.context_processors.auth',
'django.contrib.messages.context_processors.messages',
],
},
},
]
WSGI_APPLICATION = 'optimization.wsgi.application'
# Database
# https://docs.djangoproject.com/en/4.2/ref/settings/#databases
DATABASES = {
'default': {
'ENGINE': 'django.db.backends.sqlite3',
'NAME': BASE_DIR / 'db.sqlite3',
}
}
# Password validation
# https://docs.djangoproject.com/en/4.2/ref/settings/#auth-password-validators
AUTH_PASSWORD_VALIDATORS = [
{
'NAME': 'django.contrib.auth.password_validation.UserAttributeSimilarityValidator',
},
{
'NAME': 'django.contrib.auth.password_validation.MinimumLengthValidator',
},
{
'NAME': 'django.contrib.auth.password_validation.CommonPasswordValidator',
},
{
'NAME': 'django.contrib.auth.password_validation.NumericPasswordValidator',
},
]
# Internationalization
# https://docs.djangoproject.com/en/4.2/topics/i18n/
LANGUAGE_CODE = 'en-us'
TIME_ZONE = 'UTC'
USE_I18N = True
USE_TZ = True
# Static files (CSS, JavaScript, Images)
# https://docs.djangoproject.com/en/4.2/howto/static-files/
STATIC_URL = 'static/'
STATICFILES_DIRS=[
os.path.join(BASE_DIR,'static')
]
# Default primary key field type
# https://docs.djangoproject.com/en/4.2/ref/settings/#default-auto-field
DEFAULT_AUTO_FIELD = 'django.db.models.BigAutoField'
optimization/optimization/urls.py 0 → 100644
浏览文件 @ 20a22b00
"""
URL configuration for optimization project.
The `urlpatterns` list routes URLs to views. For more information please see:
https://docs.djangoproject.com/en/4.2/topics/http/urls/
Examples:
Function views
1. Add an import: from my_app import views
2. Add a URL to urlpatterns: path('', views.home, name='home')
Class-based views
1. Add an import: from other_app.views import Home
2. Add a URL to urlpatterns: path('', Home.as_view(), name='home')
Including another URLconf
1. Import the include() function: from django.urls import include, path
2. Add a URL to urlpatterns: path('blog/', include('blog.urls'))
"""
from django.contrib import admin
from django.urls import path
from caculation import views
urlpatterns = [
path('admin/', admin.site.urls),
path('',views.home),
path('page1/', views.page1)
]
optimization/optimization/wsgi.py 0 → 100644
浏览文件 @ 20a22b00
"""
WSGI config for optimization project.
It exposes the WSGI callable as a module-level variable named ``application``.
For more information on this file, see
https://docs.djangoproject.com/en/4.2/howto/deployment/wsgi/
"""
import os
from django.core.wsgi import get_wsgi_application
os.environ.setdefault('DJANGO_SETTINGS_MODULE', 'optimization.settings')
application = get_wsgi_application()
optimization/templates/home.html 0 → 100644
浏览文件 @ 20a22b00
<!DOCTYPE html>
<html lang="en">
<head>
<meta charset="UTF-8">
<title>最优化计算器</title>
</head>
<body>
<button type="button" onclick="location.href='/page1/'">最小二乘法</button>
<button type="button">按钮2</button>
</body>
</html>
\ No newline at end of file
optimization/templates/page1.html 0 → 100644
浏览文件 @ 20a22b00
<!DOCTYPE html>
<html lang="en">
<head>
<meta charset="UTF-8">
<title>最小二乘法</title>
</head>
<body>
<h1>请输入你要计算的数据：</h1>
<form action="/page1/" method="post">
{% csrf_token %}
<input type="number" step=0.001 name="x1_1" value="{{ form.x1_1|default:'' }}">x1+
<input type="number" step=0.001 name="x2_1" value="{{ form.x2_1|default:'' }}">x2=
<input type="number" step=0.001 name="c_1" value="{{ form.c_1|default:'' }}"><br>
<input type="number" step=0.001 name="x1_2" value="{{ form.x1_2|default:'' }}">x1+
<input type="number" step=0.001 name="x2_2" value="{{ form.x2_2|default:'' }}">x2=
<input type="number" step=0.001 name="c_2" value="{{ form.c_2|default:'' }}"><br>
<input type="number" step=0.001 name="x1_3" value="{{ form.x1_3|default:'' }}">x1+
<input type="number" step=0.001 name="x2_3" value="{{ form.x2_3|default:'' }}">x2=
<input type="number" step=0.001 name="c_3" value="{{ form.c_3|default:'' }}"><br>
<button type="submit">开始计算</button>
</form>
<h1>计算结果：</h1>
<p>{{ output|linebreaksbr }}</p>
</body>
</html>
\ No newline at end of file
requirements.txt
浏览文件 @ 20a22b00
Django==4.2
numpy==1.24.3
sympy==1.11.1
\ No newline at end of file
Markdown is supported
0% 或 .
You are about to add 0 people to the discussion. Proceed with caution.
先完成此消息的编辑！
想要评论请 注册 或