itkfemlinearsystemwrapperitpack.h

DTMK软件开发包,此为开源软件,是一款很好的医学图像开发资源.

   */
  void   SetEstimatedSpectralRadiusSSOR(double i) {   m_RPARM[5] = i; }

  /**
   * Get the ??
   */
  double GetEstimatedSpectralRadiusSSOR()    { return m_RPARM[5]; }

  /**
   * Set the ??
   * \param i ??
   */
  void   SetEstimatedSpectralRadiusLU(double i) {   m_RPARM[6] = i; }

  /**
   * Get the ??
   */
  double GetEstimatedSpectralRadiusLU()    { return m_RPARM[6]; }

  /**
   * Set the tolerance level
   * \param i tolerance
   */
  void   SetTolerance(double i) {   m_RPARM[7] = i; }

  /**
   * Get the tolerance level
   */
  double GetTolerance()    { return m_RPARM[7]; }

  /**
   * Set the time to convergence
   * \param i ??
   */
  void   SetTimeToConvergence(double i) {   m_RPARM[8] = i; }

  /**
   * Get the time to convergence
   */
  double GetTimeToConvergence()    { return m_RPARM[8]; }

  /**
   * Set the time for call
   * \param i ??
   */
  void   SetTimeForCall(double i) {   m_RPARM[9] = i; }

  /**
   * Get the time for call
   */
  double GetTimeForCall()    { return m_RPARM[9]; }

  /**
   * Set digits in error
   * \param i number of digits in error
   */
  void   SetDigitsInError(double i) {   m_RPARM[10] = i; }

  /**
   * Get the number of digits in the error
   */
  double GetDigitsInError()    { return m_RPARM[10]; }

  /**
   * Set the number of digits in the residual
   * \param i number of digits in the residual
   */
  void   SetDigitsInResidual(double i) {   m_RPARM[11] = i; }

  /**
   * Get the number of digits in the residual
   */
  double GetDigitsInResidual()    { return m_RPARM[11]; }

  /** 
   * Set numerical solving method to jacobian conjugate gradient 
   */
  void JacobianConjugateGradient() { m_Method = 0; }

  /** 
   * Set numerical solving method to jacobian semi iterative
   */
  void JacobianSemiIterative() { m_Method = 1; }

  /** 
   * Set numerical solving method to successive over-relaxation
   */
  void SuccessiveOverrelaxation() { m_Method = 2; }

  /** 
   * Set numerical solving method to symmetric successive over-relaxation
   * conjugate gradient 
   */
  void SymmetricSuccessiveOverrelaxationConjugateGradient() { m_Method = 3; }

  /** 
   * Set numerical solving method to symmetric successive over-relaxation
   * successive over-relaxation
   */
  void SymmetricSuccessiveOverrelaxationSuccessiveOverrelaxation() { m_Method = 4; }

  /** 
   * Set numerical solving method to reduced system conjugate gradient 
   */
  void ReducedSystemConjugateGradient() { m_Method = 5; }

  /** 
   * Set numerical solving method to reduced system semi-iteration */
  void ReducedSystemSemiIteration() { m_Method = 6; }


  /* -----------------------------------------------------------------
   * 
   * Redefine methods defined in LinearSystemWrapper 
   *
   * ----------------------------------------------------------------- 
   */

  /**
   * set maximum number of entires in a matrix
   * \param maxNonZeroValues maximum number of entries allowed in matrix
   * \note this must be called before any matrices are initialized
   */  
  virtual void SetMaximumNonZeroValuesInMatrix(unsigned int maxNonZeroValues) {m_MaximumNonZeroValues = maxNonZeroValues;}


  void ScaleMatrix(Float scale, unsigned int matrixIndex);


  /* -----------------------------------------------------------------
   * 
   * Functions required by LinearSystemWrapper 
   *
   * ----------------------------------------------------------------- 
   */

  /**
   * constructor 
   */
  LinearSystemWrapperItpack();
  
  /**
   * destructor
   */
  ~LinearSystemWrapperItpack();


  /* memory management routines */
  virtual void  InitializeMatrix(unsigned int matrixIndex);  

  virtual bool  IsMatrixInitialized(unsigned int matrixIndex);

  virtual void  DestroyMatrix(unsigned int matrixIndex);

  virtual void  InitializeVector(unsigned int vectorIndex);

  virtual bool  IsVectorInitialized(unsigned int vectorIndex);

  virtual void  DestroyVector(unsigned int vectorIndex);

  virtual void  InitializeSolution(unsigned int solutionIndex);

  virtual bool  IsSolutionInitialized(unsigned int solutionIndex);

  virtual void  DestroySolution(unsigned int solutionIndex);

  /* assembly & solving routines */
  virtual Float GetMatrixValue(unsigned int i, unsigned int j, unsigned int matrixIndex) const;

  virtual void  SetMatrixValue(unsigned int i, unsigned int j, Float value, unsigned int matrixIndex);

  virtual void  AddMatrixValue(unsigned int i, unsigned int j, Float value, unsigned int matrixIndex);

  virtual void GetColumnsOfNonZeroMatrixElementsInRow( unsigned int row, ColumnArray& cols, unsigned int matrixIndex );

  virtual Float GetVectorValue(unsigned int i, unsigned int vectorIndex) const;

  virtual void  SetVectorValue(unsigned int i, Float value, unsigned int vectorIndex);

  virtual void  AddVectorValue(unsigned int i, Float value, unsigned int vectorIndex);

  virtual Float GetSolutionValue(unsigned int i, unsigned int solutionIndex) const;

  virtual void  SetSolutionValue(unsigned int i, Float value, unsigned int solutionIndex);

  virtual void  AddSolutionValue(unsigned int i, Float value, unsigned int solutionIndex);

  virtual void  Solve(void);


  /* matrix & vector manipulation routines */
  virtual void  SwapMatrices(unsigned int matrixIndex1, unsigned int matrixIndex2);

  virtual void  SwapVectors(unsigned int vectorIndex1, unsigned int vectorIndex2);

  virtual void  SwapSolutions(unsigned int solutionIndex1, unsigned int solutionIndex2);

  virtual void  CopySolution2Vector(unsigned solutionIndex, unsigned int vectorIndex);

  virtual void  CopyVector2Solution(unsigned int vectorIndex, unsigned int solutionIndex);

  virtual void  MultiplyMatrixMatrix(unsigned int resultMatrixIndex, unsigned int leftMatrixIndex, unsigned int rightMatrixIndex);

  virtual void  MultiplyMatrixVector(unsigned int resultVectorIndex, unsigned int matrixIndex, unsigned int vectorIndex);

private:

  /** pointer to vector of matrices */
  MatrixHolder *m_Matrices;

  /** pointer to vector of force arrays */
  VectorHolder *m_Vectors;

  /** pointer to vector of solution arrays */
  VectorHolder *m_Solutions;

  /** pointer to array of unsigned int's indicating max number of entries in each matrix */
  //UnsignedIntegerArrayPtr m_MaximumNonZeroValues;
  unsigned int m_MaximumNonZeroValues;

  /** Array of pointers to available solver functions */
  ItkItpackSolverFunction    m_Methods[7];

  /** flag indicating which solver function should be used */
  integer m_Method;

  /** vector of length 12 used to initialize various parameters on input */
  integer m_IPARM[12];

  /** vector of length 12 used to initialize various parameters on input */
  doublereal m_RPARM[12];

};

/**
 * \class FEMExceptionItpackSolver
 * \brief handles errors that occur in itpack solving routines
 * \sa LinearSystemWrapperItpack
 * \sa FEMException
 */
class FEMExceptionItpackSolver : public FEMException
{
public:
  /** typedefs from f2c.h  */
  typedef long      integer;

  /**
   * Constructor. In order to construct this exception object, four parameters
   * must be provided: file, lineNumber, location and a detailed description
   * of the exception.
   */
  FEMExceptionItpackSolver(const char *file, unsigned int lineNumber, std::string location, integer errorCode);
 
  /** Virtual destructor needed for subclasses. Has to have empty throw(). */
  virtual ~FEMExceptionItpackSolver() throw() {}
  
  /** Type related information. */
  itkTypeMacro(FEMExceptionItpackSolver,FEMException);
  
};



}} // end namespace itk::fem

#endif // #ifndef __itkFEMLinearSystemWrapperItpack_h