cgls.hh

来自「COOOL：CWP面向对象最优化库（CWP Object Oriented Op」· HH 代码 · 共 84 行

//============================================================
// COOOL           version 1.1           ---     Nov,  1995
//   Center for Wave Phenomena, Colorado School of Mines
//============================================================
//
//   This code is part of a preliminary release of COOOL (CWP
// Object-Oriented Optimization Library) and associated class 
// libraries. 
//
// The COOOL library is a free software. You can do anything you want
// with it including make a fortune.  However, neither the authors,
// the Center for Wave Phenomena, nor anyone else you can think of
// makes any guarantees about anything in this package or any aspect
// of its functionality.
//
// Since you've got the source code you can also modify the
// library to suit your own purposes. We would appreciate it 
// if the headers that identify the authors are kept in the 
// source code.
//
#ifndef CGLS_HH
#define CGLS_HH

#include "QuaOptima.hh"

//=============================
// author:  H. Lydia Deng, 03/04/94, 03/13/94
//=============================
//@Man:
//@Memo: Least Squared Conjugate Gradient 
/*@Doc:
    The conjugate gradient implemented here should be used for the 
    solution of the normal equations $A^T.A.x = A^T.y$. It is coded 
    such that the product $A^T.A$ is never performed, to avoid
    numerical instabilities and non sparse matrices. This procedure 
    comes straight from the classical paper "Methods of conjugate 
    gradients for solving linear systems:, 1952, NBS J. Research
    by Hesteness and Stiefel."

*/

class LSConjugateGradient : public QuadraticOptima {

  protected:
    double		alpha;
    double		scale;
    double		scaleOld;
    Vector<double>*	search;
    Vector<double>*	modelError;
    Vector<double>*	dataError;
    int	m, n;

    void		getModelError(Vector<double>&);
    void		conjugateDirection();
    void		upDating();

 public:
    //@ManMemo: a constructor
    LSConjugateGradient(/// dimension of model space
			int n, 
			///pointer to the forward operator (matrix)
			LinearForward* p, 
			///pointer to the observed data vector
			Vector<double>* data, 
			/// maximum number of iterations in solving linear system
			int it, 
			/// the maxmimum tolerable error
			double tol, 
			///vebose or quiet
			int verb);
    //@ManMemo: a constructor
    LSConjugateGradient(int n, LinearForward* p, Vector<double>* data, int it, double tol);
    ~LSConjugateGradient();
				
    //@ManMemo: CGLS search starting from m0, and returns an optimum model
    Model<double> 	optimizer(Model<double>& m0);
    //@ManMemo: CGLS search starting from m0, and returns an optimum model
    Model<long> 	optimizer(Model<long>& m0);
};

#endif