TomSym Distribution of electrons on a sphere: Difference between revisions
(Created page with "{{Part Of Manual|title=the TomSym Manual|link=TomSym Manual}} TomSym implementation of GAMS Example (ELEC,SEQ=230) Given n electrons, find the equilibrium state dist...") |
No edit summary |
||
Line 77: | Line 77: | ||
<pre> | <pre> | ||
Problem type appears to be: con | Problem type appears to be: con | ||
Time for symbolic processing: 0. | Time for symbolic processing: 0.11941 seconds | ||
Starting numeric solver | Starting numeric solver | ||
===== * * * =================================================================== * * * | ===== * * * =================================================================== * * * | ||
TOMLAB - TOMLAB Development license 999007. Valid to 2011-12-31 | TOMLAB - TOMLAB Development license 999007. Valid to 2011-12-31 | ||
===================================================================================== | ===================================================================================== | ||
Problem: --- 1: Problem 1 f_k 243. | Problem: --- 1: Problem 1 f_k 243.812760283686030000 | ||
sum(|constr|) 0. | sum(|constr|) 0.000000003100005141 | ||
f(x_k) + sum(|constr|) 243. | f(x_k) + sum(|constr|) 243.812760286786040000 | ||
f(x_0) | f(x_0) 292.427334354468940000 | ||
Solver: snopt. EXIT=0. INFORM=1. | Solver: snopt. EXIT=0. INFORM=1. | ||
Line 91: | Line 91: | ||
Optimality conditions satisfied | Optimality conditions satisfied | ||
FuncEv | FuncEv 82 GradEv 80 ConstrEv 80 ConJacEv 80 Iter 79 MinorIter 129 | ||
CPU time: 0. | CPU time: 0.093601 sec. Elapsed time: 0.090000 sec. | ||
</pre> | </pre> |
Latest revision as of 09:34, 8 November 2011
This page is part of the TomSym Manual. See TomSym Manual. |
TomSym implementation of GAMS Example (ELEC,SEQ=230)
Given n electrons, find the equilibrium state distribution (of minimal Coulomb potential) of the electrons positioned on a conducting sphere.
This model is from the COPS benchmarking suite.
The number of electrons can be specified using the parameter np.
Dolan, E D, and More, J J, Benchmarking Optimization Software with COPS. Tech. rep., Mathematics and Computer Science Division, 2000.
Morris, J R, Deaven, D M, and Ho, K M, Genetic Algorithm Energy Minimization for Point Charges on a Sphere. Phys. Rev. B. 53 (1996), R1740-R1743.
Saff, E B, and Kuijlaars, A, Distributing Many Points on the Sphere. Math. Intelligencer 19 (1997), 5-11.
np = 25;
% Coordinates of the electrons
x = tom('x',np,1);
y = tom('y',np,1);
z = tom('z',np,1);
toms potential
% The following two methods of computing the potential energy are
% mathematically equivalent.
useSlowMethod = false; % Set to "true" to try the slower method. (Tip: Use
% a small value for np when you try this.)
if useSlowMethod
% The slow way of computing the sum
% Using nested for loops, we create a symbolic expression containing
% (np-1)*(np-2)/2 terms. Each of these terms will be symbolically
% differentiated to compute the Jacobian and Hessian matrices for use
% by the solver.
% For np=10, the autogenerated Jacobian of this function contains more
% than 15 kilobytes of Matlab code (and the Hessian more than 30 kb).
% TomSym will need several minutes to generate this code, which is
% orders of magnitude more than the fraction of a second required to
% actually solve the numeric problem.
potential = 0;
for i=1:np
for j=(i+1):np
potential = potential + 1.0/sqrt((x(i)-x(j))^2 +...
(y(i)-y(j))^2 + (z(i)-z(j))^2);
end
end
else
% Using vectorized code saves a lot of time in handling the symbolic
% expressions, because the symbolic derivative is computed for the one
% vectorized expression, not on hundreds of individual terms.
% (The autogenerated Hessian of this function contains less than 1 kb
% of Matlab code, regardless of the value of np.)
[ii,jj] = meshgrid(1:np,1:np);
[i,j] = find(ii<jj);
potential = sum(1.0./sqrt((x(i)-x(j)).^2 + (y(i)-y(j)).^2 +...
(z(i)-z(j)).^2));
end
% Ball constraint
ball = (x.^2 + y.^2 + z.^2 == 1);
% Initial guess:
theta = 2*pi*rand(np,1);
phi = pi*rand(np,1);
x0 = {x == cos(theta).*sin(phi)
y == sin(theta).*sin(phi)
z == cos(phi)};
options = struct;
options.solver = 'snopt';
solution = ezsolve(potential,ball,x0,options);
Problem type appears to be: con Time for symbolic processing: 0.11941 seconds Starting numeric solver ===== * * * =================================================================== * * * TOMLAB - TOMLAB Development license 999007. Valid to 2011-12-31 ===================================================================================== Problem: --- 1: Problem 1 f_k 243.812760283686030000 sum(|constr|) 0.000000003100005141 f(x_k) + sum(|constr|) 243.812760286786040000 f(x_0) 292.427334354468940000 Solver: snopt. EXIT=0. INFORM=1. SNOPT 7.2-5 NLP code Optimality conditions satisfied FuncEv 82 GradEv 80 ConstrEv 80 ConJacEv 80 Iter 79 MinorIter 129 CPU time: 0.093601 sec. Elapsed time: 0.090000 sec.