samer@0: package samer.mds;
samer@0: 
samer@0: import samer.core.*;
samer@0: import samer.core.types.*;
samer@0: import samer.maths.*;
samer@0: import samer.tools.*;
samer@0: 
samer@0: /**
samer@0:      This is an alternative version of MDS that uses less memory.
samer@0:      Rather than keep a record of the force in each link, it accumulates
samer@0:      the link forces on a per-object basis.
samer@0:  */
samer@0: public class MDS extends MDSBase
samer@0: {
samer@0: 	Metric		metric;		// metric function
samer@0: 	Stress		stress;		// stress function
samer@0: 
samer@0: 	public interface Metric {
samer@0: 		/** return distance between x and y, put y-x in r */
samer@0: 		double d(double [] x, double [] y, double [] r);
samer@0: 	}
samer@0: 
samer@0: 	public interface Stress {
samer@0: 		/** accumulate stress due to a link and return dStress/dcurrent */
samer@0: 		double add( double target, double current);
samer@0: 
samer@0: 		/** return normalised stress and reset for next time */
samer@0: 		double done();
samer@0: 
samer@0: 		/** return last value returned by done() */
samer@0: 		double get();
samer@0: 	}
samer@0: 
samer@0: 	public static class SamerStress implements Stress {
samer@0: 		double S=0, a=0;
samer@0: 		int			n=0;
samer@0: 
samer@0: 		public double add( double target, double current) {
samer@0: 			a += sqr(current/target-1); n++;
samer@0: 			return (current-target)/(target*target);
samer@0: 		}
samer@0: 		public double done() { S=a/n; a=0; n=0; return S; }
samer@0: 		public double get() { return S; }
samer@0: 		private final static double sqr(double t) {return t*t;}
samer@0: 	}
samer@0: 
samer@0: 	public MDS(Matrix p)
samer@0: 	{
samer@0: 		super(p);
samer@0: 
samer@0: 		Shell.push(node);
samer@0: 		metric=new Euclidean();
samer@0: 		stress=new SamerStress();
samer@0: 		Shell.pop();
samer@0: 	}
samer@0: 
samer@0: 	public void setMetric(Metric m) { metric=m; }
samer@0: 
samer@0: 	public void run()
samer@0: 	{
samer@0: 		double [][] _P=P.getArray();
samer@0: 		double	d, dS;
samer@0: 		int 		i, j, k;
samer@0: 
samer@0: 		// accumulate all forces on a per object basis
samer@0: 
samer@0: 		for (k=0; k<N; k++) Mathx.zero(F[k]); // zero out object forces
samer@0: 		for (k=0; k<M; k++) {			// for each link
samer@0: 			i=l1[k]; j=l2[k]; 				// object indices
samer@0: 			d=metric.d(_P[i],_P[j],f); 		// current distance, vector from i to j in f
samer@0: 			dS=stress.add(D[k],d);
samer@0: 
samer@0: 			if (d>0) { // ignore if current distance is zero
samer@0: 				// replace f with dS/dd * grad d
samer@0: 				Mathx.mul(f,dS/d);
samer@0: 
samer@0: 				// f contains 'force' on j from i
samer@0: 				// accumulate this force for both objects (only E dimensions used)
samer@0: 				Mathx.add(E,F[i],f);
samer@0: 				Mathx.sub(E,F[j],f);
samer@0: 			}
samer@0: 		}
samer@0: 		S.set(stress.done());
samer@0: 
samer@0: 		// now move objects at the ends of each link,
samer@0: 		// but only the first E dimensions
samer@0: 		for (k=0; k<N; k++) Mathx.muladd(E,_P[k],F[k],rate.value);
samer@0: 		P.changed();
samer@0: 	}
samer@0: }