CMinpack使用介紹

github: https://github.com/devernay/cminpack
主頁: http://devernay.github.io/cminpack/
使用手冊: http://devernay.github.io/cminpack/man.htmlhtml

## LMDIF使用說明

官方英文介紹:http://devernay.github.io/cminpack/lmdif_.htmlc++

包括函數名

lmdif,lmdif1_ - 最小化非線性函數平方和git

函數概要

include <minpack.h>
void lmdif1_(void (*fcn)(int *m, int *n, double *x, double *fvec, int *iflag),
    int *m, int *n, double *x, double *fvec,
    double *tol, int *info, int *iwa, double *wa, int *lwa);


void lmdif_(void (*fcn)(int *m, int *n, double *x,  double *fvec, int *iflag),
            int *m, int *n, double *x, double *fvec,
            double *ftol, double *xtol, double *gtol, int *maxfev, double *epsfcn, double *diag, 
            int *mode, double *factor, int *nprint, int *info, int *nfev, double *fjac,
            int *ldfjac, int *ipvt, double *qtf,
            double *wa1, double *wa2, double *wa3, double *wa4 );

詳細描述

lmdif_的目的是最小化m個n元非線性方程的平方和,使用的方法是LM算法的改進。用戶須要提供計算方程的子程序。Jacobian矩陣會經過一個前向差分(forward-difference)近似計算獲得。
lmdif1_是相同的目的,可是調用方法更簡單一些。github

語言備註

這些函數是經過FORTRAN寫的,若是從C調用,須要記住如下幾點:算法

  • 名稱重編:
    • 2.95/3.0版本的g77下,全部函數如下劃線結尾,後續版本可能會更改;
  • 使用g77編譯:
    • 即便你的程序所有用C語言寫成,你也須要使用gcc進行連接,由於這樣它會自動導入FORTRAN庫。只使用g77進行編譯是最方便的(它處理C語言也是OK的);
  • 經過引用調用:
    • 全部函數參數都是指針;
  • 列優先數組:
    • z( i , j ) = z[ ( i - 1 ) + ( j - 1 ) * n
    • fcn是用戶提供用於計算函數的子程序。在C語言當中fcn須要以下定義:
    void fcn(int m, int n, double *x, double *fvec, int *iflag) {
      /* 計算函數在x點的值,經過fvec返回。*/
    }
    iflag的值不能被fcn所修改,除非用戶想要終止lmdif/lmdif1_。在這個例子中iflag設置爲負整數。

lmdif_和lmdif1_的共同參數

m:函數個數;
n:變量個數(n<=m)
x:長度爲n的數組,設置爲初始的估計解向量。輸出的時候x內容爲最終估計的解向量。
fvec:輸出長度爲m的數組,內容爲最終輸出x計算獲得的函數解。數組

lmdif1_的參數

tol:做爲輸入,非負數。用於函數終止的條件判斷:app

  • 平方和小於tol;
  • x之間的相對偏差小於tol;

info:做爲輸出。若是用戶終止了函數的執行,info將被設置爲iflag的值(負數)(詳細見fcn的描述),不然,info的值以下幾種狀況:less

  • 0:輸入參數不合適;
  • 1:平方和的相對偏差小於tol;
  • 2:x之間的相對偏差小於tol;
  • 3:1/2兩種狀況同時符合;
  • 4: fvec is orthogonal to the columns of the Jacobian to machine precision(這個狀況是什麼暫時不是很清楚)
  • 5:調用fcn的次數達到了200*(n+1)次;
  • 6:tol設置太小,平方和沒法達到那麼小;
  • 7:tol設置太小,x的近似解沒法優化到偏差達到那麼小。

iwa:長度n的工做數組;
wa:長度lwa的工做數組;
lwa:做爲輸入,整數,不能小於mn+5n+m;
[NOTE] 這三個輸入我也不知道做用,從樣例來看不須要初始化。ide

lmdif_的參數

暫時不用這部分,跳過。函數

官方樣例解讀

源碼: https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/examples/tlmdif1_.c

/*     driver for lmdif1 example. */
#include <stdio.h>
#include <math.h>
#include <assert.h>
#include <minpack.h>
#define real __minpack_real__

// 用戶自定義的函數f()
// real -> __cminpack_real__ -> 浮點數(double)
void fcn(const int *m, const int *n, const real *x, real *fvec, int *iflag);

int main()
{
  int j, m, n, info, lwa, iwa[3], one=1;
  real tol, fnorm, x[3], fvec[15], wa[75];

  // 函數個數15; 變量數3
  m = 15;
  n = 3;

  /* the following starting values provide a rough fit. */
  // 初始位置
  // 1.e0 = 1.0e0 = 1.0
  x[0] = 1.e0;
  x[1] = 1.e0;
  x[2] = 1.e0;

  // 爲何要設置75?
  lwa = 75;

  /* set tol to the square root of the machine precision.  unless high
     precision solutions are required, this is the recommended
     setting. */
  // (建議打印一下看值是多少)

  tol = sqrt(__minpack_func__(dpmpar)(&one));

  // 須要注意指針
  __minpack_func__(lmdif1)(&fcn, &m, &n, x, fvec, &tol, &info, iwa, wa, &lwa);

  // 最終的2範數(即平方和開根號)
  fnorm = __minpack_func__(enorm)(&m, fvec);
  printf("      final l2 norm of the residuals%15.7g\n\n", (double)fnorm);
  printf("      exit parameter                %10i\n\n", info);
  printf("      final approximate solution\n");
  for (j=1; j<=n; j++) {
    printf("%s%15.7g", j%3==1?"\n ":"", (double)x[j-1]);
  }
  printf("\n");
  return 0;
}

//        The problem is to determine the values of x(1), x(2), and x(3)
//        which provide the best fit (in the least squares sense) of
//              x(1) + u(i)/(v(i)*x(2) + w(i)*x(3)),  i = 1, 15
//        to the data
//              y = (0.14,0.18,0.22,0.25,0.29,0.32,0.35,0.39,
//                   0.37,0.58,0.73,0.96,1.34,2.10,4.39),
//        where u(i) = i, v(i) = 16 - i, and w(i) = min(u(i),v(i)).  The
//        i-th component of FVEC is thus defined by
//              y(i) - (x(1) + u(i)/(v(i)*x(2) + w(i)*x(3))).

void fcn(const int *m, const int *n, const real *x, real *fvec, int *iflag)
{
  /* function fcn for lmdif1 example */

  int i;
  real tmp1,tmp2,tmp3;

  // 實際的y值
  real y[15]={1.4e-1,1.8e-1,2.2e-1,2.5e-1,2.9e-1,3.2e-1,3.5e-1,3.9e-1,
              3.7e-1,5.8e-1,7.3e-1,9.6e-1,1.34e0,2.1e0,4.39e0};
  assert(*m == 15 && *n == 3);


  if (*iflag == 0) {
    /*      insert print statements here when nprint is positive. */
    /* if the nprint parameter to lmder is positive, the function is
       called every nprint iterations with iflag=0, so that the
       function may perform special operations, such as printing
       residuals. */
    // 這段沒有很看懂,在??狀況下打印信息
    return;
  }

  /* compute residuals */

  for (i=0; i<15; i++) {
    tmp1 = i+1;
    tmp2 = 15 - i;
    tmp3 = tmp1;
    if (i >= 8) {
      tmp3 = tmp2;
    }
    fvec[i] = y[i] - (x[0] + tmp1/(x[1]*tmp2 + x[2]*tmp3));
  }
}

[NOTE] 其餘內容有待更新

相關文章
相關標籤/搜索