py_llsq_obj

View page source

Python Linear Least Squares Objective

Syntax

      llsq = cmpad.llsq_obj ( like_numpy )
      llsq . setup ( option )
      obj = llsq ( x )

Prototype

class llsq_obj :
   def __init__(self, like_numpy) :
      self.like_numpy = like_numpy
   def domain(self) :
      return self.n_arg
   def range (self) :
      return 1
   #
   def setup(self, option) :
      assert type(option) == dict
      assert type(option['n_arg']) == int
      assert type(option['n_other']) == int
      assert option['n_arg'] > 0
      assert option['n_other'] > 0

Algorithm

This is a python implementation of the llsq_obj algorithm . It is vectorized using arrays that act like numpy arrays. Different array types are used by the different AD packages.

like_numpy

This is a like_numpy class. It is used to vectorize this algorithm.

n_arg

see n_arg .

n_other

see n_other . This is the number of elements that are computed by one like_numpy operation of the llsq_obj algorithm.

Example

xam_llsq_obj.py contains an example and test of llsq_obj .

Source Code

The code below is the implementation of this function:

import cmpad
import numpy
#
# BEGIN PROTOTYPE
class llsq_obj :
   def __init__(self, like_numpy) :
      self.like_numpy = like_numpy
   def domain(self) :
      return self.n_arg
   def range (self) :
      return 1
   #
   def setup(self, option) :
      assert type(option) == dict
      assert type(option['n_arg']) == int
      assert type(option['n_other']) == int
      assert option['n_arg'] > 0
      assert option['n_other'] > 0
      # END PROTOTYPE
      #
      # n_arg, n_other
      n_arg   = option['n_arg']
      n_other = option['n_other']
      #
      # t
      if n_other == 1 :
         t = numpy.array( [ 0.0 ] )
      else :
         t = numpy.linspace(-1.0, 1.0, n_other)
      #
      # q
      q = -1.0 * (t < 0) + 1.0 * (t > 0)
      #
      # self
      self.n_arg = n_arg
      self.t     = t
      self.q     = q
   #
   def __call__(self, ax) :
      assert len(ax) == self.n_arg
      #
      # like_numpy
      like_numpy = self.like_numpy
      #
      # model
      model = like_numpy.array( [ 0.0 ] )
      ti    = numpy.ones( len(self.t) )
      for i in range(self.n_arg) :
         model = model + like_numpy.array( ax[i] ) * like_numpy.array(ti)
         ti    = ti * self.t
      #
      # squared_residual
      residual         = model - like_numpy.array( self.q )
      squared_residual = residual * residual
      #
      # objective
      objective = like_numpy.array( 0.5 * squared_residual.sum() )
      objective = objective.reshape(-1)
      #
      return objective