itkContinuousLBFGSOptimizer.cxx

Go to the documentation of this file.

/*=========================================================================

  Program:   Insight Segmentation & Registration Toolkit
  Module:    $RCSfile: itkContinuousLBFGSOptimizer.cxx,v $
  Language:  C++
  Date:      $Date: 2005/03/15 02:19:46 $
  Version:   $Revision: 1.20 $

  Copyright (c) Insight Software Consortium. All rights reserved.
  See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.

     This software is distributed WITHOUT ANY WARRANTY; without even 
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 
     PURPOSE.  See the above copyright notices for more information.

=========================================================================*/
#ifndef _itkContinuousLBFGSOptimizer_txx
#define _itkContinuousLBFGSOptimizer_txx

#include "itkContinuousLBFGSOptimizer.h"

namespace itk
{

ContinuousLBFGSOptimizer
::ContinuousLBFGSOptimizer()
{
  m_OptimizerInitialized    = false;
  m_VnlOptimizer            = 0;
  m_ErrorOccurred           = false;
  m_Trace                              = false;
  m_MaximumNumberOfFunctionEvaluations = 2000;
  m_GradientConvergenceTolerance       = 1e-5;
  m_LineSearchAccuracy                 = 0.9;
  m_DefaultStepLength                  = 1.0;
}


ContinuousLBFGSOptimizer
::~ContinuousLBFGSOptimizer()
{
  delete m_VnlOptimizer;
}

void
ContinuousLBFGSOptimizer
::PrintSelf(std::ostream& os, Indent indent) const
{
  Superclass::PrintSelf( os, indent );
  os << indent << "Trace: ";
  if ( m_Trace ) { os << "On"; } else { os << "Off"; }
  os << std::endl;
  os << indent << "MaximumNumberOfFunctionEvaluations: ";
  os << m_MaximumNumberOfFunctionEvaluations << std::endl;
  os << indent << "GradientConvergenceTolerance: ";
  os << m_GradientConvergenceTolerance << std::endl;
  os << indent << "LineSearchAccuracy: ";
  os << m_LineSearchAccuracy << std::endl;
  os << indent << "DefaultStepLength: ";
  os << m_DefaultStepLength << std::endl;
  os << indent << "ErrorOccured: ";
  os << m_ErrorOccurred << std::endl;

}

void
ContinuousLBFGSOptimizer
::SetTrace( bool flag )
{
  if ( flag == m_Trace )
    {
    return;
    }

  m_Trace = flag;
  if ( m_OptimizerInitialized )
    {
    m_VnlOptimizer->set_trace( m_Trace );
    }

  this->Modified();
}

void
ContinuousLBFGSOptimizer
::SetMaximumNumberOfFunctionEvaluations( unsigned int n )
{
  if ( n == m_MaximumNumberOfFunctionEvaluations )
    {
    return;
    }

  m_MaximumNumberOfFunctionEvaluations = n;
  if ( m_OptimizerInitialized )
    {
    m_VnlOptimizer->set_max_function_evals(
      static_cast<int>( m_MaximumNumberOfFunctionEvaluations ) );
    }

  this->Modified();
}

void
ContinuousLBFGSOptimizer
::SetGradientConvergenceTolerance( double f )
{
  if ( f == m_GradientConvergenceTolerance )
    {
    return;
    }

  m_GradientConvergenceTolerance = f;
  if ( m_OptimizerInitialized )
    {
    m_VnlOptimizer->set_g_tolerance( m_GradientConvergenceTolerance );
    }

  this->Modified();
}

void
ContinuousLBFGSOptimizer
::SetLineSearchAccuracy( double f )
{
  if ( f == m_LineSearchAccuracy )
    {
    return;
    }

  m_LineSearchAccuracy = f;
  if ( m_OptimizerInitialized )
    {
    m_VnlOptimizer->line_search_accuracy = m_LineSearchAccuracy;
    }

  this->Modified();
}

void
ContinuousLBFGSOptimizer
::SetDefaultStepLength( double f )
{
  if ( f == m_DefaultStepLength )
    {
    return;
    }

  m_DefaultStepLength = f;
  if ( m_OptimizerInitialized )
    {
    m_VnlOptimizer->default_step_length = m_DefaultStepLength;
    }

  this->Modified();
}

ContinuousLBFGSOptimizer::MeasureType
ContinuousLBFGSOptimizer
::GetValue()
{
  return this->GetNonConstCostFunctionAdaptor()->f(this->GetCurrentPosition());
}

void
ContinuousLBFGSOptimizer
::SetCostFunction( SingleValuedCostFunction * costFunction )
{
  const unsigned int numberOfParameters = 
    costFunction->GetNumberOfParameters();

  CostFunctionAdaptorType * adaptor = 
    new CostFunctionAdaptorType( numberOfParameters );
       
  adaptor->SetCostFunction( costFunction );

  if( m_OptimizerInitialized )
    { 
    delete m_VnlOptimizer;
    }
    
  this->SetCostFunctionAdaptor( adaptor );

  m_VnlOptimizer = new vnl_lbfgs( *adaptor );

  // set the optimizer parameters
  m_VnlOptimizer->set_trace( m_Trace );
  m_VnlOptimizer->set_max_function_evals(
    static_cast<int>( m_MaximumNumberOfFunctionEvaluations ) );
  m_VnlOptimizer->set_g_tolerance( m_GradientConvergenceTolerance );
  m_VnlOptimizer->line_search_accuracy = m_LineSearchAccuracy;
  m_VnlOptimizer->default_step_length  = m_DefaultStepLength;

  m_OptimizerInitialized = true;

}



void
ContinuousLBFGSOptimizer
::StartOptimization( void )
{
  
  this->InvokeEvent( StartEvent() );

  if( this->GetMaximize() )
    {
    this->GetNonConstCostFunctionAdaptor()->NegateCostFunctionOn();
    }

  ParametersType initialPosition = this->GetInitialPosition();
  ParametersType parameters( initialPosition );  

  // If the user provides the scales then we set otherwise we don't
  // for computation speed.
  // We also scale the initial parameters up if scales are defined.
  // This compensates for later scaling them down in the cost function adaptor
  // and at the end of this function.  
  if(m_ScalesInitialized)
    {
    ScalesType scales = this->GetScales();
    this->GetNonConstCostFunctionAdaptor()->SetScales(scales);
    for(unsigned int i=0;i<parameters.size();i++)
      {
      parameters[i] *= scales[i]; 
      }
    }
  
  m_ErrorOccurred = false;

  // vnl optimizers return the solution by reference 
  // in the variable provided as initial position
  m_VnlOptimizer->minimize( parameters );

  if ( parameters.size() != initialPosition.size() )
    { //there was an error in the optimization process
      m_ErrorOccurred = true;
      std::cerr << std::endl << "Error occured on in lbfgs optimizer. Continued with initial position." << std::endl;
      parameters = initialPosition;
    }
  else if(m_ScalesInitialized)
    { // we scale the parameters down if scales are defined
    ScalesType scales = this->GetScales();
    for(unsigned int i=0;i<parameters.size();i++)
      {
      parameters[i] /= scales[i]; 
      }
    }

  this->SetCurrentPosition( parameters );       
  
  this->InvokeEvent( EndEvent() );
}




vnl_lbfgs * 
ContinuousLBFGSOptimizer
::GetOptimizer()
{
  return m_VnlOptimizer;
}




} // end namespace itk

#endif

---