nmrFnJacobianSolver.h

Go to the documentation of this file.

/* -*- Mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-    */
/* ex: set filetype=cpp softtabstop=4 shiftwidth=4 tabstop=4 cindent expandtab: */

/*
  
  Author(s):    Ankur Kapoor
  Created on: 2004-10-26

  (C) Copyright 2004-2007 Johns Hopkins University (JHU), All Rights
  Reserved.

--- begin cisst license - do not edit ---

This software is provided "as is" under an open source license, with
no warranty.  The complete license can be found in license.txt and
http://www.cisst.org/cisst/license.txt.

--- end cisst license ---
*/


#ifndef _nmrFnJacobianSolver_h
#define _nmrFnJacobianSolver_h


#include <cisstVector/vctDynamicMatrix.h>

#include <cnetlib.h>

class nmrFnJacobianSolver {
    // we have this class so that we reserve memory only one
    // would help if svd of a same size matrix (or a matrix)
    // that doesnt change much is desired.

protected:
    CISSTNETLIB_INTEGER N;
    CISSTNETLIB_INTEGER Ldfjac;
    CISSTNETLIB_DOUBLE Tolerance;
    CISSTNETLIB_INTEGER Info;
    CISSTNETLIB_INTEGER Lwork;
    vctDynamicVector<CISSTNETLIB_DOUBLE> Work;

public:
    nmrFnJacobianSolver(void)
    {
        Allocate(0);
    }


    nmrFnJacobianSolver(CISSTNETLIB_INTEGER n)
    {
        Allocate(n);
    }

    
    nmrFnJacobianSolver(vctDynamicVector<CISSTNETLIB_DOUBLE> &X) {
        Allocate(X);
    }


    inline void Allocate(CISSTNETLIB_INTEGER n) {
        N = n;
    Lwork = N*(N+13)/2;
    Ldfjac = N;
    Tolerance = 1e-6;
        Work.SetSize(Lwork);
    }


    inline void Allocate(vctDynamicVector<CISSTNETLIB_DOUBLE> &X) {
        Allocate(X.size());
    }


    template <int __instanceLineF, class __elementTypeF, int __instanceLineJ, class __elementTypeJ>
    inline void Solve(nmrCallBackFunctionF<__instanceLineF, __elementTypeF> &callBackF,
                      nmrCallBackFunctionJ<__instanceLineJ, __elementTypeJ> &callBackJ,
                      vctDynamicVector<CISSTNETLIB_DOUBLE> &X, vctDynamicVector<CISSTNETLIB_DOUBLE> &F,
                      vctDynamicVector<CISSTNETLIB_DOUBLE> &J, CISSTNETLIB_DOUBLE tolerance) throw (std::runtime_error) {
        if ((N != static_cast<CISSTNETLIB_INTEGER>(X.size())) || (N != static_cast<CISSTNETLIB_INTEGER>(F.size())) || (Ldfjac*N != static_cast<CISSTNETLIB_INTEGER>(J.size()))) {
            cmnThrow(std::runtime_error("nmrFnJacobianSolver Solve: Size used for Allocate was different"));
        }
        Tolerance = tolerance;
        /* call the MINPACK C function */
        hybrj2_(callBackF.FunctionFhybrd, callBackJ.FunctionFhybrj2, &N,
                X.Pointer(), F.Pointer(), J.Pointer(),
                &Ldfjac, &Tolerance, &Info,
                Work.Pointer(), &Lwork);
    }
    
};

#endif // _nmrFnJacobianSolver_h