Winner Paulo Uribe (turbf1)

2002-05-23 00:00:00 UTC

# smallopt3

Status: Passed
Results: Average strain = 2443.605
CPU Time: 134.944
Score: 8562.81
Submitted at: 2002-05-21 11:13:02 UTC
Scored at: 2002-05-21 11:16:35 UTC

Current Rank: 723rd

Code
```function xy = solver(kd,bx)

% catch trivial case:
if size(kd,1) < 2
xy = [mean(bx(1:2)) mean(bx(3:4))];
return
end

%intialize xy
for iter = 1:10
best = Inf;
thisxy = [rand(length(kd),1)*(bx(2)-bx(1))+bx(1) rand(length(kd),1)*(bx(4)-bx(3))+bx(3)];
s=findstrain(thisxy(:),kd,bx);
if s < best
xy = thisxy;
best=s;
end
end

[s,M] = findstrain(xy(:),kd,bx);
xyold=Inf*xy;

if s==0
return
end

iter = 1;
while iter
fx=0;fy=0; % intialize force vectors
for indx = 1:length(kd)
dx=xy(:,1)-xy(indx,1);
dy=xy(:,2)-xy(indx,2);
N = sqrt(dx.^2 + dy.^2);
N(N==0)=1;
dx=dx./N;
dy=dy./N;
% now have direction vectors (dx,dy) from xy(indx,:) to each other pt
dx = dx .* M(:,indx);
dy = dy .* M(:,indx);
dx = sum(dx);
dy = sum(dy);
% now have strain vectors
xy(indx,:) = xy(indx,:) + .05*[dx dy]; % move point
% clip if outside box
xy(indx,1) = min(max(xy(indx,1),bx(1)),bx(2));
xy(indx,2) = min(max(xy(indx,2),bx(3)),bx(4));

end

iter = iter+1;
Nx = abs(xy(:,1)-xyold(:,1))/(bx(2)-bx(1));
Ny = abs(xy(:,2)-xyold(:,2))/(bx(4)-bx(3));
N = max([Nx;Ny]);
xyold=xy;
if iter > 50 | N < .0001
iter = 0;
end

end

function [s,strainMatrix]=findstrain(inp,kd,bx)
xy = reshape(inp,[length(inp)/2 2]);
% Calculate distance matrix
x=xy(:,1);
y=xy(:,2);

% catch points outside box
if any(x<bx(1) | x>bx(2) | y<bx(3) | y> bx(4))
s = Inf;
return
end

[XY1,XY2] = meshgrid(x,y);
dx = XY1 - XY1';
dy = XY2 - XY2';
dist = sqrt(dx.^2 + dy.^2);

% Calculate strain matrix
strainMatrix = dist - kd;
strainMatrix(kd < 0) = 0;
s = sum(abs(strainMatrix(:)));
return
```