# Diffing "try12" and "mytry"

 Title: try12 mytry Author: gaof DAJ Submitted: 2002-05-23 10:30:32 UTC 2002-05-23 14:57:12 UTC Status: Passed Passed Score: 7809.33 7809.29 Result: Average strain = 2032.285 Average strain = 2032.256 CPU Time: 171.107 173.089 Code: ```function xy = solver(kd,bx) npts = length(kd); if ~any(kd(:)>0), xy = zeros(npts,2); return end dx=bx(2); dy=bx(4); D=dx+1i*dy; pm=D/2; xy=zeros(1,npts)+pm; d=0.1*max(dx,dy); x=0:d:dx; y=0:d:dy; [xx,yy] = ndgrid(x,1i*y); [nx,ny]=size(xx); onesx=ones(1,nx); onesy=ones(1,ny); XX=xx(:)+yy(:); dr=abs(D); kd(find(kd>dr))=dr; S=(kd-eye(npts))>=0; nX=length(XX); xOnes=ones(nX,1); xZeros=zeros(nX,1); start=1; iones=ones(1,npts); sCut = 1; olds = 1.0e30; iterbigmax = 1; for iterbig = 1:iterbigmax [dum,p]=sort(-max(kd)); S=S(p,p); kd=kd(p,p); [v,w]=sort(p); xxx=10; if npts<50 xxx=11; end; for iter=1:xxx for iter=1:10+npts<50; for index = 1:npts I = find(S(:,index)); if ~isempty(I) st = sum(abs(abs(XX(:,ones(size(I)))-xy(xOnes,I))-kd(xZeros+index,I)),2); [null,minloc] = min(st); xy(index) = XX(minloc); end end [s,M]=findstrain(S,xy,kd); if((abs(s-olds)/(s+1)) < 0.002) | s < sCut break; end olds = s; end xy=xy.'; xy=xy(w,:); kd=kd(w,w); S=S(w,w); olds = 1.0e50; oldxy = xy; alpha = 1; for iter=1:40 [s,M]=findstrain(S,xy,kd); if( iter > 5 & abs(s-olds) / (s+1) < 0.0001 ) | s < sCut xy = oldxy; break; end if( s > olds ) xy = oldxy; break; end oldxy = xy; olds = s; for indx = 1:npts I = find(S(:,index)>0); dX=xy(I)-xy(indx); Dx=mean(dX.*(1-kd(I,index)./(abs(dX)+eps))); % Dx=mean(dX./(abs(dX)+eps) .* M(I,indx)); xy(indx) = xy(indx) + alpha*Dx; xy(indx) = min(dx,max(0,real(xy(indx))))+1i*min(dy,max(0,imag(xy(indx)))); end end if( iterbig < iterbigmax ) xy=xy.'; end end alpha = 0.1; ngrad = 40; obj = 1.0e20*ones(ngrad,1); xybig = zeros(npts,ngrad); [objnew,strainMatrix]=findstrain(S,xy,kd); obj(1) = objnew; xybig(:,1) = xy; gradient = zeros(npts,1); for ij=2:ngrad [s,sM,dM]=findstrain(S,xy,kd); gradient(:) = sum(dM./(abs(dM)+eps).*sign(sM)); xy = xy - alpha * gradient; xy = min(max(real(xy),bx(1)),bx(2)) + j * min(max(imag(xy),bx(3)),bx(4)); [objnew,strainMatrix]=findstrain(S,xy,kd); obj(ij) = objnew; xybig(:,ij) = xy; if( obj(ij) > obj(ij-1) ) xy = xybig(:,ij-1); alpha = alpha / 10; else alpha = alpha * 2; end [bestobj,index] = min(obj); if( ij > 3 & index < 2 ) | bestobj < sCut break end end [bestobj,index] = min(obj); xy = xybig(:,index); xy=[real(xy) imag(xy)]; function [s,strainMatrix,dX]=findstrain(S,xy,kd) sI=find(S); strainMatrix=S; dX=S; [X1,X2]=meshgrid(xy); dX(sI)=X1(sI)-X2(sI); strainMatrix(sI) = abs(dX(sI))-kd(sI); s = sum(abs(strainMatrix(sI)));```