how to incorporate boundary conditions as constraints in fmincon optimization?

6 次查看(过去 30 天)
Trying to optimize the six parameters in a trial non-linear wavefunction for particle in a box to get the minimum energy. I trying to use fmincon. but i am getting an error that says - "Unrecognized function or variable 'variational_pi1DB_fun'." Somehow I am not able to incorporate non-linear constarints in the optimization problem. Can someone help me please! Pleas let me know if other methods of optimization is more helpful!
clear
clc
%defining optimization variables and an optimization problem object.
a = optimvar('a','LowerBound',30,"UpperBound",30);
b = optimvar('b','LowerBound',30,"UpperBound",30);
c = optimvar('c','LowerBound',30,"UpperBound",30);
d = optimvar('d','LowerBound',30,"UpperBound",30);
e = optimvar('e','LowerBound',30,"UpperBound",30);
f = optimvar('f','LowerBound',30,"UpperBound",30);
prob = optimproblem;
L = 5;
% w=1;
% v=1.5;
% constraints
% wavefunction zero @ x = 0 and x = L;
cons1 = double(subs(variational_pi1DB_fun,x,0)) == 0;
cons2 = double(subs(variational_pi1DB_fun,x,L)) == 0;
prob.Constraints.cons1 = cons1;
prob.Constraints.cons2 = cons2;
x0.a = -0.00006;
x0.b = -19;
x0.c = -9;
x0.d = 6;
x0.e = -0.90;
x0.f = 0.043;
%new variables
t_pi1DB = a + b*x + c*x^2 + d*x^3 + e*x^4 + f*x^5;
% equation 4.2 nominator
t_pi1DB_D1 = diff(t_pi1DB,x); % 1st derivative
t_pi1DB_D2 = diff(t_pi1DB_D1,x); % 2nd derivative
t_pi1DB_DD = ((-hbar^2)./(2*m)).*t_pi1DB_D2; % applying the hamiltonian to trial wfn Psi
t_pi1DB_c = (t_pi1DB)'; % conjugate of psi
t_pi1DB_prod = t_pi1DB_c.*t_pi1DB_DD; % multiplying the Psi* from left side
% equation 4.2 denominator
t_pi1DB_prod2 = t_pi1DB_c.*t_pi1DB; % hamiltonian not applied
% integration in nominator % ground state
probability1 = int(real(t_pi1DB_prod),x,0,L); % integrate from 0 to L
disp('probability = '),disp(probability1)
probability1 = subs(probability1,{hbar,m,L},{1,1,5}); % evaluate the probability at n = 1
disp('probability1 = '),disp(probability1)
% integration in denominator
probability2 = int(real(t_pi1DB_prod2),x,0,L);
probability2 = subs(probability2,L,5); % evaluate the probability at n = 1
disp('probability2 = '),disp(probability2)
norm_cons = sqrt(probability2);
variational_pi1DB_fun = t_pi1DB/norm_cons;
% application of equation 4.2
variational_Energy = probability1/probability2;
%objective function as an expression in the optimization variables.
P = variational_Energy;
%the objective function in prob.
prob.Objective = P;
show(prob)
sol = solve(prob,x0);
disp('sol = '),disp(sol)
  3 个评论
Saurav Parmar
Saurav Parmar 2023-12-7
I tried defining constraints in % Constraints block ...it is just below the % opt variables block(top). Everything after %New varables is for the definition of Optim.problem..
Torsten
Torsten 2023-12-7
But you can't use a variable before you define it.
In the constraint section
% constraints
% wavefunction zero @ x = 0 and x = L;
cons1 = double(subs(variational_pi1DB_fun,x,0)) == 0;
cons2 = double(subs(variational_pi1DB_fun,x,L)) == 0;
variational_pi1DB_fun is not yet defined.

请先登录,再进行评论。

采纳的回答

Torsten
Torsten 2023-12-7
编辑:Torsten 2023-12-7
I think with the mixture of numerical, symbolic and optimization variables it is easier to set up the problem directly.
x0 = [ -0.00006; -19; -9; 6; -0.90; 0.043];
L = 5;
hbar = 1;
m = 1;
sol = fmincon(@(x)obj(x,L,hbar,m),x0,[],[],[],[],[],[],@(x)nonlcon(x,L,hbar,m))
Local minimum found that satisfies the constraints. Optimization completed because the objective function is non-decreasing in feasible directions, to within the value of the optimality tolerance, and constraints are satisfied to within the value of the constraint tolerance.
sol = 6×1
-0.0000 -24.1985 -0.6429 2.1927 -0.2191 -0.0000
function variational_Energy = obj(xnum,L,hbar,m)
a = xnum(1);
b = xnum(2);
c = xnum(3);
d = xnum(4);
e = xnum(5);
f = xnum(6);
syms x
%new variables
t_pi1DB = a + b*x + c*x^2 + d*x^3 + e*x^4 + f*x^5;
% equation 4.2 nominator
t_pi1DB_D1 = diff(t_pi1DB,x); % 1st derivative
t_pi1DB_D2 = diff(t_pi1DB_D1,x); % 2nd derivative
t_pi1DB_DD = ((-hbar^2)./(2*m)).*t_pi1DB_D2; % applying the hamiltonian to trial wfn Psi
t_pi1DB_c = (t_pi1DB)'; % conjugate of psi
t_pi1DB_prod = t_pi1DB_c.*t_pi1DB_DD; % multiplying the Psi* from left side
% equation 4.2 denominator
t_pi1DB_prod2 = t_pi1DB_c.*t_pi1DB; % hamiltonian not applied
% integration in nominator % ground state
probability1 = int(real(t_pi1DB_prod),x,0,L); % integrate from 0 to L
% integration in denominator
probability2 = int(real(t_pi1DB_prod2),x,0,L);
% application of equation 4.2
variational_Energy = probability1/probability2;
variational_Energy = double(variational_Energy);
end
function [C,Ceq] = nonlcon(xnum,L,hbar,m)
a = xnum(1);
b = xnum(2);
c = xnum(3);
d = xnum(4);
e = xnum(5);
f = xnum(6);
syms x
%new variables
t_pi1DB = a + b*x + c*x^2 + d*x^3 + e*x^4 + f*x^5;
t_pi1DB_c = (t_pi1DB)'; % conjugate of psi
% equation 4.2 denominator
t_pi1DB_prod2 = t_pi1DB_c.*t_pi1DB; % hamiltonian not applied
% integration in denominator
probability2 = int(real(t_pi1DB_prod2),x,0,L);
% constraints
% wavefunction zero @ x = 0 and x = L;
norm_cons = sqrt(probability2);
variational_pi1DB_fun = t_pi1DB/norm_cons;
Ceq(1) = double(subs(variational_pi1DB_fun,x,0));
Ceq(2) = double(subs(variational_pi1DB_fun,x,L));
C = [];
end
  2 个评论
Saurav Parmar
Saurav Parmar 2023-12-7
Thanks for your insightful answer...Can I ask you for some solid resources for MATLAB optimizaton routine learning/training. It looks like I am somewhat good at basic MATLAB but woefully unprepared for Optimization problems. Thanks again!!

请先登录,再进行评论。

更多回答(1 个)

Matt J
Matt J 2023-12-7
Use fcn2optimexpr to implement nonlinear constraints.
  1 个评论
Saurav Parmar
Saurav Parmar 2023-12-7
Thanks Matt for this option. It led me know that some variables had to be defined as symbolic in my problem. Ans I was completely unaware of this function. Thanks.

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Nonlinear Optimization 的更多信息

产品


版本

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by