fgridmin.m

global minimization utility
1.6K Downloads
Updated 31 Mar 2016

View License

function [i,i_] = fgridmin(y)
% Search an N-dimensional array y for points that may be proximate to a global minimum
% of the interpolated array; return the proximate points' N-D subscripts (i) in y. Also
% optionally return the subscripts (i_) for discrete local minima in y.
%
% syntax:
% i = fgridmin(y);
% [i,i_] = fgridmin(y);
%
% input arg:
%
% y: N-dimensional array, real double (See validation conditions below.)
%
% output args:
%
% i: size-[M, N] matrix of integer
%
% i_: size-[M_,N] matrix of integer (optional output)
%
% N is defined above such that size(y,N) > 1 and size(y,N+m) == 1 for all m > 0. If y is
% a scalar (case N = 0), or if size(y,m) < 3 for any m <= N, then i will be empty.
% Otherwise, i will be a size-[M,N] matrix whose rows represent N-dimensional array
% subscripts of points in y that may be proximate to a global minimum of the
% interpolated y array. (M will be zero if none is found.) The subscripts will only
% address interior points in y, i.e., 1 < i(m,j) < size(y,j) for all m, j. Also, if y
% contains any NaNs, the NaN points and points adjacent to NaNs will be excluded from
% the search. (Use NaNs to represent "missing data".)
%
% If the i_ result is requested, it will similarly contain subscripts for discrete local
% minima in y. The local minima are determined by comparing adjcent y points, including
% diagonal adjacencies. (If two adjacent points have equal y values, only one will be
% counted as a local minimum.) The i result will comprise a subset of the i_ rows.
%
% Version 04/09/2006
% Author: Kenneth C. Johnson
% software.kjinnovation.com
%
% See ZipInterp.pdf, Appendix ("Interpolator application example - global
% optimization"), on the KJ Innovation website for an application example illustrating
% the use of fgridmin.m.
%
% Get y size; define N.
size_y = size(y);
N = length(size_y);
if size_y(N)==1
N = N-1; % column vector
if size_y(N)==1
N = N-1; % scalar
end
end

% If y has no interior points return an empty subscript list.
if N==0 || any(size_y(1:N) < 3)
i = zeros(0,N);
i_ = i;
return
end

% Construct a list of index offsets for adjacent point comparisons. ("Adjacent" includes
% diagonal adjacencies.) Each index offset corresponds to a size-[1,N] subscript offset
% with elements 0, +1, or -1 (not all zero), and with the last non-zero subscript equal
% to +1. (This only includes half of the comparison offsets; the others are obtained by
% sign-inverting the offset list.)
%
% Also construct indices (n) for interior points in y.
%
offset = [];
n = 0;
stride = 1;
for j = 1:N
offset = [offset;[-offset;0;offset]+stride];
n_ = (1:size_y(j)-2)*stride;
n = repmat(n,size(n_))+repmat(n_,size(n));
n = n(:);
stride = stride*size_y(j);
end
n = n+1;

% Select indices of discrete local minima in y; restrict n to selected indices.
for m = 1:length(offset)
% Note: The first comparison below is >=, and second comparison is >, so that if two
% adjacent points have equal y values only one will be counted as a local minimum.
% Negative logic is used (~(...<...). ~(...>=...)) to filter out NaNs and points
% adjacent to NaNs. (Comparisons to NaN always yield false.)
n(~(y(n)<y(n-offset(m)))) = [];
n(~(y(n)<=y(n+offset(m)))) = [];
if isempty(n)
% no local minima
i = zeros(0,N);
i_ = i;
return
end
end

if nargout >= 2
% Convert flat indices n to N-D subscripts i_.
[i_{1:N}] = ind2sub(size_y,n);
i_ = cell2mat(i_);
end

% Determine an estimated upper bound (d) on the variation of the interpolated y from
% each y minimum within a half-step subscript range from each minimum. The d estimate
% is half the maximum difference between the y minimum and any adjacent y value, and the
% interpolated y values proximate to the minima will be assumed to be between y(n)-d and
% y(n)+d.
y_n = y(n);
d = y_n;
for m = 1:length(offset)
d = max(d,y(n-offset(m)));
d = max(d,y(n+offset(m)));
end
d = 0.5*(d-y_n);

% A conservative upper bound on the interpolated global minimum is min(y(n)+d); select
% local minima whose proximate lower bounds (y(n)-d) are within this limit.
n(y_n-d > min(y_n+d)) = [];
clear y_n d

% Convert flat indices n to N-D subscripts i.
[i{1:N}] = ind2sub(size_y,n);
i = cell2mat(i);

Cite As

Kenneth Johnson (2024). fgridmin.m (https://www.mathworks.com/matlabcentral/fileexchange/10755-fgridmin-m), MATLAB Central File Exchange. Retrieved .

MATLAB Release Compatibility
Created with R2006a
Compatible with any release
Platform Compatibility
Windows macOS Linux
Categories
Find more on Global or Multiple Starting Point Search in Help Center and MATLAB Answers

Community Treasure Hunt

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

Start Hunting!
Version Published Release Notes
1.0.0.0

Updated to add BSD License.