複素一般化線形最小二乗: GLM

LAPACKサンプルソースコード : 使用ルーチン名：ZGGGLM

ホーム > LAPACKサンプルプログラム目次 > 複素一般化線形最小二乗 > GLM

概要

本サンプルはFortran言語によりLAPACKルーチンZGGGLMを利用するサンプルプログラムです。

加重最小二乗問題を解きます。

\begin{displaymath} \min_{x} \left\Vert B^{-1} (d - A x) \right\Vert _{2} \end{displaymath}

ここで

\begin{displaymath} A = \left( \begin{array}{rrr} 0.96 - 0.81 i & -0.03 + 0.9... ...28 i & 0.20 - 0.12 i & -0.07 + 1.23 i \end{array} \right), \: \end{displaymath}


\begin{displaymath} B = \left( \begin{array}{cccc} 0.5 - 1.0 i & 0 & 0 & 0 \\... ... - 3.0 i & 0 \\ 0 & 0 & 0 & 5.0 - 4.0 i \end{array} \right) \end{displaymath}

及び

\begin{displaymath} d = \left( \begin{array}{r} 6.00 - 0.40 i \\ -5.27 + 0.90 i \\ 2.72 - 2.13 i \\ -1.30 - 2.80 i \end{array} \right). \end{displaymath}

入力データ

（本ルーチンの詳細はZGGGLM のマニュアルページを参照）
このデータをダウンロード 
ZGGGLM Example Program Data

  4             3             4                         :Values of M, N and P

( 0.96,-0.81) (-0.03, 0.96) (-0.91, 2.06)
(-0.98, 1.98) (-1.20, 0.19) (-0.66, 0.42)
( 0.62,-0.46) ( 1.01, 0.02) ( 0.63,-0.17)
( 1.08,-0.28) ( 0.20,-0.12) (-0.07, 1.23)               :End of matrix A

( 0.50,-1.00) ( 0.00, 0.00) ( 0.00, 0.00) ( 0.00, 0.00)
( 0.00, 0.00) ( 1.00,-2.00) ( 0.00, 0.00) ( 0.00, 0.00)
( 0.00, 0.00) ( 0.00, 0.00) ( 2.00,-3.00) ( 0.00, 0.00)
( 0.00, 0.00) ( 0.00, 0.00) ( 0.00, 0.00) ( 5.00,-4.00) :End of matrix B

( 6.00,-0.40)
(-5.27, 0.90)
( 2.72,-2.13)
(-1.30,-2.80)                                           :End of vector d

出力結果

（本ルーチンの詳細はZGGGLM のマニュアルページを参照）
この出力例をダウンロード 
 ZGGGLM Example Program Results

 Weighted least squares solution
 (  -0.9846,   1.9950) (   3.9929,  -4.9748) (  -3.0026,   0.9994)

 Residual vector
 ( 1.26E-04,-4.66E-04) ( 1.11E-03,-8.61E-04) ( 3.84E-03,-1.82E-03)
 ( 2.03E-03, 3.02E-03)

 Square root of the residual sum of squares
   5.79E-03

ソースコード

（本ルーチンの詳細はZGGGLM のマニュアルページを参照）

※本サンプルソースコードのご利用手順は「サンプルのコンパイル及び実行方法」をご参照下さい。
このソースコードをダウンロード 
    Program zggglm_example

!     ZGGGLM Example Program Text

!     Copyright 2017, Numerical Algorithms Group Ltd. http://www.nag.com

!     .. Use Statements ..
      Use blas_interfaces, Only: dznrm2
      Use lapack_interfaces, Only: zggglm
      Use lapack_precision, Only: dp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter :: nb = 64, nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=dp) :: rnorm
      Integer :: i, info, lda, ldb, lwork, m, n, p
!     .. Local Arrays ..
      Complex (Kind=dp), Allocatable :: a(:, :), b(:, :), d(:), work(:), x(:), &
        y(:)
!     .. Executable Statements ..
      Write (nout, *) 'ZGGGLM Example Program Results'
      Write (nout, *)
!     Skip heading in data file
      Read (nin, *)
      Read (nin, *) m, n, p
      lda = m
      ldb = m
      lwork = n + m + nb*(m+p)
      Allocate (a(lda,n), b(ldb,p), d(m), work(lwork), x(n), y(p))

!     Read A, B and D from data file

      Read (nin, *)(a(i,1:n), i=1, m)
      Read (nin, *)(b(i,1:p), i=1, m)
      Read (nin, *) d(1:m)

!     Solve the weighted least squares problem

!     minimize ||inv(B)*(d - A*x)|| (in the 2-norm)

      Call zggglm(m, n, p, a, lda, b, ldb, d, x, y, work, lwork, info)

!     Print least squares solution

      Write (nout, *) 'Weighted least squares solution'
      Write (nout, 100) x(1:n)

!     Print residual vector y = inv(B)*(d - A*x)

      Write (nout, *)
      Write (nout, *) 'Residual vector'
      Write (nout, 110) y(1:p)

!     Compute and print the square root of the residual sum of squares
      rnorm = dznrm2(p, y, 1)

      Write (nout, *)
      Write (nout, *) 'Square root of the residual sum of squares'
      Write (nout, 120) rnorm

100   Format (3(' (',F9.4,',',F9.4,')',:))
110   Format (3(' (',1P,E9.2,',',1P,E9.2,')',:))
120   Format (1X, 1P, E10.2)
    End Program


ご案内
関連情報
Privacy Policy  /  Trademarks